diff --git a/pandapower/auxiliary.py b/pandapower/auxiliary.py index 42c73ded1..c62e2dcb7 100644 --- a/pandapower/auxiliary.py +++ b/pandapower/auxiliary.py @@ -55,6 +55,13 @@ lightsim2grid_available = True except ImportError: lightsim2grid_available = False + +try: + from helmpy.core import helm + helmpy_available = True +except ImportError: + helmpy_available = False + try: import pandaplan.core.pplog as logging except ImportError: @@ -1341,7 +1348,7 @@ def _init_runpp_options(net, algorithm, calculate_voltage_angles, init, calculate_voltage_angles = True default_max_iteration = {"nr": 10, "iwamoto_nr": 10, "bfsw": 100, "gs": 10000, "fdxb": 30, - "fdbx": 30} + "fdbx": 30, "helm": 40} with_facts = len(net.svc.query("in_service & controllable")) > 0 or \ len(net.tcsc.query("in_service & controllable")) > 0 if max_iteration == "auto": @@ -1381,9 +1388,9 @@ def _init_runpp_options(net, algorithm, calculate_voltage_angles, init, logger.warning("Currently distributed_slack is implemented for 'ext_grid', 'gen' " "and 'xward' only, not for '" + "', '".join( false_slack_weight_elms) + "'.") - if algorithm != 'nr': + if algorithm != 'nr' and algorithm != 'helm': raise NotImplementedError( - 'Distributed slack is only implemented for Newton Raphson algorithm.') + 'Distributed slack is only implemented for Newton Raphson algorithm and HELM.') if tdpf: if algorithm != 'nr': diff --git a/pandapower/pf/run_helmpy_pf.py b/pandapower/pf/run_helmpy_pf.py new file mode 100644 index 000000000..0bae80c60 --- /dev/null +++ b/pandapower/pf/run_helmpy_pf.py @@ -0,0 +1,141 @@ +from time import perf_counter +import numpy as np +from pandapower.pypower.idx_bus import PD, QD, BUS_TYPE, PQ, REF, GS, BS, BUS_I, VM, VA +from pandapower.pypower.idx_gen import PG, QG, QMAX, QMIN, GEN_BUS, GEN_STATUS, VG, MBASE +from pandapower.pypower.bustypes import bustypes +from pandapower.pypower.makeSbus import makeSbus +from pandapower.pf.ppci_variables import _get_pf_variables_from_ppci, _store_results_from_pf_in_ppci +import pandas as pd + +try: + import pandaplan.core.pplog as logging +except ImportError: + import logging + +logger = logging.getLogger(__name__) + + +def convert_complex_to_polar_voltages(complex_voltage): + """Separate each voltage value in magnitude and phase angle (degrees)""" + polar_voltage = np.zeros((len(complex_voltage), 2), dtype=float) + polar_voltage[:, 0] = np.absolute(complex_voltage) + polar_voltage[:, 1] = np.angle(complex_voltage, deg=True) + return polar_voltage + + +def _runpf_helmpy_pf(ppci, options: dict, **kwargs): + """ + Runs a HELM based power flow. + + INPUT + ppci (dict) - the "internal" ppc (without out ot service elements and sorted elements) + options(dict) - options for the power flow + + """ + + from helmpy import helm + from helmpy.core.classes import CaseData, process_branches + + t0 = perf_counter() + # we cannot run DC pf before running newton with distributed slack because the slacks come pre-solved after the DC pf + # if isinstance(options["init_va_degree"], str) and options["init_va_degree"] == "dc": + # if options['distributed_slack']: + # pg_copy = ppci['gen'][:, PG].copy() + # pd_copy = ppci['bus'][:, PD].copy() + # ppci = _run_dc_pf(ppci, options["recycle"]) + # ppci['gen'][:, PG] = pg_copy + # ppci['bus'][:, PD] = pd_copy + # else: + # ppci = _run_dc_pf(ppci, options["recycle"]) + + max_coefficients = options['max_iteration'] + + enforce_Q_limits = False + if options["enforce_q_lims"]: + enforce_Q_limits = True + # else: + # ppci, success, iterations = _run_ac_pf_without_qlims_enforced(ppci, options) + # # update data matrices with solution store in ppci + # bus, gen, branch = ppci_to_pfsoln(ppci, options) + + DSB_model = False + if options['distributed_slack']: + DSB_model = True + +# -------------------------------------------- Calculation Case preparation -------------------------------------------- + buses = ppci['bus'].copy() + generators = ppci['gen'].copy() + branches = ppci['branch'].copy() + + generators[:, MBASE] = 100. + + N = len(buses) + N_generators = len(generators) + N_branches = len(branches) + case = CaseData(name='pandapower', N=N, N_generators=N_generators) + + case.N_branches = N_branches + case.Pd[:] = buses[:, PD] / 100 + case.Qd[:] = buses[:, QD] / 100 + + # weird but has to be done in this order. Otherwise, the values are wrong... + case.Shunt[:] = (buses[:, BS].copy()*1j + buses[:, GS]) + case.Yshunt[:] = np.copy(case.Shunt) + case.Shunt[:] = case.Shunt[:] / 100 + + for i in range(N): + case.Number_bus[buses[i][BUS_I]] = i + if buses[i][BUS_TYPE] == 3: + case.slack_bus = buses[i][BUS_I] + case.slack = i + + pos = 0 + for i in range(N_generators): + bus_i = case.Number_bus[generators[i][GEN_BUS]] + if bus_i != case.slack: + case.list_gen[pos] = bus_i + pos += 1 + case.Buses_type[bus_i] = 'PVLIM' + case.V[bus_i] = generators[i][VG] + case.Pg[bus_i] = generators[i][PG]/100 + case.Qgmax[bus_i] = generators[i][QMAX]/100 + case.Qgmin[bus_i] = generators[i][QMIN]/100 + + case.Buses_type[case.slack] = 'Slack' + case.Pg[case.slack] = 0 + + branches_df = pd.DataFrame(branches) + + process_branches(branches_df, N_branches, case) + + for i in range(N): + case.branches_buses[i].sort() # Variable that saves the branches + + case.Y[:] = np.copy(case.Ytrans) + for i in range(N): + if case.Yshunt[i].real != 0: + case.conduc_buses[i] = True + case.Y[i, i] += case.Yshunt[i] + if case.phase_barras[i]: + for k in range(len(case.phase_dict[i][0])): + case.Y[i, case.phase_dict[i][0][k]] += case.phase_dict[i][1][k] + + run, series_large, flag_divergence = helm(case, detailed_run_print=False, mismatch=1e-8, scale=1, + max_coefficients=max_coefficients, enforce_Q_limits=enforce_Q_limits, + results_file_name=None, save_results=False, pv_bus_model=1, + DSB_model=DSB_model, DSB_model_method=None,) + + et = perf_counter() - t0 + complex_voltage = run.V_complex_profile.copy() + polar_voltage = convert_complex_to_polar_voltages(complex_voltage) # Calculate polar voltage + + buses[:, VM] = polar_voltage[:, 0] + buses[:, VA] = polar_voltage[:, 1] + + ppci["bus"], ppci["gen"], ppci["branch"] = buses, generators, branches + ppci["success"] = not flag_divergence + ppci["internal"]["V"] = complex_voltage + ppci["iterations"] = series_large + ppci["et"] = et + + return ppci diff --git a/pandapower/powerflow.py b/pandapower/powerflow.py index dffde16e5..4231f7d32 100644 --- a/pandapower/powerflow.py +++ b/pandapower/powerflow.py @@ -14,6 +14,7 @@ from pandapower.pf.run_dc_pf import _run_dc_pf from pandapower.pf.run_newton_raphson_pf import _run_newton_raphson_pf from pandapower.pf.runpf_pypower import _runpf_pypower +from pandapower.pf.run_helmpy_pf import _runpf_helmpy_pf from pandapower.pypower.bustypes import bustypes from pandapower.pypower.idx_bus import VM from pandapower.pypower.makeYbus import makeYbus as makeYbus_pypower @@ -166,6 +167,8 @@ def _run_pf_algorithm(ppci, options, **kwargs): result = _run_newton_raphson_pf(ppci, options) elif algorithm in ['fdbx', 'fdxb', 'gs']: # algorithms existing within pypower result = _runpf_pypower(ppci, options, **kwargs)[0] + elif algorithm == 'helm': + result = _runpf_helmpy_pf(ppci, options, **kwargs) else: raise AlgorithmUnknown("Algorithm {0} is unknown!".format(algorithm)) else: diff --git a/pandapower/run.py b/pandapower/run.py index d5241e9b8..f50aa8441 100644 --- a/pandapower/run.py +++ b/pandapower/run.py @@ -90,6 +90,7 @@ def runpp(net, algorithm='nr', calculate_voltage_angles="auto", init="auto", - "gs" gauss-seidel (pypower implementation) - "fdbx" fast-decoupled (pypower implementation) - "fdxb" fast-decoupled (pypower implementation) + - "helm" holomorphic embedded loadflow method **calculate_voltage_angles** (str or bool, "auto") - consider voltage angles in loadflow calculation diff --git a/pandapower/test/loadflow/test_runpp.py b/pandapower/test/loadflow/test_runpp.py index 8a126fe3b..1a586fbb2 100644 --- a/pandapower/test/loadflow/test_runpp.py +++ b/pandapower/test/loadflow/test_runpp.py @@ -13,7 +13,7 @@ import pandapower as pp import pandapower.control -from pandapower.auxiliary import _check_connectivity, _add_ppc_options, lightsim2grid_available +from pandapower.auxiliary import _check_connectivity, _add_ppc_options, lightsim2grid_available, helmpy_available from pandapower.networks import create_cigre_network_mv, four_loads_with_branches_out, \ example_simple, simple_four_bus_system, example_multivoltage from pandapower.pd2ppc import _pd2ppc @@ -510,6 +510,40 @@ def test_bsfw_algorithm(): assert np.allclose(va_nr, va_alg) +@pytest.mark.skipif(not helmpy_available, reason="HELMpy is not installed") +def test_helm_algorithm_simple(): + import pandapower.networks as nw + net = nw.case9() + + pp.runpp(net) + vm_nr = copy.copy(net.res_bus.vm_pu) + va_nr = copy.copy(net.res_bus.va_degree) + + pp.runpp(net, algorithm='helm') + vm_alg = net.res_bus.vm_pu + va_alg = net.res_bus.va_degree + + assert np.allclose(vm_nr, vm_alg) + assert np.allclose(va_nr, va_alg) + + +@pytest.mark.skipif(not helmpy_available, reason="HELMpy is not installed") +def test_helm_algorithm_complex(): + import pandapower.networks as nw + net = nw.case9() + + pp.runpp(net) + vm_nr = copy.copy(net.res_bus.vm_pu) + va_nr = copy.copy(net.res_bus.va_degree) + + pp.runpp(net, algorithm='helm') + vm_alg = net.res_bus.vm_pu + va_alg = net.res_bus.va_degree + + assert np.allclose(vm_nr, vm_alg) + assert np.allclose(va_nr, va_alg) + + @pytest.mark.xfail(reason="unknown") def test_bsfw_algorithm_multi_net(): net = example_simple()