Improve DC initialisation

src/main/java/com/powsybl/openloadflow/OpenLoadFlowParameters.java

@@ -1676,7 +1676,7 @@ public static DcLoadFlowParameters createDcParameters(LoadFlowParameters paramet
                 .setLoadFlowModel(LoadFlowModel.DC)
                 .setShuntVoltageControl(false)
                 .setReactiveLimits(false)
-                .setHvdcAcEmulation(false) // FIXME
+                .setHvdcAcEmulation(parameters.isHvdcAcEmulation())
                 .setLowImpedanceThreshold(parametersExt.getLowImpedanceThreshold())
                 .setSvcVoltageMonitoring(false)
                 .setMaxSlackBusCount(1)

src/main/java/com/powsybl/openloadflow/dc/equations/AbstractHvdcDcEmulationFlowEquationTerm.java

@@ -0,0 +1,89 @@
+/**
+ * Copyright (c) 2022, RTE (http://www.rte-france.com)
+ * This Source Code Form is subject to the terms of the Mozilla Public
+ * License, v. 2.0. If a copy of the MPL was not distributed with this
+ * file, You can obtain one at http://mozilla.org/MPL/2.0/.
+ */
+package com.powsybl.openloadflow.dc.equations;
+
+import com.powsybl.openloadflow.equations.AbstractElementEquationTerm;
+import com.powsybl.openloadflow.equations.Variable;
+import com.powsybl.openloadflow.equations.VariableSet;
+import com.powsybl.openloadflow.network.LfBus;
+import com.powsybl.openloadflow.network.LfHvdc;
+
+import java.util.List;
+
+/**
+ * @author Anne Tilloy {@literal <anne.tilloy at rte-france.com>}
+ */
+public abstract class AbstractHvdcDcEmulationFlowEquationTerm extends AbstractElementEquationTerm<LfHvdc, DcVariableType, DcEquationType> {
+
+    protected final Variable<DcVariableType> ph1Var;
+
+    protected final Variable<DcVariableType> ph2Var;
+
+    protected final List<Variable<DcVariableType>> variables;
+
+    protected final double k;
+
+    protected final double p0;
+
+    protected final double lossFactor1;
+
+    protected final double lossFactor2;
+
+    protected final double pMaxFromCS1toCS2;
+
+    protected final double pMaxFromCS2toCS1;
+
+    protected AbstractHvdcDcEmulationFlowEquationTerm(LfHvdc hvdc, LfBus bus1, LfBus bus2, VariableSet<DcVariableType> variableSet) {
+        super(hvdc);
+        ph1Var = variableSet.getVariable(bus1.getNum(), DcVariableType.BUS_PHI);
+        ph2Var = variableSet.getVariable(bus2.getNum(), DcVariableType.BUS_PHI);
+        variables = List.of(ph1Var, ph2Var);
+        k = hvdc.getDroop() * 180 / Math.PI;
+        p0 = hvdc.getP0();
+        lossFactor1 = hvdc.getConverterStation1().getLossFactor() / 100;
+        lossFactor2 = hvdc.getConverterStation2().getLossFactor() / 100;
+        pMaxFromCS1toCS2 = hvdc.getPMaxFromCS1toCS2();
+        pMaxFromCS2toCS1 = hvdc.getPMaxFromCS2toCS1();
+    }
+
+    protected double rawP(double p0, double k, double ph1, double ph2) {
+        return p0 + k * (ph1 - ph2);
+    }
+
+    protected double boundedP(double rawP) {
+        // If there is a maximal active power
+        // it is applied at the entry of the controller VSC station
+        // on the AC side of the network.
+        if (rawP >= 0) {
+            return Math.min(rawP, pMaxFromCS1toCS2);
+        } else {
+            return Math.max(rawP, -pMaxFromCS2toCS1);
+        }
+    }
+
+    protected double ph1() {
+        return sv.get(ph1Var.getRow());
+    }
+
+    protected double ph2() {
+        return sv.get(ph2Var.getRow());
+    }
+
+    protected double getVscLossMultiplier() {
+        return (1 - lossFactor1) * (1 - lossFactor2);
+    }
+
+    @Override
+    public List<Variable<DcVariableType>> getVariables() {
+        return variables;
+    }
+
+    @Override
+    public boolean hasRhs() {
+        return false;
+    }
+}

src/main/java/com/powsybl/openloadflow/dc/equations/DcEquationSystemCreator.java

@@ -135,11 +135,37 @@ private void createBranches(EquationSystem<DcVariableType, DcEquationType> equat
         }
     }
 
+    private void createHvdcs(EquationSystem<DcVariableType, DcEquationType> equationSystem) {
+        for (LfHvdc hvdc : network.getHvdcs()) {
+            EquationTerm<DcVariableType, DcEquationType> p1 = null;
+            EquationTerm<DcVariableType, DcEquationType> p2 = null;
+            if (hvdc.getBus1() != null && hvdc.getBus2() != null && hvdc.isAcEmulation()) {
+                p1 = new HvdcDcEmulationSide1ActiveFlowEquationTerm(hvdc, hvdc.getBus1(), hvdc.getBus2(), equationSystem.getVariableSet());
+                p2 = new HvdcDcEmulationSide2ActiveFlowEquationTerm(hvdc, hvdc.getBus1(), hvdc.getBus2(), equationSystem.getVariableSet());
+            }
+
+            if (p1 != null) {
+                equationSystem.getEquation(hvdc.getBus1().getNum(), DcEquationType.BUS_TARGET_P)
+                        .orElseThrow()
+                        .addTerm(p1);
+                hvdc.setP1(p1);
+            }
+
+            if (p2 != null) {
+                equationSystem.getEquation(hvdc.getBus2().getNum(), DcEquationType.BUS_TARGET_P)
+                        .orElseThrow()
+                        .addTerm(p2);
+                hvdc.setP2(p2);
+            }
+        }
+    }
+
     public EquationSystem<DcVariableType, DcEquationType> create(boolean withListener) {
         EquationSystem<DcVariableType, DcEquationType> equationSystem = new EquationSystem<>();
 
         createBuses(equationSystem);
         createBranches(equationSystem);
+        createHvdcs(equationSystem);
 
         EquationSystemPostProcessor.findAll().forEach(pp -> pp.onCreate(equationSystem));
 

src/main/java/com/powsybl/openloadflow/dc/equations/HvdcDcEmulationSide1ActiveFlowEquationTerm.java

@@ -0,0 +1,73 @@
+/**
+ * Copyright (c) 2022, RTE (http://www.rte-france.com)
+ * This Source Code Form is subject to the terms of the Mozilla Public
+ * License, v. 2.0. If a copy of the MPL was not distributed with this
+ * file, You can obtain one at http://mozilla.org/MPL/2.0/.
+ */
+package com.powsybl.openloadflow.dc.equations;
+
+import com.powsybl.openloadflow.equations.Variable;
+import com.powsybl.openloadflow.equations.VariableSet;
+import com.powsybl.openloadflow.network.LfBus;
+import com.powsybl.openloadflow.network.LfHvdc;
+
+import java.util.Objects;
+
+/**
+ * @author Anne Tilloy {@literal <anne.tilloy at rte-france.com>}
+ */
+public class HvdcDcEmulationSide1ActiveFlowEquationTerm extends AbstractHvdcDcEmulationFlowEquationTerm {
+
+    public HvdcDcEmulationSide1ActiveFlowEquationTerm(LfHvdc hvdc, LfBus bus1, LfBus bus2, VariableSet<DcVariableType> variableSet) {
+        super(hvdc, bus1, bus2, variableSet);
+    }
+
+    private double p1(double ph1, double ph2) {
+        double rawP = rawP(p0, k, ph1, ph2);
+        double boundedP = boundedP(rawP);
+        return (isController(rawP) ? 1 : getVscLossMultiplier()) * boundedP;
+    }
+
+    private static boolean isController(double rawP) {
+        return rawP >= 0;
+    }
+
+    private boolean isInOperatingRange(double rawP) {
+        return rawP < pMaxFromCS1toCS2 && rawP > -pMaxFromCS2toCS1;
+    }
+
+    protected double dp1dph1(double ph1, double ph2) {
+        double rawP = rawP(p0, k, ph1, ph2);
+        if (isInOperatingRange(rawP)) {
+            return (isController(rawP) ? 1 : getVscLossMultiplier()) * k;
+        } else {
+            return 0;
+        }
+    }
+
+    protected double dp1dph2(double ph1, double ph2) {
+        return -dp1dph1(ph1, ph2);
+    }
+
+    @Override
+    public double eval() {
+        return p1(ph1(), ph2());
+    }
+
+    @Override
+    public double der(Variable<DcVariableType> variable) {
+        Objects.requireNonNull(variable);
+        if (variable.equals(ph1Var)) {
+            return dp1dph1(ph1(), ph2());
+        } else if (variable.equals(ph2Var)) {
+            return dp1dph2(ph1(), ph2());
+        } else {
+            throw new IllegalStateException("Unknown variable: " + variable);
+        }
+    }
+
+    @Override
+    protected String getName() {
+        return "ac_emulation_p_1";
+    }
+}

src/main/java/com/powsybl/openloadflow/dc/equations/HvdcDcEmulationSide2ActiveFlowEquationTerm.java

@@ -0,0 +1,73 @@
+/**
+ * Copyright (c) 2022, RTE (http://www.rte-france.com)
+ * This Source Code Form is subject to the terms of the Mozilla Public
+ * License, v. 2.0. If a copy of the MPL was not distributed with this
+ * file, You can obtain one at http://mozilla.org/MPL/2.0/.
+ */
+package com.powsybl.openloadflow.dc.equations;
+
+import com.powsybl.openloadflow.equations.Variable;
+import com.powsybl.openloadflow.equations.VariableSet;
+import com.powsybl.openloadflow.network.LfBus;
+import com.powsybl.openloadflow.network.LfHvdc;
+
+import java.util.Objects;
+
+/**
+ * @author Anne Tilloy {@literal <anne.tilloy at rte-france.com>}
+ */
+public class HvdcDcEmulationSide2ActiveFlowEquationTerm extends AbstractHvdcDcEmulationFlowEquationTerm {
+
+    public HvdcDcEmulationSide2ActiveFlowEquationTerm(LfHvdc hvdc, LfBus bus1, LfBus bus2, VariableSet<DcVariableType> variableSet) {
+        super(hvdc, bus1, bus2, variableSet);
+    }
+
+    private double p2(double ph1, double ph2) {
+        double rawP = rawP(p0, k, ph1, ph2);
+        double boundedP = boundedP(rawP);
+        return -(isController(rawP) ? 1 : getVscLossMultiplier()) * boundedP;
+    }
+
+    private boolean isController(double rawP) {
+        return rawP < 0;
+    }
+
+    private boolean isInOperatingRange(double rawP) {
+        return rawP < pMaxFromCS2toCS1 && rawP > -pMaxFromCS1toCS2;
+    }
+
+    private double dp2dph1(double ph1, double ph2) {
+        double rawP = rawP(p0, k, ph1, ph2);
+        if (isInOperatingRange(rawP)) {
+            return -(isController(rawP) ? 1 : getVscLossMultiplier()) * k;
+        } else {
+            return 0;
+        }
+    }
+
+    private double dp2dph2(double ph1, double ph2) {
+        return -dp2dph1(ph1, ph2);
+    }
+
+    @Override
+    public double eval() {
+        return p2(ph1(), ph2());
+    }
+
+    @Override
+    public double der(Variable<DcVariableType> variable) {
+        Objects.requireNonNull(variable);
+        if (variable.equals(ph1Var)) {
+            return dp2dph1(ph1(), ph2());
+        } else if (variable.equals(ph2Var)) {
+            return dp2dph2(ph1(), ph2());
+        } else {
+            throw new IllegalStateException("Unknown variable: " + variable);
+        }
+    }
+
+    @Override
+    protected String getName() {
+        return "ac_emulation_p_2";
+    }
+}

src/main/java/com/powsybl/openloadflow/sa/DcSecurityAnalysis.java

@@ -55,6 +55,7 @@ protected DcLoadFlowParameters createParameters(LoadFlowParameters lfParameters,
         var dcParameters = OpenLoadFlowParameters.createDcParameters(network, lfParameters,
                 lfParametersExt, matrixFactory, connectivityFactory, false);
         dcParameters.getNetworkParameters()
+                .setHvdcAcEmulation(false) // we only support it for DC load flow (initialization)
                 .setBreakers(breakers)
                 .setCacheEnabled(false) // force not caching as not supported in secu analysis
                 .setReferenceBusSelector(ReferenceBusSelector.DEFAULT_SELECTOR); // not supported yet

src/test/java/com/powsybl/openloadflow/OpenLoadFlowProviderTest.java

@@ -50,7 +50,7 @@ void test() {
     void testDcParameters() {
         Network network = Mockito.mock(Network.class);
         DcLoadFlowParameters dcParameters = OpenLoadFlowParameters.createDcParameters(network, new LoadFlowParameters().setReadSlackBus(true), new OpenLoadFlowParameters(), new DenseMatrixFactory(), new EvenShiloachGraphDecrementalConnectivityFactory<>(), true);
-        assertEquals("DcLoadFlowParameters(networkParameters=LfNetworkParameters(slackBusSelector=NetworkSlackBusSelector, connectivityFactory=EvenShiloachGraphDecrementalConnectivityFactory, generatorVoltageRemoteControl=false, minImpedance=false, twtSplitShuntAdmittance=false, breakers=false, plausibleActivePowerLimit=5000.0, computeMainConnectedComponentOnly=true, countriesToBalance=[], distributedOnConformLoad=false, phaseControl=false, transformerVoltageControl=false, voltagePerReactivePowerControl=false, generatorReactivePowerRemoteControl=false, transformerReactivePowerControl=false, loadFlowModel=DC, reactiveLimits=false, hvdcAcEmulation=false, minPlausibleTargetVoltage=0.8, maxPlausibleTargetVoltage=1.2, loaderPostProcessorSelection=[], reactiveRangeCheckMode=MAX, lowImpedanceThreshold=1.0E-8, svcVoltageMonitoring=false, maxSlackBusCount=1, debugDir=null, secondaryVoltageControl=false, cacheEnabled=false, asymmetrical=false, minNominalVoltageTargetVoltageCheck=20.0, linePerUnitMode=IMPEDANCE, useLoadModel=false, simulateAutomationSystems=false, referenceBusSelector=ReferenceBusFirstSlackSelector, voltageTargetPriorities=[GENERATOR, TRANSFORMER, SHUNT]), equationSystemCreationParameters=DcEquationSystemCreationParameters(updateFlows=true, forcePhaseControlOffAndAddAngle1Var=true, useTransformerRatio=true, dcApproximationType=IGNORE_R), matrixFactory=DenseMatrixFactory, distributedSlack=true, balanceType=PROPORTIONAL_TO_GENERATION_P_MAX, setVToNan=true, maxOuterLoopIterations=20)",
+        assertEquals("DcLoadFlowParameters(networkParameters=LfNetworkParameters(slackBusSelector=NetworkSlackBusSelector, connectivityFactory=EvenShiloachGraphDecrementalConnectivityFactory, generatorVoltageRemoteControl=false, minImpedance=false, twtSplitShuntAdmittance=false, breakers=false, plausibleActivePowerLimit=5000.0, computeMainConnectedComponentOnly=true, countriesToBalance=[], distributedOnConformLoad=false, phaseControl=false, transformerVoltageControl=false, voltagePerReactivePowerControl=false, generatorReactivePowerRemoteControl=false, transformerReactivePowerControl=false, loadFlowModel=DC, reactiveLimits=false, hvdcAcEmulation=true, minPlausibleTargetVoltage=0.8, maxPlausibleTargetVoltage=1.2, loaderPostProcessorSelection=[], reactiveRangeCheckMode=MAX, lowImpedanceThreshold=1.0E-8, svcVoltageMonitoring=false, maxSlackBusCount=1, debugDir=null, secondaryVoltageControl=false, cacheEnabled=false, asymmetrical=false, minNominalVoltageTargetVoltageCheck=20.0, linePerUnitMode=IMPEDANCE, useLoadModel=false, simulateAutomationSystems=false, referenceBusSelector=ReferenceBusFirstSlackSelector, voltageTargetPriorities=[GENERATOR, TRANSFORMER, SHUNT]), equationSystemCreationParameters=DcEquationSystemCreationParameters(updateFlows=true, forcePhaseControlOffAndAddAngle1Var=true, useTransformerRatio=true, dcApproximationType=IGNORE_R), matrixFactory=DenseMatrixFactory, distributedSlack=true, balanceType=PROPORTIONAL_TO_GENERATION_P_MAX, setVToNan=true, maxOuterLoopIterations=20)",
                 dcParameters.toString());
     }
 

src/test/java/com/powsybl/openloadflow/ac/AcLoadFlowVscTest.java

@@ -208,7 +208,7 @@ void testHvdcAcEmulationNonSupported() {
     }
 
     @Test
-    void testHvdcAcEmulationNonSupported2() {
+    void testDcLoadFlowWithHvdcAcEmulation() {
         Network network = HvdcNetworkFactory.createWithHvdcInAcEmulation();
         network.getHvdcLine("hvdc34").newExtension(HvdcAngleDroopActivePowerControlAdder.class)
                 .withDroop(180)
@@ -227,10 +227,10 @@ void testHvdcAcEmulationNonSupported2() {
         assertTrue(result.isFullyConverged());
 
         VscConverterStation cs3 = network.getVscConverterStation("cs3");
-        assertActivePowerEquals(-1.956, cs3.getTerminal());
+        assertActivePowerEquals(-0.09, cs3.getTerminal());
 
         VscConverterStation cs4 = network.getVscConverterStation("cs4");
-        assertActivePowerEquals(2.0, cs4.getTerminal());
+        assertActivePowerEquals(0.092, cs4.getTerminal());
     }
 
     @Test