refactor: change pf alg

rspower/lib/dsbus_dv.txt

@@ -102,5 +102,5 @@ fn dSbus_dV(Ybus, V, vcart) {
         dSbus_dV1 = c(0,1) * diagV * conj(diagIbus - Ybus * diagV);     // dSbus/dVa
         dSbus_dV2 = diagV * conj(Ybus * diagVnorm) + conj(diagIbus) * diagVnorm;    // dSbus/dVm
     }
-    return horzcat(dSbus_dV1, dSbus_dV2);
+    return dSbus_dV1+dSbus_dV2;
 }
\ No newline at end of file

rspower/lib/make_jac.txt

@@ -17,7 +17,5 @@ fn make_jac(baseMVA, bus, branch, gen) {
     gen_vg = get_multi(slice(gen, [0], [VG-1,VG]), on_k);
     V = set(V, index, gen_vg ./ abs(voltage_g) .* voltage_g);
     dsbus_dv = dSbus_dV(Ybus, V, 0.);
-    j_real = real(dsbus_dv);
-    j_imag = imag(dsbus_dv);
-    return vertcat(j_real, j_imag);
+    return dsbus_dv;
 }
\ No newline at end of file

rspower/lib/newtonpf.txt

@@ -53,6 +53,21 @@ fn newtonpf(Ybus, Sbus, V0, ref, pv, pq, mpopt) {
     j5 = j4 + 1;    j6 = j4 + npq;      // j5:j6 - V mag of pq buses
 
     // evaluate F(x0)
-    mis = V .* cosnj(Ybus * V) - Sbus(Vm);
-    return V;
+    mis = V .* conj(Ybus * V);
+    F = mis;
+    // check tolerance
+    normF = norm_max(F);
+    // do Newton iterations
+    while ~~converged && i < max_it {
+        // update iteration counter
+        i = i + 1;
+        // evaluate Jacobian
+        J = dSbus_dV(Ybus, V);
+        // compute update step
+        dx = linsolve(J, -F);
+        if normF < tol {
+            converged = 1;
+        }
+    }
+    return dx;
 }
\ No newline at end of file

rspower/lib/runpf.txt

@@ -37,7 +37,7 @@ fn runpf() {
     loop {
         // compute bus power injections
         Sbus = make_sbus(baseMVA, bus, gen, 1, 1, 1);
-        V = newtonpf(Ybus, Sbus, V_init, ref, pv, pq, [1e-8, 100, 0]);
+        V = newtonpf(Ybus, Sbus, V_init, ref, pv, pq, [1e-6, 100]);
         if ~~is_empty(V) && qlim {
             // if V is empty, it means the power flow did not converge
             break;