Refactor equation vector

src/main/java/com/powsybl/openloadflow/ac/nr/NewtonRaphson.java

@@ -13,10 +13,15 @@
 import com.powsybl.openloadflow.network.LfElement;
 import com.powsybl.openloadflow.network.LfNetwork;
 import com.powsybl.openloadflow.network.util.VoltageInitializer;
+import org.apache.commons.lang3.tuple.Pair;
 import org.slf4j.Logger;
 import org.slf4j.LoggerFactory;
 
+import java.util.Comparator;
+import java.util.List;
+import java.util.Map;
 import java.util.Objects;
+import java.util.stream.Collectors;
 
 /**
  * @author Geoffroy Jamgotchian <geoffroy.jamgotchian at rte-france.com>
@@ -37,38 +42,52 @@ public class NewtonRaphson {
 
     private final TargetVector<AcVariableType, AcEquationType> targetVector;
 
+    private final EquationVector<AcVariableType, AcEquationType> equationVector;
+
     public NewtonRaphson(LfNetwork network, NewtonRaphsonParameters parameters,
-                         EquationSystem<AcVariableType, AcEquationType> equationSystem, JacobianMatrix<AcVariableType, AcEquationType> j,
-                         TargetVector<AcVariableType, AcEquationType> targetVector) {
+                         EquationSystem<AcVariableType, AcEquationType> equationSystem,
+                         JacobianMatrix<AcVariableType, AcEquationType> j,
+                         TargetVector<AcVariableType, AcEquationType> targetVector,
+                         EquationVector<AcVariableType, AcEquationType> equationVector) {
         this.network = Objects.requireNonNull(network);
         this.parameters = Objects.requireNonNull(parameters);
         this.equationSystem = Objects.requireNonNull(equationSystem);
         this.j = Objects.requireNonNull(j);
         this.targetVector = Objects.requireNonNull(targetVector);
+        this.equationVector = Objects.requireNonNull(equationVector);
+    }
+
+    public static List<Pair<Equation<AcVariableType, AcEquationType>, Double>> findLargestMismatches(EquationSystem<AcVariableType, AcEquationType> equationSystem, double[] mismatch, int count) {
+        return equationSystem.getIndex().getSortedEquationsToSolve().stream()
+                .map(equation -> Pair.of(equation, mismatch[equation.getColumn()]))
+                .filter(e -> Math.abs(e.getValue()) > Math.pow(10, -7))
+                .sorted(Comparator.comparingDouble((Map.Entry<Equation<AcVariableType, AcEquationType>, Double> e) -> Math.abs(e.getValue())).reversed())
+                .limit(count)
+                .collect(Collectors.toList());
     }
 
-    private NewtonRaphsonStatus runIteration(double[] fx) {
+    private NewtonRaphsonStatus runIteration() {
         LOGGER.debug("Start iteration {}", iteration);
 
         try {
             // solve f(x) = j * dx
             try {
-                j.solveTransposed(fx);
+                j.solveTransposed(equationVector.getArray());
             } catch (MatrixException e) {
                 LOGGER.error(e.toString(), e);
                 return NewtonRaphsonStatus.SOLVER_FAILED;
             }
+            // f(x) now contains dx
 
-            // update x
-            equationSystem.getStateVector().minus(fx);
+            // update x and f(x) will be automatically updated
+            equationSystem.getStateVector().minus(equationVector.getArray());
 
-            // recalculate f(x) with new x
-            equationSystem.updateEquationVector(fx);
-
-            Vectors.minus(fx, targetVector.toArray());
+            // substract targets from f(x)
+            equationVector.minus(targetVector);
+            // f(x) now contains equation mismatches
 
             if (LOGGER.isTraceEnabled()) {
-                equationSystem.findLargestMismatches(fx, 5)
+                findLargestMismatches(equationSystem, equationVector.getArray(), 5)
                         .forEach(e -> {
                             Equation<AcVariableType, AcEquationType> equation = e.getKey();
                             String elementId = equation.getElement(network).map(LfElement::getId).orElse("?");
@@ -77,7 +96,7 @@ private NewtonRaphsonStatus runIteration(double[] fx) {
             }
 
             // test stopping criteria and log norm(fx)
-            NewtonRaphsonStoppingCriteria.TestResult testResult = parameters.getStoppingCriteria().test(fx);
+            NewtonRaphsonStoppingCriteria.TestResult testResult = parameters.getStoppingCriteria().test(equationVector.getArray());
 
             LOGGER.debug("|f(x)|={}", testResult.getNorm());
 
@@ -168,15 +187,12 @@ public NewtonRaphsonResult run(VoltageInitializer voltageInitializer) {
         // initialize state vector
         initStateVector(network, equationSystem, voltageInitializer);
 
-        // initialize mismatch vector (difference between equation values and targets)
-        double[] fx = equationSystem.createEquationVector();
-
-        Vectors.minus(fx, targetVector.toArray());
+        Vectors.minus(equationVector.getArray(), targetVector.getArray());
 
         // start iterations
         NewtonRaphsonStatus status = NewtonRaphsonStatus.NO_CALCULATION;
         while (iteration <= parameters.getMaxIteration()) {
-            NewtonRaphsonStatus newStatus = runIteration(fx);
+            NewtonRaphsonStatus newStatus = runIteration();
             if (newStatus != null) {
                 status = newStatus;
                 break;

src/main/java/com/powsybl/openloadflow/ac/outerloop/AcLoadFlowContext.java

@@ -10,6 +10,7 @@
 import com.powsybl.openloadflow.ac.equations.AcEquationType;
 import com.powsybl.openloadflow.ac.equations.AcVariableType;
 import com.powsybl.openloadflow.equations.EquationSystem;
+import com.powsybl.openloadflow.equations.EquationVector;
 import com.powsybl.openloadflow.equations.JacobianMatrix;
 import com.powsybl.openloadflow.equations.TargetVector;
 import com.powsybl.openloadflow.network.LfNetwork;
@@ -31,6 +32,8 @@ public class AcLoadFlowContext implements AutoCloseable {
 
     private AcTargetVector targetVector;
 
+    private EquationVector<AcVariableType, AcEquationType> equationVector;
+
     public AcLoadFlowContext(LfNetwork network, AcLoadFlowParameters parameters) {
         this.network = Objects.requireNonNull(network);
         this.parameters = Objects.requireNonNull(parameters);
@@ -65,6 +68,13 @@ public TargetVector<AcVariableType, AcEquationType> getTargetVector() {
         return targetVector;
     }
 
+    public EquationVector<AcVariableType, AcEquationType> getEquationVector() {
+        if (equationVector == null) {
+            equationVector = new EquationVector<>(getEquationSystem());
+        }
+        return equationVector;
+    }
+
     @Override
     public void close() {
         if (jacobianMatrix != null) {

src/main/java/com/powsybl/openloadflow/ac/outerloop/AcloadFlowEngine.java

@@ -88,8 +88,12 @@ public AcLoadFlowResult run() {
         voltageInitializer.prepare(context.getNetwork());
 
         RunningContext runningContext = new RunningContext();
-        NewtonRaphson newtonRaphson = new NewtonRaphson(context.getNetwork(), context.getParameters().getNewtonRaphsonParameters(),
-                context.getEquationSystem(), context.getJacobianMatrix(), context.getTargetVector());
+        NewtonRaphson newtonRaphson = new NewtonRaphson(context.getNetwork(),
+                                                        context.getParameters().getNewtonRaphsonParameters(),
+                                                        context.getEquationSystem(),
+                                                        context.getJacobianMatrix(),
+                                                        context.getTargetVector(),
+                                                        context.getEquationVector());
 
         List<OuterLoop> outerLoops = context.getParameters().getOuterLoops();
         List<Pair<OuterLoop, OuterLoopContextImpl>> outerLoopsAndContexts = outerLoops.stream()

src/main/java/com/powsybl/openloadflow/equations/AbstractVector.java

@@ -0,0 +1,94 @@
+/**
+ * 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.equations;
+
+import java.util.Objects;
+
+/**
+ * @author Geoffroy Jamgotchian <geoffroy.jamgotchian at rte-france.com>
+ */
+public abstract class AbstractVector<V extends Enum<V> & Quantity, E extends Enum<E> & Quantity> implements Vector {
+
+    protected final EquationSystem<V, E> equationSystem;
+
+    private double[] array;
+
+    private enum Status {
+        VALID,
+        VECTOR_INVALID, // structure has changed
+        VALUES_INVALID // same structure but values have to be updated
+    }
+
+    private Status status = Status.VECTOR_INVALID;
+
+    private final EquationSystemIndexListener<V, E> equationSystemIndexListener = new EquationSystemIndexListener<>() {
+
+        @Override
+        public void onEquationChange(Equation<V, E> equation, ChangeType changeType) {
+            invalidateVector();
+        }
+
+        @Override
+        public void onVariableChange(Variable<V> variable, ChangeType changeType) {
+            // nothing to do
+        }
+
+        @Override
+        public void onEquationTermChange(EquationTerm<V, E> term) {
+            // nothing to do
+        }
+    };
+
+    protected AbstractVector(EquationSystem<V, E> equationSystem) {
+        this.equationSystem = Objects.requireNonNull(equationSystem);
+        equationSystem.getIndex().addListener(equationSystemIndexListener);
+    }
+
+    protected void invalidateValues() {
+        if (status == Status.VALID) {
+            status = Status.VALUES_INVALID;
+        }
+    }
+
+    protected void invalidateVector() {
+        status = Status.VECTOR_INVALID;
+    }
+
+    protected void validate() {
+        status = Status.VALID;
+    }
+
+    @Override
+    public double[] getArray() {
+        switch (status) {
+            case VECTOR_INVALID:
+                array = createArray();
+                validate();
+                break;
+
+            case VALUES_INVALID:
+                updateArray(array);
+                validate();
+                break;
+
+            default:
+                // nothing to do
+                break;
+        }
+        return array;
+    }
+
+    protected abstract double[] createArray();
+
+    protected abstract void updateArray(double[] array);
+
+    @Override
+    public void minus(Vector other) {
+        Objects.requireNonNull(other);
+        Vectors.minus(getArray(), other.getArray());
+    }
+}

src/main/java/com/powsybl/openloadflow/equations/EquationSystem.java

@@ -6,31 +6,23 @@
  */
 package com.powsybl.openloadflow.equations;
 
-import com.google.common.base.Stopwatch;
 import com.powsybl.commons.PowsyblException;
 import com.powsybl.openloadflow.network.ElementType;
 import com.powsybl.openloadflow.network.LfNetwork;
 import org.apache.commons.lang3.tuple.Pair;
-import org.slf4j.Logger;
-import org.slf4j.LoggerFactory;
 
 import java.io.IOException;
 import java.io.StringWriter;
 import java.io.UncheckedIOException;
 import java.io.Writer;
 import java.util.*;
-import java.util.concurrent.TimeUnit;
 import java.util.stream.Collectors;
 
-import static com.powsybl.openloadflow.util.Markers.PERFORMANCE_MARKER;
-
 /**
  * @author Geoffroy Jamgotchian <geoffroy.jamgotchian at rte-france.com>
  */
 public class EquationSystem<V extends Enum<V> & Quantity, E extends Enum<E> & Quantity> {
 
-    private static final Logger LOGGER = LoggerFactory.getLogger(EquationSystem.class);
-
     private final boolean indexTerms;
 
     private final Map<Pair<Integer, E>, Equation<V, E>> equations = new HashMap<>();
@@ -172,27 +164,6 @@ public List<String> getColumnNames(LfNetwork network) {
                 .collect(Collectors.toList());
     }
 
-    public double[] createEquationVector() {
-        double[] fx = new double[index.getSortedEquationsToSolve().size()];
-        updateEquationVector(fx);
-        return fx;
-    }
-
-    public void updateEquationVector(double[] fx) {
-        if (fx.length != index.getSortedEquationsToSolve().size()) {
-            throw new IllegalArgumentException("Bad equation vector length: " + fx.length);
-        }
-
-        Stopwatch stopwatch = Stopwatch.createStarted();
-
-        Arrays.fill(fx, 0);
-        for (Equation<V, E> equation : index.getSortedEquationsToSolve()) {
-            fx[equation.getColumn()] = equation.eval();
-        }
-
-        LOGGER.debug(PERFORMANCE_MARKER, "Equation vector updated in {} us", stopwatch.elapsed(TimeUnit.MICROSECONDS));
-    }
-
     public void addListener(EquationSystemListener<V, E> listener) {
         Objects.requireNonNull(listener);
         listeners.add(listener);
@@ -247,13 +218,4 @@ public String writeToString(boolean writeInactiveEquations) {
     public String writeToString() {
         return writeToString(false);
     }
-
-    public List<Pair<Equation<V, E>, Double>> findLargestMismatches(double[] mismatch, int count) {
-        return index.getSortedEquationsToSolve().stream()
-                .map(equation -> Pair.of(equation, mismatch[equation.getColumn()]))
-                .filter(e -> Math.abs(e.getValue()) > Math.pow(10, -7))
-                .sorted(Comparator.comparingDouble((Map.Entry<Equation<V, E>, Double> e) -> Math.abs(e.getValue())).reversed())
-                .limit(count)
-                .collect(Collectors.toList());
-    }
 }

src/main/java/com/powsybl/openloadflow/equations/EquationVector.java

@@ -0,0 +1,60 @@
+/**
+ * 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.equations;
+
+import com.google.common.base.Stopwatch;
+import org.slf4j.Logger;
+import org.slf4j.LoggerFactory;
+
+import java.util.Arrays;
+import java.util.concurrent.TimeUnit;
+
+import static com.powsybl.openloadflow.util.Markers.PERFORMANCE_MARKER;
+
+/**
+ * @author Geoffroy Jamgotchian <geoffroy.jamgotchian at rte-france.com>
+ */
+public class EquationVector<V extends Enum<V> & Quantity, E extends Enum<E> & Quantity> extends AbstractVector<V, E>
+        implements StateVectorListener {
+
+    private static final Logger LOGGER = LoggerFactory.getLogger(EquationVector.class);
+
+    public EquationVector(EquationSystem<V, E> equationSystem) {
+        super(equationSystem);
+        equationSystem.getStateVector().addListener(this);
+    }
+
+    @Override
+    public void onStateUpdate() {
+        invalidateValues();
+    }
+
+    @Override
+    protected double[] createArray() {
+        double[] array = new double[equationSystem.getIndex().getSortedEquationsToSolve().size()];
+        updateArray(array);
+        return array;
+    }
+
+    @Override
+    protected void updateArray(double[] array) {
+        Stopwatch stopwatch = Stopwatch.createStarted();
+
+        var equations = equationSystem.getIndex().getSortedEquationsToSolve();
+
+        if (array.length != equations.size()) {
+            throw new IllegalArgumentException("Bad equation vector length: " + array.length);
+        }
+
+        Arrays.fill(array, 0); // necessary?
+        for (Equation<V, E> equation : equations) {
+            array[equation.getColumn()] = equation.eval();
+        }
+
+        LOGGER.debug(PERFORMANCE_MARKER, "Equation vector updated in {} us", stopwatch.elapsed(TimeUnit.MICROSECONDS));
+    }
+}