[PoC] Robust optimization with Marmot

ra-optimisation/search-tree-rao/src/main/java/com/powsybl/openrao/searchtreerao/commons/SensitivityComputer.java

@@ -63,6 +63,7 @@ public static final class SensitivityComputerBuilder {
         private Set<FlowCnec> loopFlowCnecs;
         private AppliedRemedialActions appliedRemedialActions;
         private Instant outageInstant;
+        private Set<String> extraInjectionIds;
 
         public SensitivityComputerBuilder withToolProvider(ToolProvider toolProvider) {
             this.toolProvider = toolProvider;
@@ -114,6 +115,11 @@ public SensitivityComputerBuilder withOutageInstant(Instant outageInstant) {
             return this;
         }
 
+        public SensitivityComputerBuilder withExtraInjections(Set<String> extraInjectionIds) {
+            this.extraInjectionIds = extraInjectionIds;
+            return this;
+        }
+
         public SensitivityComputer build() {
             Objects.requireNonNull(toolProvider);
             Objects.requireNonNull(flowCnecs);
@@ -128,7 +134,8 @@ public SensitivityComputer build() {
                     computePtdfs,
                     computeLoopFlows,
                     appliedRemedialActions,
-                    outageInstant);
+                    outageInstant,
+                    extraInjectionIds);
             BranchResultAdapterImpl.BranchResultAdpaterBuilder builder = BranchResultAdapterImpl.create();
             if (loopFlowComputation != null) {
                 builder.withCommercialFlowsResults(loopFlowComputation, loopFlowCnecs);

ra-optimisation/search-tree-rao/src/main/java/com/powsybl/openrao/searchtreerao/commons/ToolProvider.java

@@ -21,7 +21,9 @@
 import com.powsybl.openrao.raoapi.parameters.RaoParameters;
 import com.powsybl.openrao.raoapi.parameters.LoopFlowParameters;
 import com.powsybl.openrao.raoapi.parameters.RelativeMarginsParameters;
+import com.powsybl.openrao.searchtreerao.marmot.Marmot;
 import com.powsybl.openrao.sensitivityanalysis.AppliedRemedialActions;
+import com.powsybl.openrao.sensitivityanalysis.InjectionSensitivityProvider;
 import com.powsybl.openrao.sensitivityanalysis.SystematicSensitivityInterface;
 import com.powsybl.glsk.commons.ZonalData;
 import com.powsybl.glsk.commons.ZonalDataImpl;
@@ -83,15 +85,16 @@ public SystematicSensitivityInterface getSystematicSensitivityInterface(Set<Flow
                                                                             Set<RangeAction<?>> rangeActions,
                                                                             boolean computePtdfs,
                                                                             boolean computeLoopFlows, Instant outageInstant) {
-        return getSystematicSensitivityInterface(cnecs, rangeActions, computePtdfs, computeLoopFlows, null, outageInstant);
+        return getSystematicSensitivityInterface(cnecs, rangeActions, computePtdfs, computeLoopFlows, null, outageInstant, null);
     }
 
     public SystematicSensitivityInterface getSystematicSensitivityInterface(Set<FlowCnec> cnecs,
                                                                             Set<RangeAction<?>> rangeActions,
                                                                             boolean computePtdfs,
                                                                             boolean computeLoopFlows,
                                                                             AppliedRemedialActions appliedRemedialActions,
-                                                                            Instant outageInstant) {
+                                                                            Instant outageInstant,
+                                                                            Set<String> extraInjectionIds) {
 
         SystematicSensitivityInterface.SystematicSensitivityInterfaceBuilder builder = SystematicSensitivityInterface.builder()
             .withSensitivityProviderName(getSensitivityProvider(raoParameters))
@@ -100,6 +103,16 @@ public SystematicSensitivityInterface getSystematicSensitivityInterface(Set<Flow
             .withAppliedRemedialActions(appliedRemedialActions)
             .withOutageInstant(outageInstant);
 
+        if (extraInjectionIds != null && !extraInjectionIds.isEmpty()) {
+            builder.withSensitivityProvider(new InjectionSensitivityProvider(cnecs, extraInjectionIds, Collections.singleton(Unit.MEGAWATT)));
+        } else if (Marmot.getScenarioRepo() != null) {
+            // TODO remove this and make it cleaner
+            Set<String> injectionIds = Marmot.getScenarioRepo().getInjectionIds();
+            if (!injectionIds.isEmpty()) {
+                builder.withSensitivityProvider(new InjectionSensitivityProvider(cnecs, injectionIds, Collections.singleton(Unit.MEGAWATT)));
+            }
+        }
+
         if (!getSensitivityWithLoadFlowParameters(raoParameters).getLoadFlowParameters().isDc()) {
             builder.withLoadflow(cnecs, Collections.singleton(Unit.AMPERE));
         }

ra-optimisation/search-tree-rao/src/main/java/com/powsybl/openrao/searchtreerao/linearoptimisation/algorithms/ShiftedFlowResult.java

@@ -0,0 +1,81 @@
+/*
+ * Copyright (c) 2025, 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.openrao.searchtreerao.linearoptimisation.algorithms;
+
+import com.powsybl.iidm.network.TwoSides;
+import com.powsybl.openrao.commons.OpenRaoException;
+import com.powsybl.openrao.commons.Unit;
+import com.powsybl.openrao.data.crac.api.Instant;
+import com.powsybl.openrao.data.crac.api.State;
+import com.powsybl.openrao.data.crac.api.cnec.FlowCnec;
+import com.powsybl.openrao.data.raoresult.api.ComputationStatus;
+import com.powsybl.openrao.searchtreerao.result.api.FlowResult;
+import com.powsybl.openrao.searchtreerao.result.api.SensitivityResult;
+
+import java.util.Map;
+
+/**
+ * @author Peter Mitri {@literal <peter.mitri at rte-france.com>}
+ */
+public class ShiftedFlowResult implements FlowResult {
+    private FlowResult initFlowResult;
+    private Map<String, Double> shiftedInjections;
+    private SensitivityResult sensitivityResult;
+
+    public ShiftedFlowResult(FlowResult initFlowResult, Map<String, Double> shiftedInjections, SensitivityResult sensitivityResult) {
+        this.initFlowResult = initFlowResult;
+        this.shiftedInjections = shiftedInjections;
+        this.sensitivityResult = sensitivityResult;
+    }
+
+    @Override
+    public double getFlow(FlowCnec flowCnec, TwoSides side, Unit unit) {
+        double flow = initFlowResult.getFlow(flowCnec, side, unit);
+        for (String shiftedInjection : this.shiftedInjections.keySet()) {
+            double sensi = sensitivityResult.getSensitivityValue(flowCnec, side, shiftedInjection, unit);
+            flow += sensi * shiftedInjections.get(shiftedInjection);
+        }
+        return flow;
+    }
+
+    @Override
+    public double getFlow(FlowCnec flowCnec, TwoSides side, Unit unit, Instant optimizedInstant) {
+        throw new OpenRaoException("Not implemented");
+    }
+
+    @Override
+    public double getMargin(FlowCnec flowCnec, Unit unit) {
+        throw new OpenRaoException("Not implemented");
+    }
+
+    @Override
+    public double getCommercialFlow(FlowCnec flowCnec, TwoSides side, Unit unit) {
+        // throw new OpenRaoException("Not implemented");
+        return initFlowResult.getCommercialFlow(flowCnec, side, unit);
+    }
+
+    @Override
+    public double getPtdfZonalSum(FlowCnec flowCnec, TwoSides side) {
+        return initFlowResult.getPtdfZonalSum(flowCnec, side);
+    }
+
+    @Override
+    public Map<FlowCnec, Map<TwoSides, Double>> getPtdfZonalSums() {
+        return initFlowResult.getPtdfZonalSums();
+    }
+
+    @Override
+    public ComputationStatus getComputationStatus() {
+        return initFlowResult.getComputationStatus();
+    }
+
+    @Override
+    public ComputationStatus getComputationStatus(State state) {
+        return initFlowResult.getComputationStatus(state);
+    }
+}

ra-optimisation/search-tree-rao/src/main/java/com/powsybl/openrao/searchtreerao/linearoptimisation/algorithms/fillers/MaxLoopFlowFiller.java

@@ -28,6 +28,8 @@
 import java.util.Set;
 import java.util.TreeSet;
 
+import static com.powsybl.openrao.commons.Unit.MEGAWATT;
+
 /**
  * @author Pengbo Wang {@literal <pengbo.wang at rte-international.com>}
  */
@@ -116,6 +118,7 @@ private void buildLoopFlowConstraintsAndUpdateObjectiveFunction(LinearProblem li
                     side,
                     Optional.ofNullable(timestamp)
                 );
+                double referenceFlow = flowResult.getFlow(cnec, side, MEGAWATT);
 
                 // build constraint which defines the loopFlow :
                 // - MaxLoopFlow + commercialFlow <= flowVariable + loopflowViolationVariable <= POSITIVE_INF
@@ -132,6 +135,7 @@ private void buildLoopFlowConstraintsAndUpdateObjectiveFunction(LinearProblem li
                 );
                 positiveLoopflowViolationConstraint.setCoefficient(flowVariable, 1);
                 positiveLoopflowViolationConstraint.setCoefficient(loopflowViolationVariable, 1.0);
+                positiveLoopflowViolationConstraint.setIsLazy(referenceFlow >= -loopFlowUpperBound + flowResult.getCommercialFlow(cnec, side, Unit.MEGAWATT));
 
                 OpenRaoMPConstraint negativeLoopflowViolationConstraint = linearProblem.addMaxLoopFlowConstraint(
                     -linearProblem.infinity(),
@@ -143,6 +147,7 @@ private void buildLoopFlowConstraintsAndUpdateObjectiveFunction(LinearProblem li
                 );
                 negativeLoopflowViolationConstraint.setCoefficient(flowVariable, 1);
                 negativeLoopflowViolationConstraint.setCoefficient(loopflowViolationVariable, -1);
+                negativeLoopflowViolationConstraint.setIsLazy(referenceFlow <= loopFlowUpperBound + flowResult.getCommercialFlow(cnec, side, Unit.MEGAWATT));
 
                 //update objective function with loopflowViolationCost
                 linearProblem.getObjective().setCoefficient(loopflowViolationVariable, loopFlowViolationCost / cnec.getMonitoredSides().size());

ra-optimisation/search-tree-rao/src/main/java/com/powsybl/openrao/searchtreerao/linearoptimisation/algorithms/fillers/MaxMinMarginFiller.java

@@ -20,10 +20,7 @@
 import com.powsybl.openrao.searchtreerao.result.api.SensitivityResult;
 
 import java.time.OffsetDateTime;
-import java.util.Comparator;
-import java.util.Optional;
-import java.util.Set;
-import java.util.TreeSet;
+import java.util.*;
 
 import static com.powsybl.openrao.commons.Unit.MEGAWATT;
 
@@ -61,7 +58,7 @@ public void fill(LinearProblem linearProblem, FlowResult flowResult, Sensitivity
         }
 
         // build constraints
-        buildMinimumMarginConstraints(linearProblem, validFlowCnecs);
+        buildMinimumMarginConstraints(linearProblem, validFlowCnecs, flowResult);
         if (costOptimization) {
             addMinMarginShiftedViolationConstraint(linearProblem);
         }
@@ -121,30 +118,42 @@ private void buildMinimumMarginVariable(LinearProblem linearProblem, Set<FlowCne
      * MM <= (fmax[c] - F[c]) * 1000 / (Unom * sqrt(3))     (ABOVE_THRESHOLD)
      * MM <= (F[c] - fmin[c]) * 1000 / (Unom * sqrt(3))     (BELOW_THRESHOLD)
      */
-    private void buildMinimumMarginConstraints(LinearProblem linearProblem, Set<FlowCnec> validFlowCnecs) {
+    private void buildMinimumMarginConstraints(LinearProblem linearProblem, Set<FlowCnec> validFlowCnecs, FlowResult flowResult) {
         OpenRaoMPVariable minimumMarginVariable = linearProblem.getMinimumMarginVariable(Optional.ofNullable(timestamp));
 
+        List<OpenRaoMPConstraint> allMarginConstraints = new ArrayList<>();
         validFlowCnecs.forEach(cnec -> cnec.getMonitoredSides().forEach(side -> {
+
             OpenRaoMPVariable flowVariable = linearProblem.getFlowVariable(cnec, side, Optional.ofNullable(timestamp));
 
             Optional<Double> minFlow;
             Optional<Double> maxFlow;
             minFlow = cnec.getLowerBound(side, MEGAWATT);
             maxFlow = cnec.getUpperBound(side, MEGAWATT);
             double unitConversionCoefficient = RaoUtil.getFlowUnitMultiplier(cnec, side, unit, MEGAWATT);
+            double referenceFlow = flowResult.getFlow(cnec, side, unit) * unitConversionCoefficient;
 
             if (minFlow.isPresent()) {
                 OpenRaoMPConstraint minimumMarginNegative = linearProblem.addMinimumMarginConstraint(-linearProblem.infinity(), -minFlow.get(), cnec, side, LinearProblem.MarginExtension.BELOW_THRESHOLD, Optional.ofNullable(timestamp));
                 minimumMarginNegative.setCoefficient(minimumMarginVariable, unitConversionCoefficient);
                 minimumMarginNegative.setCoefficient(flowVariable, -1);
+                minimumMarginNegative.setIsLazy(referenceFlow >= minFlow.get());
+                allMarginConstraints.add(minimumMarginNegative);
             }
 
             if (maxFlow.isPresent()) {
                 OpenRaoMPConstraint minimumMarginPositive = linearProblem.addMinimumMarginConstraint(-linearProblem.infinity(), maxFlow.get(), cnec, side, LinearProblem.MarginExtension.ABOVE_THRESHOLD, Optional.ofNullable(timestamp));
                 minimumMarginPositive.setCoefficient(minimumMarginVariable, unitConversionCoefficient);
                 minimumMarginPositive.setCoefficient(flowVariable, 1);
+                minimumMarginPositive.setIsLazy(referenceFlow <= maxFlow.get());
+                allMarginConstraints.add(minimumMarginPositive);
             }
         }));
+
+        // Make sure at least one constraint is not lazy, otherwise XPRESS finds the problem to be dual infeasible
+        if (allMarginConstraints.stream().allMatch(OpenRaoMPConstraint::isLazy)) {
+            allMarginConstraints.stream().findAny().ifPresent(constraint -> constraint.setIsLazy(false));
+        }
     }
 
     /**

ra-optimisation/search-tree-rao/src/main/java/com/powsybl/openrao/searchtreerao/linearoptimisation/algorithms/fillers/MnecFiller.java

@@ -59,7 +59,7 @@ public MnecFiller(FlowResult initialFlowResult,
     public void fill(LinearProblem linearProblem, FlowResult flowResult, SensitivityResult sensitivityResult, RangeActionActivationResult rangeActionActivationResult) {
         Set<FlowCnec> validMonitoredCnecs = FillersUtil.getFlowCnecsComputationStatusOk(monitoredCnecs, sensitivityResult);
         buildMarginViolationVariable(linearProblem, validMonitoredCnecs);
-        buildMnecMarginConstraints(linearProblem, validMonitoredCnecs);
+        buildMnecMarginConstraints(linearProblem, validMonitoredCnecs, flowResult);
         fillObjectiveWithMnecPenaltyCost(linearProblem, validMonitoredCnecs);
     }
 
@@ -74,19 +74,21 @@ private void buildMarginViolationVariable(LinearProblem linearProblem, Set<FlowC
         ));
     }
 
-    private void buildMnecMarginConstraints(LinearProblem linearProblem, Set<FlowCnec> validMonitoredCnecs) {
+    private void buildMnecMarginConstraints(LinearProblem linearProblem, Set<FlowCnec> validMonitoredCnecs, FlowResult flowResult) {
         validMonitoredCnecs.forEach(mnec -> mnec.getMonitoredSides().forEach(side -> {
                 double mnecInitialFlowInMW = initialFlowResult.getFlow(mnec, side, unit) * RaoUtil.getFlowUnitMultiplier(mnec, side, unit, MEGAWATT);
 
                 OpenRaoMPVariable flowVariable = linearProblem.getFlowVariable(mnec, side, Optional.ofNullable(timestamp));
                 OpenRaoMPVariable mnecViolationVariable = linearProblem.getMnecViolationVariable(mnec, side, Optional.ofNullable(timestamp));
+                double referenceFlow = flowResult.getFlow(mnec, side, unit) * RaoUtil.getFlowUnitMultiplier(mnec, side, unit, MEGAWATT);
 
                 Optional<Double> maxFlow = mnec.getUpperBound(side, MEGAWATT);
                 if (maxFlow.isPresent()) {
                     double ub = Math.max(maxFlow.get(), mnecInitialFlowInMW + mnecAcceptableMarginDecrease) - mnecConstraintAdjustmentCoefficient;
                     OpenRaoMPConstraint maxConstraint = linearProblem.addMnecFlowConstraint(-linearProblem.infinity(), ub, mnec, side, LinearProblem.MarginExtension.BELOW_THRESHOLD, Optional.ofNullable(timestamp));
                     maxConstraint.setCoefficient(flowVariable, 1);
                     maxConstraint.setCoefficient(mnecViolationVariable, -1);
+                    maxConstraint.setIsLazy(referenceFlow <= ub);
                 }
 
                 Optional<Double> minFlow = mnec.getLowerBound(side, MEGAWATT);
@@ -95,6 +97,7 @@ private void buildMnecMarginConstraints(LinearProblem linearProblem, Set<FlowCne
                     OpenRaoMPConstraint maxConstraint = linearProblem.addMnecFlowConstraint(lb, linearProblem.infinity(), mnec, side, LinearProblem.MarginExtension.ABOVE_THRESHOLD, Optional.ofNullable(timestamp));
                     maxConstraint.setCoefficient(flowVariable, 1);
                     maxConstraint.setCoefficient(mnecViolationVariable, 1);
+                    maxConstraint.setIsLazy(referenceFlow >= lb);
                 }
             }
         ));