Fixing bug on PTDF calculation with distributed slack

src/ptdf_calculations.jl

@@ -153,10 +153,10 @@ function _calculate_PTDF_matrix_KLU(
         @info "Distributed bus"
         copyto!(PTDFm_t, BA)
         PTDFm_t[valid_ix, :] = KLU.solve!(K, PTDFm_t[valid_ix, :])
-        PTDFm_t[ref_bus_positions, :] .= 0.0
+        PTDFm_t[collect(ref_bus_positions), :] .= 0.0
         slack_array = dist_slack / sum(dist_slack)
         slack_array = reshape(slack_array, 1, buscount)
-        return PTDFm_t - (PTDFm_t * slack_array) * ones(1, buscount)
+        return PTDFm_t - ones(buscount, 1) * (slack_array * PTDFm_t)
     else
         error("Distributed bus specification doesn't match the number of buses.")
     end
@@ -246,10 +246,16 @@ function _calculate_PTDF_matrix_DENSE(
         return PTDFm_t
     elseif length(dist_slack) == buscount
         @info "Distributed bus"
+        PTDFm_t[valid_ixs, :] = gemm(
+            'N',
+            'N',
+            getri!(ABA, bipiv),
+            BA[valid_ixs, :],
+        )
         slack_array = dist_slack / sum(dist_slack)
-        slack_array = reshape(slack_array, buscount, 1)
+        slack_array = reshape(slack_array, 1, buscount)
         return PTDFm_t -
-               gemm('N', 'N', gemm('N', 'N', PTDFm, slack_array), ones(1, buscount))
+               gemm('N', 'N', ones(buscount, 1), gemm('N', 'N', slack_array, PTDFm_t))
     else
         error("Distributed bus specification doesn't match the number of buses.")
     end

test/test_ptdf.jl

@@ -192,3 +192,22 @@ end
     # Test Throw error when isolated buses are connected to an available branch
     @test_throws IS.ConflictingInputsError ptdf_2 = PTDF(sys_2)
 end
+
+@test "PTDF matrices with distributed slack" begin
+    sys5 = PSB.build_system(PSB.PSITestSystems, "c_sys5")
+
+    bus_number = length(PNM.get_buses(sys))
+
+    dist_slack = 1/bus_number*ones(bus_number)
+    slack_array = dist_slack / sum(dist_slack)
+    slack_array = reshape(slack_array, 1, buscount)
+
+    P5_1 = PTDF(sys5; slack_array=slack_array, linear_solver = "KLU")
+    P5_2 = PTDF(sys5; slack_array=slack_array, linear_solver = "Dense")
+    P5_3 = PTDF(sys5; slack_array=slack_array, linear_solver = "MKLPardiso")
+
+    @assert isapprox(P5_1.data, P5_2.data, atol=1e-5)
+    @assert isapprox(P5_1.data, P5_3.data, atol=1e-5)
+    @assert isapprox(P5_2.data, P5_3.data, atol=1e-5)
+
+end