diff --git a/examples/small_signal_and_emt/2-bus_2-src_swing/inputs/buses.csv b/examples/small_signal_and_emt/2-bus_2-src_swing/inputs/buses.csv new file mode 100644 index 0000000..b50acd6 --- /dev/null +++ b/examples/small_signal_and_emt/2-bus_2-src_swing/inputs/buses.csv @@ -0,0 +1,3 @@ +name,base_power_MVA,base_voltage_kV,base_frequency_Hz,minimum_voltage_pu,maximum_voltage_pu +lima,1.00E+02,2.30E+02,60,1,1 +santiago,1.00E+02,2.30E+02,60,0.95,1.3 \ No newline at end of file diff --git a/examples/small_signal_and_emt/2-bus_2-src_swing/inputs/lines.csv b/examples/small_signal_and_emt/2-bus_2-src_swing/inputs/lines.csv new file mode 100644 index 0000000..74a86b2 --- /dev/null +++ b/examples/small_signal_and_emt/2-bus_2-src_swing/inputs/lines.csv @@ -0,0 +1,2 @@ +name,from_bus,to_bus,base_power_MVA,base_voltage_kV,base_frequency_Hz,r_pu,x_pu,g_pu,b_pu +tx_1,lima,santiago,1.00E+02,2.30E+02,60,0.01,0.5,0.05,0.066666667 \ No newline at end of file diff --git a/examples/small_signal_and_emt/2-bus_2-src_swing/inputs/loads.csv b/examples/small_signal_and_emt/2-bus_2-src_swing/inputs/loads.csv new file mode 100644 index 0000000..ed9aa45 --- /dev/null +++ b/examples/small_signal_and_emt/2-bus_2-src_swing/inputs/loads.csv @@ -0,0 +1,3 @@ +bus,timepoint,load_MW,load_MVAR +lima,t_1,0,0 +santiago,t_1,0,0 \ No newline at end of file diff --git a/examples/small_signal_and_emt/2-bus_2-src_swing/inputs/sources_with_swing.csv b/examples/small_signal_and_emt/2-bus_2-src_swing/inputs/sources_with_swing.csv new file mode 100644 index 0000000..2269f3c --- /dev/null +++ b/examples/small_signal_and_emt/2-bus_2-src_swing/inputs/sources_with_swing.csv @@ -0,0 +1,3 @@ +name,bus,minimum_active_power_MW,maximum_active_power_MW,minimum_reactive_power_MVAR,maximum_reactive_power_MVAR,base_power_MVA,base_voltage_kV,base_frequency_Hz,cost_variable_USDperMWh,r_pu,x_pu,inertia_constant_s,damping_pu +gen1,lima,-200,200,-500,500,1.00E+02,2.30E+02,60,0,0.01,0.5,3,1 +gen2,santiago,100,100,74,75,1.00E+02,2.30E+02,60,0,0.05,0.2,3,1 \ No newline at end of file diff --git a/examples/small_signal_and_emt/2-bus_2-src_swing/inputs/timepoints.csv b/examples/small_signal_and_emt/2-bus_2-src_swing/inputs/timepoints.csv new file mode 100644 index 0000000..fb3202e --- /dev/null +++ b/examples/small_signal_and_emt/2-bus_2-src_swing/inputs/timepoints.csv @@ -0,0 +1,2 @@ +name +t_1 \ No newline at end of file diff --git a/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/ac_power_flow/active_power_balance_by_bus.csv b/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/ac_power_flow/active_power_balance_by_bus.csv new file mode 100644 index 0000000..89c366a --- /dev/null +++ b/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/ac_power_flow/active_power_balance_by_bus.csv @@ -0,0 +1,3 @@ +bus,timepoint,generator_dispatch_MW,load_shedding_MW,load_MW,net_line_leaving_flow_MW +lima,t_1,-85.86698271936963,-9.97494096261266e-9,0.0,-85.86698272918828 +santiago,t_1,100.0,-9.974940980027277e-9,0.0,99.99999998988655 diff --git a/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/ac_power_flow/bus_voltage.csv b/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/ac_power_flow/bus_voltage.csv new file mode 100644 index 0000000..0c979a0 --- /dev/null +++ b/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/ac_power_flow/bus_voltage.csv @@ -0,0 +1,3 @@ +id,bus,timepoint,voltage_magnitude_pu,voltage_angle_deg +0,lima,t_1,1.0,0.0 +1,santiago,t_1,1.2780887210409115,20.643650769134787 diff --git a/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/ac_power_flow/costs_summary.csv b/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/ac_power_flow/costs_summary.csv new file mode 100644 index 0000000..2d7a436 --- /dev/null +++ b/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/ac_power_flow/costs_summary.csv @@ -0,0 +1,2 @@ +component,cost +total_cost_USD,-0.00003989976382117789 diff --git a/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/ac_power_flow/generator_dispatch.csv b/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/ac_power_flow/generator_dispatch.csv new file mode 100644 index 0000000..24fe03e --- /dev/null +++ b/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/ac_power_flow/generator_dispatch.csv @@ -0,0 +1,3 @@ +id,type,generator,timepoint,active_power_MW,reactive_power_MVAR +0,sources_with_swing,gen1,t_1,-85.86698271936963,-44.05416612599554 +1,sources_with_swing,gen2,t_1,100.0,74.77059585403816 diff --git a/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/ac_power_flow/line_flows.csv b/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/ac_power_flow/line_flows.csv new file mode 100644 index 0000000..5459b8a --- /dev/null +++ b/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/ac_power_flow/line_flows.csv @@ -0,0 +1,2 @@ +line,from_bus,to_bus,existing_capacity_MW,active_power_from_bus_MW,reactive_power_from_bus_MVAR,active_power_to_bus_MW,reactive_power_to_bus_MVAR,active_power_loss_MW,reactive_power_loss_MVAR +tx_1,lima,santiago,inf,-85.86698272918831,-44.0541661359936,99.99999998988658,74.77059584444507,14.133017260698267,30.71642970845147 diff --git a/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/ac_power_flow/load_shedding.csv b/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/ac_power_flow/load_shedding.csv new file mode 100644 index 0000000..19fbe0a --- /dev/null +++ b/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/ac_power_flow/load_shedding.csv @@ -0,0 +1,3 @@ +bus,timepoint,active_load_shedding_MW,reactive_load_shedding_MVAR +lima,t_1,-9.97494096261266e-9,-9.974940964279889e-9 +santiago,t_1,-9.974940980027277e-9,-9.97494091425806e-9 diff --git a/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/ac_power_flow/reactive_power_balance_by_bus.csv b/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/ac_power_flow/reactive_power_balance_by_bus.csv new file mode 100644 index 0000000..af27197 --- /dev/null +++ b/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/ac_power_flow/reactive_power_balance_by_bus.csv @@ -0,0 +1,3 @@ +bus,timepoint,generator_dispatch_MVAR,load_shedding_MVAR,load_MVAR,net_line_leaving_flow_MVAR +lima,t_1,-44.05416612599554,-9.974940964279889e-9,0.0,-44.054166135993555 +santiago,t_1,74.77059585403816,-9.97494091425806e-9,0.0,74.77059584444504 diff --git a/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/ac_power_flow/solver_status.csv b/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/ac_power_flow/solver_status.csv new file mode 100644 index 0000000..96d02af --- /dev/null +++ b/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/ac_power_flow/solver_status.csv @@ -0,0 +1,5 @@ +attribute,value +solver_name,ipopt +solver_status,ok +termination_condition,optimal +time_spent_seconds,0.021149158477783203 diff --git a/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/component_connection_matrices/F.csv b/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/component_connection_matrices/F.csv new file mode 100644 index 0000000..aca6ead --- /dev/null +++ b/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/component_connection_matrices/F.csv @@ -0,0 +1,19 @@ +Index,"('sources_with_swing_0', 'i_bus_D')","('sources_with_swing_0', 'i_bus_Q')","('sources_with_swing_1', 'i_bus_D')","('sources_with_swing_1', 'i_bus_Q')","('shunt_parallel_rc_0', 'v_bus_D')","('shunt_parallel_rc_0', 'v_bus_Q')","('shunt_parallel_rc_1', 'v_bus_D')","('shunt_parallel_rc_1', 'v_bus_Q')","('branch_series_rl_0', 'i_br_D')","('branch_series_rl_0', 'i_br_Q')" +"('sources_with_swing_0', 'p_m')",0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0 +"('sources_with_swing_0', 'v_ref_d')",0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0 +"('sources_with_swing_0', 'v_ref_q')",0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0 +"('sources_with_swing_0', 'v_bus_D')",0.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0,0.0,0.0 +"('sources_with_swing_0', 'v_bus_Q')",0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0,0.0 +"('sources_with_swing_1', 'p_m')",0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0 +"('sources_with_swing_1', 'v_ref_d')",0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0 +"('sources_with_swing_1', 'v_ref_q')",0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0 +"('sources_with_swing_1', 'v_bus_D')",0.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0 +"('sources_with_swing_1', 'v_bus_Q')",0.0,0.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0 +"('shunt_parallel_rc_0', 'i_bus_D')",1.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,-1.0,0.0 +"('shunt_parallel_rc_0', 'i_bus_Q')",0.0,1.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,-1.0 +"('shunt_parallel_rc_1', 'i_bus_D')",0.0,0.0,1.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0 +"('shunt_parallel_rc_1', 'i_bus_Q')",0.0,0.0,0.0,1.0,0.0,0.0,0.0,0.0,0.0,1.0 +"('branch_series_rl_0', 'v_from_bus_D')",0.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0,0.0,0.0 +"('branch_series_rl_0', 'v_from_bus_Q')",0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0,0.0 +"('branch_series_rl_0', 'v_to_bus_D')",0.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0 +"('branch_series_rl_0', 'v_to_bus_Q')",0.0,0.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0 diff --git a/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/component_connection_matrices/G.csv b/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/component_connection_matrices/G.csv new file mode 100644 index 0000000..4656042 --- /dev/null +++ b/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/component_connection_matrices/G.csv @@ -0,0 +1,19 @@ +Index,"('sources_with_swing_0', 'p_m')","('sources_with_swing_0', 'v_ref_d')","('sources_with_swing_0', 'v_ref_q')","('sources_with_swing_1', 'p_m')","('sources_with_swing_1', 'v_ref_d')","('sources_with_swing_1', 'v_ref_q')" +"('sources_with_swing_0', 'p_m')",1.0,0.0,0.0,0.0,0.0,0.0 +"('sources_with_swing_0', 'v_ref_d')",0.0,1.0,0.0,0.0,0.0,0.0 +"('sources_with_swing_0', 'v_ref_q')",0.0,0.0,1.0,0.0,0.0,0.0 +"('sources_with_swing_0', 'v_bus_D')",0.0,0.0,0.0,0.0,0.0,0.0 +"('sources_with_swing_0', 'v_bus_Q')",0.0,0.0,0.0,0.0,0.0,0.0 +"('sources_with_swing_1', 'p_m')",0.0,0.0,0.0,1.0,0.0,0.0 +"('sources_with_swing_1', 'v_ref_d')",0.0,0.0,0.0,0.0,1.0,0.0 +"('sources_with_swing_1', 'v_ref_q')",0.0,0.0,0.0,0.0,0.0,1.0 +"('sources_with_swing_1', 'v_bus_D')",0.0,0.0,0.0,0.0,0.0,0.0 +"('sources_with_swing_1', 'v_bus_Q')",0.0,0.0,0.0,0.0,0.0,0.0 +"('shunt_parallel_rc_0', 'i_bus_D')",0.0,0.0,0.0,0.0,0.0,0.0 +"('shunt_parallel_rc_0', 'i_bus_Q')",0.0,0.0,0.0,0.0,0.0,0.0 +"('shunt_parallel_rc_1', 'i_bus_D')",0.0,0.0,0.0,0.0,0.0,0.0 +"('shunt_parallel_rc_1', 'i_bus_Q')",0.0,0.0,0.0,0.0,0.0,0.0 +"('branch_series_rl_0', 'v_from_bus_D')",0.0,0.0,0.0,0.0,0.0,0.0 +"('branch_series_rl_0', 'v_from_bus_Q')",0.0,0.0,0.0,0.0,0.0,0.0 +"('branch_series_rl_0', 'v_to_bus_D')",0.0,0.0,0.0,0.0,0.0,0.0 +"('branch_series_rl_0', 'v_to_bus_Q')",0.0,0.0,0.0,0.0,0.0,0.0 diff --git a/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/component_connection_matrices/H.csv b/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/component_connection_matrices/H.csv new file mode 100644 index 0000000..2d6a91f --- /dev/null +++ b/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/component_connection_matrices/H.csv @@ -0,0 +1,11 @@ +Index,"('sources_with_swing_0', 'i_bus_D')","('sources_with_swing_0', 'i_bus_Q')","('sources_with_swing_1', 'i_bus_D')","('sources_with_swing_1', 'i_bus_Q')","('shunt_parallel_rc_0', 'v_bus_D')","('shunt_parallel_rc_0', 'v_bus_Q')","('shunt_parallel_rc_1', 'v_bus_D')","('shunt_parallel_rc_1', 'v_bus_Q')","('branch_series_rl_0', 'i_br_D')","('branch_series_rl_0', 'i_br_Q')" +"('sources_with_swing_0', 'i_bus_D')",1.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0 +"('sources_with_swing_0', 'i_bus_Q')",0.0,1.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0 +"('sources_with_swing_1', 'i_bus_D')",0.0,0.0,1.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0 +"('sources_with_swing_1', 'i_bus_Q')",0.0,0.0,0.0,1.0,0.0,0.0,0.0,0.0,0.0,0.0 +"('shunt_parallel_rc_0', 'v_bus_D')",0.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0,0.0,0.0 +"('shunt_parallel_rc_0', 'v_bus_Q')",0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0,0.0 +"('shunt_parallel_rc_1', 'v_bus_D')",0.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0,0.0 +"('shunt_parallel_rc_1', 'v_bus_Q')",0.0,0.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0,0.0 +"('branch_series_rl_0', 'i_br_D')",0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,1.0,0.0 +"('branch_series_rl_0', 'i_br_Q')",0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,1.0 diff --git a/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/component_connection_matrices/L.csv b/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/component_connection_matrices/L.csv new file mode 100644 index 0000000..94ebce2 --- /dev/null +++ b/examples/small_signal_and_emt/2-bus_2-src_swing/outputs/component_connection_matrices/L.csv @@ -0,0 +1,11 @@ +Index,"('sources_with_swing_0', 'p_m')","('sources_with_swing_0', 'v_ref_d')","('sources_with_swing_0', 'v_ref_q')","('sources_with_swing_1', 'p_m')","('sources_with_swing_1', 'v_ref_d')","('sources_with_swing_1', 'v_ref_q')" +"('sources_with_swing_0', 'i_bus_D')",0.0,0.0,0.0,0.0,0.0,0.0 +"('sources_with_swing_0', 'i_bus_Q')",0.0,0.0,0.0,0.0,0.0,0.0 +"('sources_with_swing_1', 'i_bus_D')",0.0,0.0,0.0,0.0,0.0,0.0 +"('sources_with_swing_1', 'i_bus_Q')",0.0,0.0,0.0,0.0,0.0,0.0 +"('shunt_parallel_rc_0', 'v_bus_D')",0.0,0.0,0.0,0.0,0.0,0.0 +"('shunt_parallel_rc_0', 'v_bus_Q')",0.0,0.0,0.0,0.0,0.0,0.0 +"('shunt_parallel_rc_1', 'v_bus_D')",0.0,0.0,0.0,0.0,0.0,0.0 +"('shunt_parallel_rc_1', 'v_bus_Q')",0.0,0.0,0.0,0.0,0.0,0.0 +"('branch_series_rl_0', 'i_br_D')",0.0,0.0,0.0,0.0,0.0,0.0 +"('branch_series_rl_0', 'i_br_Q')",0.0,0.0,0.0,0.0,0.0,0.0 diff --git a/examples/small_signal_and_emt/2-bus_2-src_swing/run.py b/examples/small_signal_and_emt/2-bus_2-src_swing/run.py new file mode 100644 index 0000000..9175f6c --- /dev/null +++ b/examples/small_signal_and_emt/2-bus_2-src_swing/run.py @@ -0,0 +1,75 @@ +""" +Simulates two sources with swing dynamics connected via a transmission line. + +First, we compute the system-wide small-signal model using STING. +This small-signal model also contains EMT initial conditions. + +You should obtain the following eigenvalues: +shape: (14, 5) +┌──────────┬───────────┬──────────────────────┬──────────────────┬───────────────────────┐ +│ real ┆ imag ┆ natural_frequency_hz ┆ damping_ratio_pu ┆ time_constant_seconds │ +│ --- ┆ --- ┆ --- ┆ --- ┆ --- │ +│ f64 ┆ f64 ┆ f64 ┆ f64 ┆ f64 │ +╞══════════╪═══════════╪══════════════════════╪══════════════════╪═══════════════════════╡ +│ 0.0 ┆ 0.0 ┆ 0.0 ┆ -1.0 ┆ -2.6276e7 │ +│ -0.063 ┆ 8.892 ┆ 1.415 ┆ 0.007 ┆ 15.8731 │ +│ -0.063 ┆ -8.892 ┆ 1.415 ┆ 0.007 ┆ 15.8731 │ +│ -0.163 ┆ 0.0 ┆ 0.026 ┆ 1.0 ┆ 6.1304 │ +│ -21.987 ┆ 376.796 ┆ 60.071 ┆ 0.058 ┆ 0.0455 │ +│ -21.987 ┆ -376.796 ┆ 60.071 ┆ 0.058 ┆ 0.0455 │ +│ -159.587 ┆ 2149.013 ┆ 342.968 ┆ 0.074 ┆ 0.0063 │ +│ -159.587 ┆ -2149.013 ┆ 342.968 ┆ 0.074 ┆ 0.0063 │ +│ -159.589 ┆ 2902.987 ┆ 462.722 ┆ 0.055 ┆ 0.0063 │ +│ -159.589 ┆ -2902.987 ┆ 462.722 ┆ 0.055 ┆ 0.0063 │ +│ -166.836 ┆ 3750.796 ┆ 597.548 ┆ 0.044 ┆ 0.006 │ +│ -166.836 ┆ -3750.796 ┆ 597.548 ┆ 0.044 ┆ 0.006 │ +│ -166.837 ┆ 4504.775 ┆ 717.449 ┆ 0.037 ┆ 0.006 │ +│ -166.837 ┆ -4504.775 ┆ 717.449 ┆ 0.037 ┆ 0.006 │ +└──────────┴───────────┴──────────────────────┴──────────────────┴───────────────────────┘ +""" + +# Import Python standard and third-party packages +from pathlib import Path +# Import sting package +from sting import main +from sting.system.core import System +import os +import numpy as np +import polars as pl + +# Specify path of the case study directory +case_dir = Path(__file__).resolve().parent + +# Step function inputs to simulate +def step1(t): + return 0.01 if t >= 0.5 else 0.0 + +def step2(t): + return 0.0 + +inputs = { + 'sources_with_swing_0': { + 'v_ref_d': step1 + }, + 'sources_with_swing_1': { + 'v_ref_d': step2 + } + } +t_max = 1.0 # Simulation length in seconds + +# Construct system and small-signal model +sys, ssm = main.run_ssm(case_directory=case_dir) +ssm.simulate_ssm(t_max=t_max, inputs=inputs) + +# Run EMT simulation +main.run_emt(case_directory=case_dir, inputs=inputs, t_max=t_max) + +emt_dir = os.path.join(case_dir, "outputs", "simulation_emt") +ssm_dir = os.path.join(case_dir, "outputs", "small_signal_model") + +ans = dict() +for component in sys: + if hasattr(component, "compare_ssm_emt"): + ans |= getattr(component, "compare_ssm_emt")(emt_dir, ssm_dir) + +mprint("\nok") \ No newline at end of file diff --git a/sting/generator/core.py b/sting/generator/core.py index 07a7945..7eb4e2e 100644 --- a/sting/generator/core.py +++ b/sting/generator/core.py @@ -46,6 +46,8 @@ class Generator(Component): base_power_MVA: float = None base_voltage_kV: float = None base_frequency_Hz: float = None + damping_pu: float = None + inertia_constant_s: float = None technology: str = None site: str = None cap_existing_power_MW: float = None diff --git a/sting/generator/source_with_swing.py b/sting/generator/source_with_swing.py new file mode 100644 index 0000000..be13333 --- /dev/null +++ b/sting/generator/source_with_swing.py @@ -0,0 +1,289 @@ +""" +This module implements a voltage source that incorporates: +- A variable frequency and swing dynamics. +- Equilibrium assumptions: + - The source frequency equals to the nominal bus frequency. + - Mechanical input power equals to the electrical output power of the source. +- Series RL branch: It is in series with the voltage source. +""" +# ------------- +# Import python packages +# -------------- +import numpy as np +from scipy.linalg import block_diag +from dataclasses import dataclass +from typing import NamedTuple +import plotly.graph_objects as go +from plotly.subplots import make_subplots +import os +import polars as pl +# ------------- +# Import sting code +# ------------- +from sting.utils.dynamical_systems import StateSpaceModel, DynamicalVariables +from sting.utils.transformations import dq02abc, abc2dq0 +from sting.generator.core import Generator + +# ------------- +# Sub-classes +# ------------- +class InitialConditionsEMT(NamedTuple): + v_bus_D: float + v_bus_Q: float + v_int_d: float + v_int_q: float + i_bus_d: float + i_bus_q: float + i_bus_D: float + i_bus_Q: float + angle_ref: float + p_source: float + +class VariablesEMT(NamedTuple): + x: DynamicalVariables + u: DynamicalVariables + y: DynamicalVariables + +# ------------- +# Main class +# ------------- +@dataclass(slots=True, kw_only=True, eq=False) +class SourceWithSwing(Generator): + r_pu: float + x_pu: float + inertia_constant_s: float + damping_pu: float + + emt_init: InitialConditionsEMT = None + + def _build_small_signal_model(self): + + r = self.r_pu + x = self.x_pu + h = self.inertia_constant_s + d = self.damping_pu + + wb = 2 * np.pi * self.base_frequency_Hz + angle_ref = self.emt_init.angle_ref * np.pi / 180 + cosphi = np.cos(self.emt_init.angle_ref * np.pi / 180) + sinphi = np.sin(self.emt_init.angle_ref * np.pi / 180) + + v_bus_D, v_bus_Q = self.emt_init.v_bus_D, self.emt_init.v_bus_Q + v_int_d, v_int_q = self.emt_init.v_int_d, self.emt_init.v_int_q + + i_bus_D, i_bus_Q = self.emt_init.i_bus_D, self.emt_init.i_bus_Q + i_bus_d, i_bus_q = self.emt_init.i_bus_d, self.emt_init.i_bus_q + + p_source = self.emt_init.p_source + + # Roation matrix (turn off code formatters for matrices) + # fmt: off + R = np.array( + [[cosphi, -sinphi], + [sinphi, cosphi]]) + + # Define state-space matrices + A = wb * np.array( + [[0, 1/wb, 0, 0 ], + [0, -d / (2.0 * h * wb), -v_int_d/(2*h), -v_int_q/(2*h) ], + [1/x * (sinphi * v_bus_D - cosphi * v_bus_Q), i_bus_q/wb, -r/x, 1 ], + [1/x * (cosphi * v_bus_D + sinphi * v_bus_Q), -i_bus_d/wb, -1, -r/x ]]) + + B = wb * np.array( + [[0, 0, 0, 0, 0 ], + [1/(2*h), -i_bus_d/(2*h), -i_bus_q/(2*h), 0, 0 ], + [0, 1/x, 0, -cosphi/x, -sinphi/x ], + [0, 0, 1/x, sinphi/x, -cosphi/x ]]) + # B = B @ block_diag(np.eye(2), R.T) + # fmt: on + C = np.hstack((np.array([[-i_bus_Q, 0], [i_bus_D, 0]]), R)) + + D = np.zeros((2, 5)) + + # Inputs + u = DynamicalVariables( + name=["p_m", "v_ref_d", "v_ref_q", "v_bus_D", "v_bus_Q"], + component=f"{self.type_}_{self.id}", + type=["device", "device", "device", "grid", "grid"], + init=[0.0, 0.0, 0.0, 0.0, 0.0], + ) + + # Outputs + y = DynamicalVariables( + name=["i_bus_D", "i_bus_Q"], + component=f"{self.type_}_{self.id}", + init=[0.0, 0.0], + ) + + # States + x = DynamicalVariables( + name=["delta", "omega", "i_bus_d", "i_bus_q"], + component=f"{self.type_}_{self.id}", + init=[0.0, 0.0, 0.0, 0.0], + ) + + self.ssm = StateSpaceModel(A=A, B=B, C=C, D=D, u=u, y=y, x=x) + + def _calculate_emt_initial_conditions(self): + vmag_bus = self.power_flow_variables.vmag_bus + vphase_bus = self.power_flow_variables.vphase_bus + p_bus = self.power_flow_variables.p_bus + q_bus = self.power_flow_variables.q_bus + + v_bus_DQ = vmag_bus * np.exp(vphase_bus * 1j * np.pi / 180) + i_bus_DQ = ((p_bus + 1j * q_bus) / v_bus_DQ).conjugate() + + v_int_DQ = v_bus_DQ + i_bus_DQ * (self.r_pu + 1j * self.x_pu) + angle_ref = np.angle(v_int_DQ, deg=True) + + v_int_dq = v_int_DQ * np.exp(-angle_ref * np.pi / 180 * 1j) + i_bus_dq = i_bus_DQ * np.exp(-angle_ref * np.pi / 180 * 1j) + + p_source = v_int_dq * i_bus_dq.conjugate() + + self.emt_init = InitialConditionsEMT( + v_bus_D=v_bus_DQ.real, + v_bus_Q=v_bus_DQ.imag, + v_int_d=v_int_dq.real, + v_int_q=v_int_dq.imag, + i_bus_d=i_bus_dq.real, + i_bus_q=i_bus_dq.imag, + i_bus_D=i_bus_DQ.real, + i_bus_Q=i_bus_DQ.imag, + angle_ref=angle_ref, + p_source=p_source.real, + ) + + def define_variables_emt(self): + + # States + # ------ + + # Initial conditions + i_bus_d, i_bus_q = self.emt_init.i_bus_d, self.emt_init.i_bus_q + angle_ref = self.emt_init.angle_ref * np.pi / 180 + i_bus_a, i_bus_b, i_bus_c = dq02abc(i_bus_d, i_bus_q, 0, angle_ref) + + wb = 2 * np.pi * self.base_frequency_Hz + + x = DynamicalVariables( + name=["i_bus_a", "i_bus_b", "i_bus_c", "angle_ref", "omega"], + component=f"{self.type_}_{self.id}", + init=[i_bus_a, i_bus_b, i_bus_c, angle_ref, wb], + ) + + # Inputs + # ------ + + # Initial conditions + v_ref_d, v_ref_q = self.emt_init.v_int_d, self.emt_init.v_int_q + v_bus_D, v_bus_Q = self.emt_init.v_bus_D, self.emt_init.v_bus_Q + v_bus_a, v_bus_b, v_bus_c = dq02abc(v_bus_D, v_bus_Q, 0, 0) + p_m = self.emt_init.p_source + + u = DynamicalVariables( + name=["p_m", "v_ref_d", "v_ref_q", "v_bus_a", "v_bus_b", "v_bus_c"], + component=f"{self.type_}_{self.id}", + type=["device", "device", "device", "grid", "grid", "grid"], + init=[p_m, v_ref_d, v_ref_q, v_bus_a, v_bus_b, v_bus_c], + ) + + # Outputs + y = DynamicalVariables( + name=["i_bus_a", "i_bus_b", "i_bus_c"], + component=f"{self.type_}_{self.id}", + ) + + self.variables_emt = VariablesEMT(x=x, u=u, y=y) + + def get_derivative_state_emt(self): + + # Get state values + i_bus_a, i_bus_b, i_bus_c, angle_ref, omega = self.variables_emt.x.value + + # Get input values + p_m, v_ref_d, v_ref_q, v_bus_a, v_bus_b, v_bus_c = self.variables_emt.u.value + + v_ref_a, v_ref_b, v_ref_c = dq02abc(v_ref_d, v_ref_q, 0, angle_ref) + i_bus_d, i_bus_q, _ = abc2dq0(i_bus_a, i_bus_b, i_bus_c, angle_ref) + + # Get parameters + r = self.r_pu + x = self.x_pu + h = self.inertia_constant_s + d = self.damping_pu + + wb = 2 * np.pi * self.base_frequency_Hz + + p_e = v_ref_d * i_bus_d + v_ref_q * i_bus_q + + # Differential equations + d_i_bus_a = wb / x * (v_ref_a - v_bus_a - r * i_bus_a) + d_i_bus_b = wb / x * (v_ref_b - v_bus_b - r * i_bus_b) + d_i_bus_c = wb / x * (v_ref_c - v_bus_c - r * i_bus_c) + d_angle_ref = omega + d_omega = wb / (2.0 * h) * (p_m - p_e - d * (omega - wb) / wb) + + return [d_i_bus_a, d_i_bus_b, d_i_bus_c, d_angle_ref, d_omega] + + def get_output_emt(self): + + i_bus_a, i_bus_b, i_bus_c, angle_ref, omega = self.variables_emt.x.value + + return [i_bus_a, i_bus_b, i_bus_c] + + def plot_results_emt(self): + """ + Plot EMT simulation results + """ + + i_bus_a, i_bus_b, i_bus_c, angle_ref, omega = self.variables_emt.x.value + i_bus_d, i_bus_q, _ = zip(*map(abc2dq0, i_bus_a, i_bus_b, i_bus_c, angle_ref)) + t = self.variables_emt.x.time + + wb = 2.0 * np.pi * self.base_frequency_Hz + angle_ref_init = self.emt_init.angle_ref * np.pi / 180.0 + + delta_dev = angle_ref - angle_ref_init - wb * t + omega_dev = omega - wb + i_bus_d_dev = i_bus_d - self.emt_init.i_bus_d + i_bus_q_dev = i_bus_q - self.emt_init.i_bus_q + + results = DynamicalVariables( + name=["delta", "omega", "i_bus_d", "i_bus_q"], + component=f"{self.type_}_{self.id}", + value=[delta_dev, omega_dev, i_bus_d_dev, i_bus_q_dev], + time=t, + ) + return results + + + def compare_ssm_emt(self, emt_directory, ssm_directory): + # Read the SSM and EMT states + emt = pl.read_csv(os.path.join(emt_directory, f"{self.type_}_{self.id}_states.csv")) + ssm = pl.read_csv(os.path.join(ssm_directory, f"{self.type_}_{self.id}_states.csv")) + + # Transform EMT abc states to dq0 states + + i_a, i_b, i_c, angle_ref, omega = [c.to_numpy() for c in emt.select("i_bus_a", "i_bus_b", "i_bus_c", "angle_ref", "omega")] + i_emt_d, i_emt_q, _ = zip(*map(abc2dq0, i_a, i_b, i_c, angle_ref)) + t = emt["time"].to_numpy() + + wb = 2.0 * np.pi * self.base_frequency_Hz + angle_ref_init = self.emt_init.angle_ref * np.pi / 180.0 + + delta_emt_dev = (angle_ref - angle_ref_init - wb * t) + omega_emt_dev = omega - wb + i_emt_d_dev = i_emt_d - self.emt_init.i_bus_d + i_emt_q_dev = i_emt_q - self.emt_init.i_bus_q + + # Unpack the SSM dq states + delta_ssm, omega_ssm, i_ssm_d, i_ssm_q = [c.to_numpy() for c in ssm.select("delta", "omega", "i_bus_d", "i_bus_q")] + + # Return deltas + return { + f"({self.type_}_{self.id}, delta)": (delta_emt_dev, delta_ssm), + f"({self.type_}_{self.id}, omega)": (omega_emt_dev, omega_ssm), + f"({self.type_}_{self.id}, i_bus_d)": (i_emt_d_dev, i_ssm_d), + f"({self.type_}_{self.id}, i_bus_q)": (i_emt_q_dev, i_ssm_q) + } \ No newline at end of file diff --git a/sting/system/core.py b/sting/system/core.py index 071724e..5f8cfd5 100644 --- a/sting/system/core.py +++ b/sting/system/core.py @@ -25,6 +25,7 @@ from sting.storage.core import Storage from sting.generator.infinite_source import InfiniteSource from sting.generator.switching_load import SwitchingLoad +from sting.generator.source_with_swing import SourceWithSwing from sting.generator.gfli_a import GFLIa from sting.generator.gfmi_c import GFMIc from sting.generator.gfmi_d import GFMId @@ -57,6 +58,7 @@ class System: storage: list[Storage] = None infinite_sources: list[InfiniteSource] = None switching_loads: list[SwitchingLoad] = None + sources_with_swing: list[SourceWithSwing] = None gfmi_c: list[GFMIc] = None gfmi_d: list[GFMId] = None gfmi_e: list[GFMIe] = None