Respect for slack voltage angle in distributed power flow calculation

CHANGELOG.md

@@ -50,6 +50,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 - Fixed some unreachable code [#167](https://github.com/ie3-institute/simona/issues/167)
 - Fix treatment of non-InitializeTrigger triggers in initialization within SimScheduler [#237](https://github.com/ie3-institute/simona/issues/237)
 - Fix breaking SIMONA caused by introducing temperature dependant load profiles in PSDM [#255](https://github.com/ie3-institute/simona/issues/255)
+- Respect for voltage angle in DBFS slack voltage exchange protocol [#69](https://github.com/ie3-institute/simona/issues/69)
 
 ### Removed
 - Remove workaround for tscfg tmp directory [#178](https://github.com/ie3-institute/simona/issues/178)

src/main/scala/edu/ie3/simona/agent/grid/DBFSAlgorithm.scala

@@ -433,21 +433,23 @@ trait DBFSAlgorithm extends PowerFlowSupport with GridResultsSupport {
       )
 
       val gridModel = gridAgentBaseData.gridEnv.gridModel
-      val powerFlowResult = newtonRaphsonPF(
+
+      val (operatingPoint, slackNodeVoltages) = composeOperatingPoint(
+        gridModel.gridComponents.nodes,
+        gridModel.gridComponents.transformers,
+        gridModel.gridComponents.transformers3w,
+        gridModel.nodeUuidToIndexMap,
+        gridAgentBaseData.receivedValueStore,
+        gridModel.mainRefSystem
+      )
+
+      newtonRaphsonPF(
         gridModel,
         gridAgentBaseData.powerFlowParams.maxIterations,
-        composeOperatingPoint(
-          gridModel.gridComponents.nodes,
-          gridModel.gridComponents.transformers,
-          gridModel.gridComponents.transformers3w,
-          gridModel.nodeUuidToIndexMap,
-          gridAgentBaseData.receivedValueStore,
-          gridModel.mainRefSystem
-        )
-      )(gridAgentBaseData.powerFlowParams.epsilon)
-
-      // if res is valid, ask our assets (if any) for updated power values based on the newly determined nodal voltages
-      powerFlowResult match {
+        operatingPoint,
+        slackNodeVoltages
+      )(gridAgentBaseData.powerFlowParams.epsilon) match {
+        // if res is valid, ask our assets (if any) for updated power values based on the newly determined nodal voltages
         case validPowerFlowResult: ValidNewtonRaphsonPFResult =>
           log.debug(
             "{}",
@@ -579,16 +581,21 @@ trait DBFSAlgorithm extends PowerFlowSupport with GridResultsSupport {
           s"$actorName Unable to get results from previous sweep ${gridAgentBaseData.currentSweepNo - 1}!"
         )
       )
-      newtonRaphsonPF(
-        gridModel,
-        gridAgentBaseData.powerFlowParams.maxIterations,
+
+      val (operatingPoint, slackNodeVoltages) =
         composeOperatingPointWithUpdatedSlackVoltages(
           receivedSlackValues,
           previousSweepData.sweepData,
           gridModel.gridComponents.transformers,
           gridModel.gridComponents.transformers3w,
           gridModel.mainRefSystem
         )
+
+      newtonRaphsonPF(
+        gridModel,
+        gridAgentBaseData.powerFlowParams.maxIterations,
+        operatingPoint,
+        slackNodeVoltages
       )(gridAgentBaseData.powerFlowParams.epsilon) match {
         case validPowerFlowResult: ValidNewtonRaphsonPFResult =>
           log.debug(
@@ -658,7 +665,10 @@ trait DBFSAlgorithm extends PowerFlowSupport with GridResultsSupport {
 
         case failedNewtonRaphsonPFResult: FailedNewtonRaphsonPFResult =>
           val powerFlowDoneData =
-            PowerFlowDoneData(gridAgentBaseData, failedNewtonRaphsonPFResult)
+            PowerFlowDoneData(
+              gridAgentBaseData,
+              failedNewtonRaphsonPFResult
+            )
           log.warning(
             "Power flow with updated slack voltage did finally not converge!"
           )
@@ -705,7 +715,7 @@ trait DBFSAlgorithm extends PowerFlowSupport with GridResultsSupport {
 
       /* This is the highest grid agent, therefore no data is received for the slack node. Suppress, that it is looked
        * up in the empty store. */
-      val operationPoint = composeOperatingPoint(
+      val (operationPoint, slackNodeVoltages) = composeOperatingPoint(
         gridModel.gridComponents.nodes,
         gridModel.gridComponents.transformers,
         gridModel.gridComponents.transformers3w,
@@ -729,7 +739,8 @@ trait DBFSAlgorithm extends PowerFlowSupport with GridResultsSupport {
          newtonRaphsonPF(
            gridModel,
            gridAgentBaseData.powerFlowParams.maxIterations,
-           operationPoint
+           operationPoint,
+           slackNodeVoltages
          )(gridAgentBaseData.powerFlowParams.epsilon) match {
            case validPowerFlowResult: ValidNewtonRaphsonPFResult =>
              log.debug(

src/main/scala/edu/ie3/simona/agent/grid/PowerFlowSupport.scala

@@ -6,7 +6,6 @@
 
 package edu.ie3.simona.agent.grid
 
-import java.util.UUID
 import breeze.math.Complex
 import edu.ie3.powerflow.NewtonRaphsonPF
 import edu.ie3.powerflow.model.NodeData.{PresetData, StateData}
@@ -16,22 +15,16 @@ import edu.ie3.powerflow.model.StartData.WithForcedStartVoltages
 import edu.ie3.powerflow.model.enums.NodeType
 import edu.ie3.simona.agent.grid.ReceivedValues.ReceivedSlackValues
 import edu.ie3.simona.exceptions.agent.DBFSAlgorithmException
-import edu.ie3.simona.model.grid.{
-  GridModel,
-  NodeModel,
-  RefSystem,
-  Transformer3wModel,
-  TransformerModel
-}
+import edu.ie3.simona.model.grid._
 import edu.ie3.simona.ontology.messages.PowerMessage.ProvidePowerMessage
 import edu.ie3.simona.ontology.messages.VoltageMessage.ProvideSlackVoltageMessage
 import edu.ie3.util.quantities.PowerSystemUnits
 import tech.units.indriya.ComparableQuantity
+import tech.units.indriya.quantity.Quantities
 
+import java.util.UUID
 import javax.measure.Quantity
 import javax.measure.quantity.{Dimensionless, ElectricPotential}
-import tech.units.indriya.quantity.Quantities
-
 import scala.collection.mutable
 import scala.util.{Failure, Success, Try}
 
@@ -67,7 +60,8 @@ trait PowerFlowSupport {
     *   receivedValuesStore)
     * @return
     *   current operating point of the grid to be used with
-    *   [[edu.ie3.powerflow.NewtonRaphsonPF.calculate()]]
+    *   [[edu.ie3.powerflow.NewtonRaphsonPF.calculate()]] as well as the complex
+    *   slack node target voltages
     */
   protected def composeOperatingPoint(
       nodes: Set[NodeModel],
@@ -78,14 +72,13 @@ trait PowerFlowSupport {
       gridMainRefSystem: RefSystem,
       targetVoltageFromReceivedData: Boolean = true,
       ignoreTargetVoltage: Boolean = false
-  ): Array[PresetData] = {
-
-    val mainRefSystemPowerUnit = gridMainRefSystem.nominalPower.getUnit
-
+  ): (Array[PresetData], WithForcedStartVoltages) = {
     nodes.toArray.map { nodeModel =>
       // note: currently we only support pq nodes as we not distinguish between pq/pv nodes -
       // when slack emulators or pv-node assets are added this needs to be considered here
       val nodeType = if (nodeModel.isSlack) NodeType.SL else NodeType.PQ
+
+      /* Determine the operating point for this given node */
       val nodeIdx = nodeUuidToIndexMap.getOrElse(
         nodeModel.uuid,
         throw new RuntimeException(
@@ -97,6 +90,7 @@ trait PowerFlowSupport {
         receivedValuesStore.nodeToReceivedPower
           .get(nodeModel.uuid) match {
           case Some(actorRefsWithPower) =>
+            val powerUnit = gridMainRefSystem.nominalPower.getUnit
             val (p, q) = actorRefsWithPower
               .map { case (_, powerMsg) => powerMsg }
               .collect {
@@ -113,8 +107,8 @@ trait PowerFlowSupport {
               }
               .foldLeft(
                 (
-                  Quantities.getQuantity(0, mainRefSystemPowerUnit),
-                  Quantities.getQuantity(0, mainRefSystemPowerUnit)
+                  Quantities.getQuantity(0, powerUnit),
+                  Quantities.getQuantity(0, powerUnit)
                 )
               ) { case ((pSum, qSum), powerMessage) =>
                 (pSum.add(powerMessage.p), qSum.add(powerMessage.q))
@@ -127,8 +121,8 @@ trait PowerFlowSupport {
           case None => new Complex(0, 0)
         }
 
-      val targetVoltageInPu =
-        if (targetVoltageFromReceivedData && nodeType == NodeType.SL) {
+      val targetVoltage =
+        if (targetVoltageFromReceivedData && nodeModel.isSlack) {
           /* If the preset voltage is meant to be determined by means of received data and the node is a slack node
            * (only then there is received data), look it up and transform it */
           val receivedSlackVoltage =
@@ -148,17 +142,37 @@ trait PowerFlowSupport {
             gridMainRefSystem
           )
         } else {
-          /* Either, the received data shall not be considered or the node is not a slack node: Depending on if the
-           * node's is meant to be neglected or not, return a dummy value. */
-          Complex.one * (if (!ignoreTargetVoltage)
-                           nodeModel.vTarget
-                             .to(PowerSystemUnits.PU)
-                             .getValue
-                             .doubleValue()
-                         else 1.0)
+          // Either the received data shall not be considered or the node is not a slack node
+          Complex.one *
+            (if (!ignoreTargetVoltage)
+               nodeModel.vTarget
+                 .to(PowerSystemUnits.PU)
+                 .getValue
+                 .doubleValue()
+             else 1.0)
         }
 
-      PresetData(nodeIdx, nodeType, apparentPower, targetVoltageInPu.abs)
+      val optStateData = Option.when(nodeModel.isSlack)(
+        StateData(
+          nodeIdx,
+          NodeType.SL,
+          targetVoltage,
+          apparentPower
+        )
+      )
+
+      (
+        PresetData(
+          nodeIdx,
+          nodeType,
+          apparentPower,
+          targetVoltage.abs
+        ),
+        optStateData
+      )
+    }.unzip match {
+      case (operatingPoint, stateData) =>
+        (operatingPoint, WithForcedStartVoltages(stateData.flatten))
     }
   }
 
@@ -181,51 +195,63 @@ trait PowerFlowSupport {
     *   instance of [[RefSystem]] of the grid under investigation
     * @return
     *   current operating point of the grid to be used with
-    *   [[edu.ie3.powerflow.NewtonRaphsonPF.calculate()]]
+    *   [[edu.ie3.powerflow.NewtonRaphsonPF.calculate()]] as well as the complex
+    *   slack node target voltages
     */
   protected def composeOperatingPointWithUpdatedSlackVoltages(
       receivedSlackValues: ReceivedSlackValues,
       sweepDataValues: Vector[SweepValueStore.SweepValueStoreData],
       transformers2w: Set[TransformerModel],
       transformers3w: Set[Transformer3wModel],
       gridMainRefSystem: RefSystem
-  ): Array[PresetData] = {
-    sweepDataValues
-      .map(sweepValueStoreData => {
-
-        val nodeStateData = sweepValueStoreData.stateData
-        val targetVoltage = if (nodeStateData.nodeType == NodeType.SL) {
-          val receivedSlackVoltage = receivedSlackValues.values
-            .map { case (_, slackVoltageMsg) => slackVoltageMsg }
-            .find(_.nodeUuid == sweepValueStoreData.nodeUuid)
-            .getOrElse(
-              throw new RuntimeException(
-                s"Unable to find node with uuid " +
-                  s"${sweepValueStoreData.nodeUuid} in received slack voltage values!"
-              )
+  ): (Array[PresetData], WithForcedStartVoltages) =
+    sweepDataValues.map { sweepValueStoreData =>
+      val nodeStateData = sweepValueStoreData.stateData
+      val targetVoltage = if (nodeStateData.nodeType == NodeType.SL) {
+        val receivedSlackVoltage = receivedSlackValues.values
+          .map { case (_, slackVoltageMsg) => slackVoltageMsg }
+          .find(_.nodeUuid == sweepValueStoreData.nodeUuid)
+          .getOrElse(
+            throw new RuntimeException(
+              s"Unable to find node with uuid " +
+                s"${sweepValueStoreData.nodeUuid} in received slack voltage values!"
             )
-
-          transformVoltage(
-            receivedSlackVoltage,
-            sweepValueStoreData.nodeUuid,
-            transformers2w,
-            transformers3w,
-            gridMainRefSystem
           )
-        } else
-          new Complex(1, 0)
 
-        // note: target voltage will be ignored for slack node if provided
+        transformVoltage(
+          receivedSlackVoltage,
+          sweepValueStoreData.nodeUuid,
+          transformers2w,
+          transformers3w,
+          gridMainRefSystem
+        )
+      } else
+        Complex.one
+
+      // note: target voltage will be ignored for slack node if provided
+      (
         PresetData(
           nodeStateData.index,
           nodeStateData.nodeType,
           nodeStateData.power,
           targetVoltage.abs
+        ),
+        Option.when(nodeStateData.nodeType == NodeType.SL)(
+          StateData(
+            nodeStateData.index,
+            nodeStateData.nodeType,
+            targetVoltage,
+            nodeStateData.power
+          )
         )
-      })
-      .toArray
-
-  }
+      )
+    }.unzip match {
+      case (operatingPoint, stateData) =>
+        (
+          operatingPoint.toArray,
+          WithForcedStartVoltages(stateData.flatten.toArray)
+        )
+    }
 
   /** A debug method that composes a string with voltage information (in p.u.)
     * from a [[ValidNewtonRaphsonPFResult]]
@@ -409,6 +435,8 @@ trait PowerFlowSupport {
     *   Maximum permissible iterations
     * @param operatingPoint
     *   Current operation point of the grid
+    * @param slackVoltages
+    *   Complex target voltages of the slack nodes
     * @param epsilons
     *   Ascending ordered list of convergence thresholds for relaxation
     * @return
@@ -417,7 +445,8 @@ trait PowerFlowSupport {
   protected final def newtonRaphsonPF(
       gridModel: GridModel,
       maxIterations: Int,
-      operatingPoint: Array[PresetData]
+      operatingPoint: Array[PresetData],
+      slackVoltages: WithForcedStartVoltages
   )(epsilons: Vector[Double]): PowerFlowResult = {
     epsilons.headOption match {
       case Some(epsilon) =>
@@ -426,26 +455,15 @@ trait PowerFlowSupport {
             gridModel.nodeUuidToIndexMap,
             gridModel.gridComponents
           )
-        // / add WithForcedVoltageVector for slackNode
-        // / as we know that there is at least one slack in every grid available here, we can do a direct call on element zero
-        // // NOTE: currently only the first slack node is taken, needs to be adapted when several slacks are present
-        val forcedSlackNodeVoltage =
-          operatingPoint.find(_.nodeType == NodeType.SL) match {
-            case Some(slackNodeData) =>
-              WithForcedStartVoltages(Array(StateData(slackNodeData)))
-            case None =>
-              throw new DBFSAlgorithmException(
-                s"Unable to find a slack node in grid ${gridModel.subnetNo}."
-              )
-          }
+
         // / execute
         val powerFlow =
           NewtonRaphsonPF(epsilon, maxIterations, admittanceMatrix)
 
         Try {
           powerFlow.calculate(
             operatingPoint,
-            Some(forcedSlackNodeVoltage)
+            Some(slackVoltages)
           )
         }.map {
           case _: PowerFlowResult.FailedPowerFlowResult if epsilons.size > 1 =>
@@ -456,7 +474,12 @@ trait PowerFlowSupport {
               epsilon,
               epsilonsLeft.headOption.getOrElse("")
             )
-            newtonRaphsonPF(gridModel, maxIterations, operatingPoint)(
+            newtonRaphsonPF(
+              gridModel,
+              maxIterations,
+              operatingPoint,
+              slackVoltages
+            )(
               epsilonsLeft
             )
           case result =>