Skip to content

Commit

Permalink
tests passing for symmetric case
Browse files Browse the repository at this point in the history
  • Loading branch information
@leovsch
leovsch committed Nov 6, 2024
1 parent 135e882 commit 1722b64
Show file tree
Hide file tree
Showing 2 changed files with 49 additions and 198 deletions.
22 changes: 7 additions & 15 deletions power_grid_model_c/power_grid_model/include/power_grid_model/component/asym_line.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -36,7 +36,6 @@ class AsymLine : public Branch {
}

y_series = inv(z_series);
std::cout << "AsymLine class y_series: " << std::endl << y_series << std::endl;
y_shunt = 2 * pi * system_frequency * c_matrix * 1.0i;
}

Expand Down Expand Up @@ -65,29 +64,23 @@ class AsymLine : public Branch {
}

ComplexTensor<asymmetric_t> compute_z_series_from_input(const power_grid_model::AsymLineInput& asym_line_input) {
ComplexTensor<asymmetric_t> z_series;
ComplexTensor<asymmetric_t> z_series_abc;
if (is_nan(asym_line_input.r_na) && is_nan(asym_line_input.x_na)) {
ComplexTensor<asymmetric_t> r_matrix = ComplexTensor<asymmetric_t>(asym_line_input.r_aa, asym_line_input.r_ba, asym_line_input.r_bb, asym_line_input.r_ca, asym_line_input.r_cb, asym_line_input.r_cc);
ComplexTensor<asymmetric_t> x_matrix = ComplexTensor<asymmetric_t>(asym_line_input.x_aa, asym_line_input.x_ba, asym_line_input.x_bb, asym_line_input.x_ca, asym_line_input.x_cb, asym_line_input.x_cc);
z_series = r_matrix + x_matrix * 1.0i;
z_series_abc = r_matrix + x_matrix * 1.0i;
}
else {
ComplexTensor4 r_matrix = ComplexTensor4(asym_line_input.r_aa, asym_line_input.r_ba, asym_line_input.r_bb, asym_line_input.r_ca, asym_line_input.r_cb, asym_line_input.r_cc, asym_line_input.r_na, asym_line_input.r_nb, asym_line_input.r_nc, asym_line_input.r_nn);
ComplexTensor4 x_matrix = ComplexTensor4(asym_line_input.x_aa, asym_line_input.x_ba, asym_line_input.x_bb, asym_line_input.x_ca, asym_line_input.x_cb, asym_line_input.x_cc, asym_line_input.x_na, asym_line_input.x_nb, asym_line_input.x_nc, asym_line_input.x_nn);

std::cout << "AsymLine class r_matrix: " << std::endl << r_matrix << std::endl;
std::cout << "AsymLine class x_matrix: " << std::endl << x_matrix << std::endl;
ComplexTensor4 y = r_matrix + 1.0i * x_matrix;
std::cout << "AsymLine class y_matrix: " << std::endl << y << std::endl;

ComplexTensor<asymmetric_t> Y_reduced = this->kron_reduction(y);
std::cout << "AsymLine class Y_reduced: " << std::endl << Y_reduced << std::endl;

DoubleComplex a = std::pow(e, 1.0i * (2.0 / 3.0) * pi);
ComplexTensor<asymmetric_t> a_matrix = ComplexTensor<asymmetric_t>(1, 1, pow(a, 2), 1, a, pow(a, 2));
ComplexTensor<asymmetric_t> a_matrix_inv = a_matrix.inverse();
z_series = a_matrix_inv * Y_reduced * a_matrix;
z_series_abc = kron_reduction(y);
}
DoubleComplex a = std::pow(e, 1.0i * (2.0 / 3.0) * pi);
ComplexTensor<asymmetric_t> a_matrix = ComplexTensor<asymmetric_t>(1, 1, pow(a, 2), 1, a, pow(a, 2));
ComplexTensor<asymmetric_t> a_matrix_inv = (1.0/3.0) * ComplexTensor<asymmetric_t>(1, 1, a, 1, pow(a, 2), a);
ComplexTensor<asymmetric_t> z_series = (a_matrix_inv.matrix() * z_series_abc.matrix() * a_matrix.matrix()).array();
return z_series;
}

Expand Down Expand Up @@ -117,7 +110,6 @@ class AsymLine : public Branch {
BranchCalcParam<symmetric_t> sym_calc_param() const override {
DoubleComplex y1_series_ = average_of_diagonal_of_matrix(y_series) - average_of_off_diagonal_of_matrix(y_series);
DoubleComplex y1_shunt_ = average_of_diagonal_of_matrix(y_shunt) - average_of_off_diagonal_of_matrix(y_shunt);
std::cout << "AsymLine class: " << y1_series_ << "," << y1_shunt_ << std::endl;;
return calc_param_y_sym(y1_series_, y1_shunt_, 1.0);
}

Expand Down
225 changes: 42 additions & 183 deletions tests/cpp_unit_tests/test_asym_line.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -45,16 +45,13 @@ TEST_CASE("Test asym line") {
double const base_y = base_i * base_i / base_power_1p;
Branch& branch = asym_line;
ComplexTensor<asymmetric_t> const y_series = ComplexTensor<asymmetric_t>(0.66303418-0.34266364i, -0.02971114+0.03783535i, 0.04762194+0.00681293i, 0.04762194+0.00681293i, 2.48612768-0.46271628i, 0.05942228-0.07567069i, -0.02971114+0.03783535i, -0.09524388-0.01362585i, 2.48612768-0.46271628i);
std::cout << "AsymLine test y_series: " << std::endl << y_series << std::endl;

ComplexTensor<asymmetric_t> const y_shunt = 2 * pi * system_frequency * ComplexTensor<asymmetric_t>(input.c0 + input.c1, -input.c1) * 1.0i;

DoubleComplex const y1_series = (y_series(0,0) + y_series(1,1) + y_series(2,2)) / 3.0 - (y_series(0,2) + y_series(1,1) + y_series(2,0)) / 3.0;
DoubleComplex const y1_shunt = (y_shunt(0,0) + y_shunt(1,1) + y_shunt(2,2)) / 3.0 - (y_shunt(0,2) + y_shunt(1,1) + y_shunt(2,0)) / 3.0;

// symmetric
std::cout << "AsymLine test: " << y1_series << "," << y1_shunt << std::endl;

DoubleComplex const yff1 = y1_series + 0.5 * y1_shunt;
DoubleComplex const yft1 = -y1_series;
DoubleComplex const ys1 = 0.5 * y1_shunt + 1.0 / (1.0 / y1_series + 2.0 / y1_shunt);
Expand Down Expand Up @@ -109,190 +106,52 @@ TEST_CASE("Test asym line") {
CHECK(cabs(param.ytt() - yff1) < numerical_tolerance);
CHECK(cabs(param.ytf() - yft1) < numerical_tolerance);
CHECK(cabs(param.yft() - yft1) < numerical_tolerance);
// // to connected
// CHECK(branch.update(BranchUpdate{1, false, na_IntS}).topo);
// param = branch.calc_param<symmetric_t>();
// CHECK(cabs(param.yff() - 0.0) < numerical_tolerance);
// CHECK(cabs(param.ytt() - ys1) < numerical_tolerance);
// CHECK(cabs(param.ytf() - 0.0) < numerical_tolerance);
// CHECK(cabs(param.yft() - 0.0) < numerical_tolerance);
// // not connected
// CHECK(branch.set_status(na_IntS, false));
// param = branch.calc_param<symmetric_t>();
// CHECK(cabs(param.yff() - 0.0) < numerical_tolerance);
// CHECK(cabs(param.ytt() - 0.0) < numerical_tolerance);
// CHECK(cabs(param.ytf() - 0.0) < numerical_tolerance);
// CHECK(cabs(param.yft() - 0.0) < numerical_tolerance);
// // not changing
// CHECK(!branch.set_status(false, false));
// // from connected
// CHECK(branch.set_status(true, na_IntS));
// param = branch.calc_param<symmetric_t>();
// CHECK(cabs(param.yff() - ys1) < numerical_tolerance);
// CHECK(cabs(param.ytt() - 0.0) < numerical_tolerance);
// CHECK(cabs(param.ytf() - 0.0) < numerical_tolerance);
// CHECK(cabs(param.yft() - 0.0) < numerical_tolerance);
// to connected
CHECK(branch.update(BranchUpdate{1, false, na_IntS}).topo);
param = branch.calc_param<symmetric_t>();
CHECK(cabs(param.yff() - 0.0) < numerical_tolerance);
CHECK(cabs(param.ytt() - ys1) < numerical_tolerance);
CHECK(cabs(param.ytf() - 0.0) < numerical_tolerance);
CHECK(cabs(param.yft() - 0.0) < numerical_tolerance);
// not connected
CHECK(branch.set_status(na_IntS, false));
param = branch.calc_param<symmetric_t>();
CHECK(cabs(param.yff() - 0.0) < numerical_tolerance);
CHECK(cabs(param.ytt() - 0.0) < numerical_tolerance);
CHECK(cabs(param.ytf() - 0.0) < numerical_tolerance);
CHECK(cabs(param.yft() - 0.0) < numerical_tolerance);
// not changing
CHECK(!branch.set_status(false, false));
// from connected
CHECK(branch.set_status(true, na_IntS));
param = branch.calc_param<symmetric_t>();
CHECK(cabs(param.yff() - ys1) < numerical_tolerance);
CHECK(cabs(param.ytt() - 0.0) < numerical_tolerance);
CHECK(cabs(param.ytf() - 0.0) < numerical_tolerance);
CHECK(cabs(param.yft() - 0.0) < numerical_tolerance);
}

SUBCASE("Asymmetric parameters") {
// double connected
// BranchCalcParam<asymmetric_t> param = asym_line.calc_param<asymmetric_t>();
// CHECK((cabs(param.yff() - yffa) < numerical_tolerance).all());
// CHECK((cabs(param.ytt() - yffa) < numerical_tolerance).all());
// CHECK((cabs(param.ytf() - yfta) < numerical_tolerance).all());
// CHECK((cabs(param.yft() - yfta) < numerical_tolerance).all());
// // no source
// param = branch.calc_param<asymmetric_t>(false);
// CHECK((cabs(param.yff() - 0.0) < numerical_tolerance).all());
// CHECK((cabs(param.ytt() - 0.0) < numerical_tolerance).all());
// CHECK((cabs(param.ytf() - 0.0) < numerical_tolerance).all());
// CHECK((cabs(param.yft() - 0.0) < numerical_tolerance).all());
// // from connected
// CHECK(branch.set_status(na_IntS, false));
// param = asym_line.calc_param<asymmetric_t>();
// CHECK((cabs(param.yff() - ysa) < numerical_tolerance).all());
// CHECK((cabs(param.ytt() - 0.0) < numerical_tolerance).all());
// CHECK((cabs(param.ytf() - 0.0) < numerical_tolerance).all());
// CHECK((cabs(param.yft() - 0.0) < numerical_tolerance).all());
BranchCalcParam<asymmetric_t> param = asym_line.calc_param<asymmetric_t>();
CHECK((cabs(param.yff() - yffa) < numerical_tolerance).all());
CHECK((cabs(param.ytt() - yffa) < numerical_tolerance).all());
CHECK((cabs(param.ytf() - yfta) < numerical_tolerance).all());
CHECK((cabs(param.yft() - yfta) < numerical_tolerance).all());
// no source
param = branch.calc_param<asymmetric_t>(false);
CHECK((cabs(param.yff() - 0.0) < numerical_tolerance).all());
CHECK((cabs(param.ytt() - 0.0) < numerical_tolerance).all());
CHECK((cabs(param.ytf() - 0.0) < numerical_tolerance).all());
CHECK((cabs(param.yft() - 0.0) < numerical_tolerance).all());
// from connected
CHECK(branch.set_status(na_IntS, false));
param = asym_line.calc_param<asymmetric_t>();
CHECK((cabs(param.yff() - ysa) < numerical_tolerance).all());
CHECK((cabs(param.ytt() - 0.0) < numerical_tolerance).all());
CHECK((cabs(param.ytf() - 0.0) < numerical_tolerance).all());
CHECK((cabs(param.yft() - 0.0) < numerical_tolerance).all());
}

// SUBCASE("Symmetric results") {
// BranchOutput<symmetric_t> output = branch.get_output<symmetric_t>(1.0, 0.9);
// CHECK(output.id == 1);
// CHECK(output.energized);
// CHECK(output.loading == doctest::Approx(loading));
// CHECK(output.i_from == doctest::Approx(cabs(i1f)));
// CHECK(output.i_to == doctest::Approx(cabs(i1t)));
// CHECK(output.s_from == doctest::Approx(cabs(s_f)));
// CHECK(output.s_to == doctest::Approx(cabs(s_t)));
// CHECK(output.p_from == doctest::Approx(real(s_f)));
// CHECK(output.p_to == doctest::Approx(real(s_t)));
// CHECK(output.q_from == doctest::Approx(imag(s_f)));
// CHECK(output.q_to == doctest::Approx(imag(s_t)));
// }

// SUBCASE("Symmetric results with direct power and current output") {
// BranchSolverOutput<symmetric_t> branch_solver_output{};
// branch_solver_output.i_f = 1.0 - 2.0i;
// branch_solver_output.i_t = 2.0 - 1.0i;
// branch_solver_output.s_f = 1.0 - 1.5i;
// branch_solver_output.s_t = 1.5 - 1.5i;
// BranchOutput<symmetric_t> output = branch.get_output<symmetric_t>(branch_solver_output);
// CHECK(output.id == 1);
// CHECK(output.energized);
// CHECK(output.loading == doctest::Approx(cabs(2.0 - 1.0i) * base_i / input.i_n));
// CHECK(output.i_from == doctest::Approx(cabs(1.0 - 2.0i) * base_i));
// CHECK(output.i_to == doctest::Approx(cabs(2.0 - 1.0i) * base_i));
// CHECK(output.s_from == doctest::Approx(cabs(1.0 - 1.5i) * base_power<symmetric_t>));
// CHECK(output.s_to == doctest::Approx(cabs(1.5 - 1.5i) * base_power<symmetric_t>));
// CHECK(output.p_from == doctest::Approx(1.0 * base_power<symmetric_t>));
// CHECK(output.p_to == doctest::Approx(1.5 * base_power<symmetric_t>));
// CHECK(output.q_from == doctest::Approx(-1.5 * base_power<symmetric_t>));
// CHECK(output.q_to == doctest::Approx(-1.5 * base_power<symmetric_t>));
// }

// SUBCASE("No source results") {
// BranchOutput<asymmetric_t> output = branch.get_null_output<asymmetric_t>();
// CHECK(output.id == 1);
// CHECK(!output.energized);
// CHECK(output.loading == 0.0);
// CHECK(output.i_from(0) == 0.0);
// CHECK(output.i_to(1) == 0.0);
// CHECK(output.s_from(2) == 0.0);
// CHECK(output.s_to(0) == 0.0);
// CHECK(output.p_from(1) == 0.0);
// CHECK(output.p_to(2) == 0.0);
// CHECK(output.q_from(0) == 0.0);
// CHECK(output.q_to(1) == 0.0);
// }

// SUBCASE("No source short circuit results") {
// BranchShortCircuitOutput output = branch.get_null_sc_output();
// CHECK(output.id == 1);
// CHECK(!output.energized);
// CHECK(output.i_from(0) == 0.0);
// CHECK(output.i_to(1) == 0.0);
// CHECK(output.i_from_angle(0) == 0.0);
// CHECK(output.i_to_angle(1) == 0.0);
// }

// SUBCASE("Asymmetric results") {
// BranchOutput<asymmetric_t> output = branch.get_output<asymmetric_t>(uaf, uat);
// CHECK(output.id == 1);
// CHECK(output.energized);
// CHECK(output.loading == doctest::Approx(loading));
// CHECK(output.i_from(0) == doctest::Approx(cabs(i1f)));
// CHECK(output.i_to(1) == doctest::Approx(cabs(i1t)));
// CHECK(output.s_from(2) == doctest::Approx(cabs(s_f) / 3.0));
// CHECK(output.s_to(0) == doctest::Approx(cabs(s_t) / 3.0));
// CHECK(output.p_from(1) == doctest::Approx(real(s_f) / 3.0));
// CHECK(output.p_to(2) == doctest::Approx(real(s_t) / 3.0));
// CHECK(output.q_from(0) == doctest::Approx(imag(s_f) / 3.0));
// CHECK(output.q_to(1) == doctest::Approx(imag(s_t) / 3.0));
// }

// SUBCASE("Asym short circuit results") {
// BranchShortCircuitOutput asym_output = branch.get_sc_output(if_sc_asym, it_sc_asym);
// CHECK(asym_output.id == 1);
// CHECK(asym_output.energized);
// CHECK(asym_output.i_from(1) == doctest::Approx(cabs(if_sc) * base_i));
// CHECK(asym_output.i_from(2) == doctest::Approx(cabs(if_sc) * base_i));
// CHECK(asym_output.i_to(0) == doctest::Approx(cabs(it_sc) * base_i));
// CHECK(asym_output.i_to(1) == doctest::Approx(cabs(it_sc) * base_i));
// CHECK(asym_output.i_from_angle(0) == doctest::Approx(pi / 4));
// CHECK(asym_output.i_from_angle(2) == doctest::Approx(pi / 4 + deg_120));
// CHECK(asym_output.i_to_angle(1) == doctest::Approx(pi / 3 - deg_120));
// CHECK(asym_output.i_to_angle(2) == doctest::Approx(pi / 3 + deg_120));
// CHECK(asym_output.id == 1);
// }

// SUBCASE("Sym short circuit results") {
// BranchShortCircuitOutput sym_output = branch.get_sc_output(if_sc, it_sc);
// BranchShortCircuitOutput asym_output = branch.get_sc_output(if_sc_asym, it_sc_asym);
// CHECK(sym_output.energized == asym_output.energized);
// CHECK(sym_output.i_from(1) == doctest::Approx(asym_output.i_from(1)));
// CHECK(sym_output.i_from(2) == doctest::Approx(asym_output.i_from(2)));
// CHECK(sym_output.i_to(0) == doctest::Approx(asym_output.i_to(0)));
// CHECK(sym_output.i_to(1) == doctest::Approx(asym_output.i_to(1)));
// CHECK(sym_output.i_from_angle(0) == doctest::Approx(asym_output.i_from_angle(0)));
// CHECK(sym_output.i_from_angle(2) == doctest::Approx(asym_output.i_from_angle(2)));
// CHECK(sym_output.i_to_angle(1) == doctest::Approx(asym_output.i_to_angle(1)));
// CHECK(sym_output.i_to_angle(2) == doctest::Approx(asym_output.i_to_angle(2)));
// }

// SUBCASE("Update inverse") {
// BranchUpdate branch_update{1, na_IntS, na_IntS};
// auto expected = branch_update;

// SUBCASE("Identical") {
// // default values
// }

// SUBCASE("From status") {
// SUBCASE("same") { branch_update.from_status = static_cast<IntS>(line.from_status()); }
// SUBCASE("different") { branch_update.from_status = IntS{0}; }
// expected.from_status = static_cast<IntS>(line.from_status());
// }

// SUBCASE("To status") {
// SUBCASE("same") { branch_update.to_status = static_cast<IntS>(line.to_status()); }
// SUBCASE("different") { branch_update.to_status = IntS{0}; }
// expected.to_status = static_cast<IntS>(line.to_status());
// }

// SUBCASE("multiple") {
// branch_update.from_status = IntS{0};
// branch_update.to_status = IntS{0};
// expected.from_status = static_cast<IntS>(line.from_status());
// expected.to_status = static_cast<IntS>(line.to_status());
// }

// auto const inv = line.inverse(branch_update);

// CHECK(inv.id == expected.id);
// CHECK(inv.from_status == expected.from_status);
// CHECK(inv.to_status == expected.to_status);
// }
}

} // namespace power_grid_model

0 comments on commit 1722b64

Please sign in to comment.