Todo added

PSD/complex_current_and_nodes.py

@@ -38,6 +38,33 @@ class SwingEquationNode(NodeType):
             index=j)
 
 
+class DroopEquationNode(NodeType):
+    def __init__(self, number_of_variables, v_set = 1., infeed = 1., damping = 0.1, H = 1.5):
+        super(DroopEquationNode, self).__init__(number_of_variables)
+        self.v_set = v_set
+        self.v_set_squared = v_set ** 2
+        self.infeed = infeed
+        self.damping = damping
+        self.H_inv = 1./H
+
+    def node_dynamics(self, v, omega, i, t, dv, domega, index):
+        # TODO Jonathan: Implement the Schiffer equations here:
+        # The current is assumed to be negative when it is flowing away from the node
+        dv[index] = 1.j * omega[index] * v[index] - (v[index].real * v[index].real + v[index].imag * v[index].imag - self.v_set_squared) * v[index]
+        domega[index] = self.H_inv * (self.infeed - self.damping * omega[index] - (v[index] * i[index].conjugate()).real)
+
+    def node_dynamics_string(self, j="{index}"):
+        return """
+    dv[{index}] = 1.j * omega[{index}] * v[{index}] - (v[{index}].real * v[{index}].real + v[{index}].imag * v[{index}].imag - {v_set_squared}) * v[{index}]
+    domega[{index}] = {H_inv:.64f} * ({infeed} - {damping} * omega[{index}] - (v[{index}] * i[{index}].conjugate()).real)"""\
+            .format(
+            v_set_squared=self.v_set_squared,
+            H_inv=self.H_inv,
+            infeed=self.infeed,
+            damping=self.damping,
+            index=j)
+
+
 @njit(complex128[:](complex128[:], int32[:], int32[:], complex128[:]))
 def numba_sp_dot(data, indptr, indices, v):
     res = np.zeros(len(indptr) - 1, dtype=np.complex128)
@@ -48,6 +75,7 @@ def numba_sp_dot(data, indptr, indices, v):
             index += 1
     return res
 
+
 @njit
 def numba_sp_complex_currents(data, indptr, indices, v, i):
     index = 0
@@ -132,7 +160,7 @@ if __name__ == "__main__":
 
     node_list = list()
     node_list.append(SwingEquationNode(3, infeed=1.))
-    node_list.append(SwingEquationNode(3, infeed=-1.))
+    node_list.append(SwingEquationNode(3, infeed=-1))
 
     Y = -8.j * np.ones((2, 2), dtype=np.complex128)
     Y[0, 0] *= -1.
@@ -165,6 +193,9 @@ if __name__ == "__main__":
 
     from scipy.integrate import odeint
 
+    node_list[0].infeed += 0.1
+    rhs = define_network_rhs(node_list, Y)
+
     states = odeint(rhs, ic, times)
 
     import matplotlib.pyplot as plt