Separate internal and external residuals in PowerElectronics models

CHANGELOG.md

@@ -54,6 +54,7 @@
 - Added SEXS-PTI Exciter Model
 - Added 200 Bus Synthetic Illinois Case
 - Added node objects to `PowerElectronics` module & updated all examples to make use of them.
+- Separated internal and external residuals of `PowerElectronics` models.
 
 ## v0.1
 

GridKit/Model/PowerElectronics/Capacitor/Capacitor.cpp

@@ -57,15 +57,19 @@ namespace GridKit
    *
    */
   template <class ScalarT, typename IdxT>
-  int Capacitor<ScalarT, IdxT>::evaluateResidual()
+  int Capacitor<ScalarT, IdxT>::evaluateInternalResidual()
+  {
+    f_[2] = -C_ * yp_[2] + y_[0] - y_[1] - y_[2];
+    return 0;
+  }
+
+  template <class ScalarT, typename IdxT>
+  int Capacitor<ScalarT, IdxT>::evaluateExternalResidual()
   {
     // input
     f_[0] = C_ * yp_[2];
     // output
     f_[1] = -C_ * yp_[2];
-
-    // internal
-    f_[2] = -C_ * yp_[2] + y_[0] - y_[1] - y_[2];
     return 0;
   }
 

GridKit/Model/PowerElectronics/Capacitor/Capacitor.hpp

@@ -47,7 +47,8 @@ namespace GridKit
 
     int initialize();
     int tagDifferentiable();
-    int evaluateResidual();
+    int evaluateInternalResidual() final;
+    int evaluateExternalResidual() final;
     int evaluateJacobian();
     int evaluateIntegrand();
 

GridKit/Model/PowerElectronics/CircuitComponent.hpp

@@ -24,14 +24,6 @@ namespace GridKit
 
     CircuitComponent() = default;
 
-    ~CircuitComponent()
-    {
-      if (connection_nodes_ != nullptr)
-      {
-        delete[] connection_nodes_;
-      }
-    };
-
     /**
      * @note Cannot be marked final, since it is overriden to recurse in the system model.
      */
@@ -70,12 +62,7 @@ namespace GridKit
      */
     int setExternalConnectionNodes(IdxT local_index, IdxT global_index)
     {
-      if (connection_nodes_ == nullptr)
-      {
-        connection_nodes_ = new IdxT[static_cast<size_t>(size_)];
-      }
-
-      connection_nodes_[local_index] = global_index;
+      connection_nodes_[static_cast<size_t>(local_index)] = global_index;
       return 0;
     }
 
@@ -108,6 +95,8 @@ namespace GridKit
       jacobian_coo_cols_   = std::make_unique<IdxT[]>(static_cast<size_t>(nnz_));
       jacobian_coo_values_ = std::make_unique<RealT[]>(static_cast<size_t>(nnz_));
 
+      connection_nodes_ = std::make_unique<IdxT[]>(static_cast<size_t>(size_));
+
       y_.resize(static_cast<size_t>(size_));
       yp_.resize(static_cast<size_t>(size_));
       f_.resize(static_cast<size_t>(size_));
@@ -145,6 +134,38 @@ namespace GridKit
       return jacobian_coo_values_.get();
     }
 
+    /**
+     * @brief Evaluating the residual of a CircuitComponent should be done by evaluating the
+     * internal residuals and external residuals. CircuitComponents should overload those
+     * functions for their residuals (and the system will call those function instead of this one),
+     * so there is no reason to overload this functionality.
+     *
+     * @return An error code, or 0 is successful.
+     */
+    int evaluateResidual() final
+    {
+      if (int err_code = evaluateInternalResidual())
+        return err_code;
+
+      return evaluateExternalResidual();
+    }
+
+    /**
+     * @brief Evaluate all of the residuals of internal variables of the component,
+     * modifying \ref f_.
+     *
+     * @return An error code, or 0 if successful.
+     */
+    virtual int evaluateInternalResidual() = 0;
+
+    /**
+     * @brief Evaluate all of the residuals of external variables of the component,
+     * modifying \ref f_
+     *
+     * @return An error code, or 0 if successful.
+     */
+    virtual int evaluateExternalResidual() = 0;
+
   protected:
     /**
      * @brief Reset the Jacobian so it can be constructed. Helper method for \ref setJacValues().
@@ -371,11 +392,11 @@ namespace GridKit
     }
 
   protected:
-    size_t         n_extern_;
-    size_t         n_intern_;
-    std::set<IdxT> extern_indices_;
+    size_t                  n_extern_;
+    size_t                  n_intern_;
+    std::set<IdxT>          extern_indices_;
     ///@todo may want to replace the mapping of connection_nodes to Node objects instead of IdxT. Allows for container free setup
-    IdxT*          connection_nodes_ = nullptr;
+    std::unique_ptr<IdxT[]> connection_nodes_;
 
   protected:
     IdxT size_{0};

GridKit/Model/PowerElectronics/DistributedGenerator/DistributedGenerator.cpp

@@ -89,26 +89,9 @@ namespace GridKit
    *
    */
   template <class ScalarT, typename IdxT>
-  int DistributedGenerator<ScalarT, IdxT>::evaluateResidual()
+  int DistributedGenerator<ScalarT, IdxT>::evaluateInternalResidual()
   {
-    //   ### Externals Componenets ###
-
     ScalarT omega = wb_ - mp_ * y_[4];
-    // ref common ref motor angle
-    /// @todo fix boolian conditional, unclear result
-    if (refframe_)
-    {
-      f_[0] = omega - y_[0];
-    }
-    else
-    {
-      f_[0] = 0.0;
-    }
-
-    // output
-    // current transformed to common frame
-    f_[1] = std::cos(y_[3]) * y_[14] - std::sin(y_[3]) * y_[15];
-    f_[2] = std::sin(y_[3]) * y_[14] + std::cos(y_[3]) * y_[15];
 
     // Take incoming voltages to current rotator reference frame
     ScalarT vbd_in = std::cos(y_[3]) * y_[1] + std::sin(y_[3]) * y_[2];
@@ -152,6 +135,27 @@ namespace GridKit
     return 0;
   }
 
+  template <class ScalarT, typename IdxT>
+  int DistributedGenerator<ScalarT, IdxT>::evaluateExternalResidual()
+  {
+    ScalarT omega = wb_ - mp_ * y_[4];
+    // ref common ref motor angle
+    if (refframe_)
+    {
+      f_[0] = omega - y_[0];
+    }
+    else
+    {
+      f_[0] = 0.0;
+    }
+
+    // output
+    // current transformed to common frame
+    f_[1] = std::cos(y_[3]) * y_[14] - std::sin(y_[3]) * y_[15];
+    f_[2] = std::sin(y_[3]) * y_[14] + std::cos(y_[3]) * y_[15];
+    return 0;
+  }
+
   /**
    * @brief Compute the jacobian of the DistributedGenerator for iteration. dF/dy - \alpha dF/dy'
    *

GridKit/Model/PowerElectronics/DistributedGenerator/DistributedGenerator.hpp

@@ -74,7 +74,8 @@ namespace GridKit
     int initialize();
     int allocate() final;
     int tagDifferentiable();
-    int evaluateResidual();
+    int evaluateInternalResidual() final;
+    int evaluateExternalResidual() final;
     int evaluateJacobian();
     int evaluateIntegrand();
     int initializeAdjoint();

GridKit/Model/PowerElectronics/InductionMotor/InductionMotor.cpp

@@ -68,19 +68,24 @@ namespace GridKit
    *
    */
   template <class ScalarT, typename IdxT>
-  int InductionMotor<ScalarT, IdxT>::evaluateResidual()
+  int InductionMotor<ScalarT, IdxT>::evaluateInternalResidual()
   {
+    f_[5] = (1.0 / 3.0) * (2.0 * y_[0] - y_[1] - y_[2]) - Rs_ * y_[5] - (Lls_ + Lms_) * yp_[5] - Lms_ * yp_[6];
+    f_[6] = (1.0 / std::sqrt(3.0)) * (-y_[1] + y_[2]) - Rs_ * y_[6] - (Lls_ + Lms_) * yp_[6] - Lms_ * yp_[5];
+    f_[7] = (y_[0] + y_[1] + y_[2]) / 3.0 - Rs_ * y_[7] - Lls_ * yp_[7];
+    f_[8] = Rr_ * y_[8] + (Llr_ + Lms_) * yp_[8] + Lms_ * yp_[5] - (P_ / 2.0) * y_[3] * ((Llr_ + Lms_) * y_[9] + Lms_ * y_[6]);
+    f_[9] = Rr_ * y_[9] + (Llr_ + Lms_) * yp_[9] + Lms_ * yp_[6] + (P_ / 2.0) * y_[3] * ((Llr_ + Lms_) * y_[8] + Lms_ * y_[5]);
+    return 0;
+  }
 
+  template <class ScalarT, typename IdxT>
+  int InductionMotor<ScalarT, IdxT>::evaluateExternalResidual()
+  {
     f_[0] = y_[5] + y_[7];
     f_[1] = (-1.0 / 2.0) * y_[5] - (std::sqrt(3.0) / 2.0) * y_[6] + y_[7];
     f_[2] = (-1.0 / 2.0) * y_[5] + (std::sqrt(3.0) / 2.0) * y_[6] + y_[7];
     f_[3] = RJ_ * yp_[3] - (3.0 / 4.0) * P_ * Lms_ * (y_[5] * y_[9] - y_[6] * y_[8]);
     f_[4] = yp_[4] - y_[3];
-    f_[5] = (1.0 / 3.0) * (2.0 * y_[0] - y_[1] - y_[2]) - Rs_ * y_[5] - (Lls_ + Lms_) * yp_[5] - Lms_ * yp_[6];
-    f_[6] = (1.0 / std::sqrt(3.0)) * (-y_[1] + y_[2]) - Rs_ * y_[6] - (Lls_ + Lms_) * yp_[6] - Lms_ * yp_[5];
-    f_[7] = (y_[0] + y_[1] + y_[2]) / 3.0 - Rs_ * y_[7] - Lls_ * yp_[7];
-    f_[8] = Rr_ * y_[8] + (Llr_ + Lms_) * yp_[8] + Lms_ * yp_[5] - (P_ / 2.0) * y_[3] * ((Llr_ + Lms_) * y_[9] + Lms_ * y_[6]);
-    f_[9] = Rr_ * y_[9] + (Llr_ + Lms_) * yp_[9] + Lms_ * yp_[6] + (P_ / 2.0) * y_[3] * ((Llr_ + Lms_) * y_[8] + Lms_ * y_[5]);
     return 0;
   }
 

GridKit/Model/PowerElectronics/InductionMotor/InductionMotor.hpp

@@ -47,7 +47,8 @@ namespace GridKit
 
     int initialize();
     int tagDifferentiable();
-    int evaluateResidual();
+    int evaluateInternalResidual() final;
+    int evaluateExternalResidual() final;
     int evaluateJacobian();
     int evaluateIntegrand();
 

GridKit/Model/PowerElectronics/Inductor/Inductor.cpp

@@ -56,14 +56,19 @@ namespace GridKit
    *
    */
   template <class ScalarT, typename IdxT>
-  int Inductor<ScalarT, IdxT>::evaluateResidual()
+  int Inductor<ScalarT, IdxT>::evaluateInternalResidual()
+  {
+    f_[2] = -L_ * yp_[2] + y_[1] - y_[0];
+    return 0;
+  }
+
+  template <class ScalarT, typename IdxT>
+  int Inductor<ScalarT, IdxT>::evaluateExternalResidual()
   {
     // input
     f_[0] = -y_[2];
     // output
     f_[1] = y_[2];
-    // internal
-    f_[2] = -L_ * yp_[2] + y_[1] - y_[0];
     return 0;
   }
 

GridKit/Model/PowerElectronics/Inductor/Inductor.hpp

@@ -50,7 +50,8 @@ namespace GridKit
     int initialize();
     int allocate() final;
     int tagDifferentiable();
-    int evaluateResidual();
+    int evaluateInternalResidual() final;
+    int evaluateExternalResidual() final;
     int evaluateJacobian();
     int evaluateIntegrand();
 

GridKit/Model/PowerElectronics/LinearTransformer/LinearTransformer.cpp

@@ -67,15 +67,20 @@ namespace GridKit
 
   /**
    * @brief Computes the component resisdual
-   *
    */
   template <class ScalarT, typename IdxT>
-  int LinearTransformer<ScalarT, IdxT>::evaluateResidual()
+  int LinearTransformer<ScalarT, IdxT>::evaluateInternalResidual()
+  {
+    f_[2] = y_[0] - R0_ * y_[2] - L0_ * yp_[2] - M_ * yp_[3];
+    f_[3] = y_[1] - R1_ * y_[3] - M_ * yp_[2] - L1_ * yp_[3];
+    return 0;
+  }
+
+  template <class ScalarT, typename IdxT>
+  int LinearTransformer<ScalarT, IdxT>::evaluateExternalResidual()
   {
     f_[0] = y_[2];
     f_[1] = y_[3];
-    f_[2] = y_[0] - R0_ * y_[2] - L0_ * yp_[2] - M_ * yp_[3];
-    f_[2] = y_[1] - R1_ * y_[3] - M_ * yp_[2] - L1_ * yp_[3];
     return 0;
   }
 

GridKit/Model/PowerElectronics/LinearTransformer/LinearTransformer.hpp

@@ -47,7 +47,8 @@ namespace GridKit
 
     int initialize();
     int tagDifferentiable();
-    int evaluateResidual();
+    int evaluateInternalResidual() final;
+    int evaluateExternalResidual() final;
     int evaluateJacobian();
     int evaluateIntegrand();
 

GridKit/Model/PowerElectronics/MicrogridBusDQ/MicrogridBusDQ.cpp

@@ -56,16 +56,22 @@ namespace GridKit
     return 0;
   }
 
+  template <class ScalarT, typename IdxT>
+  int MicrogridBusDQ<ScalarT, IdxT>::evaluateInternalResidual()
+  {
+    return 0;
+  }
+
   /**
-   * @brief Evaluate residual of microgrid line
+   * @brief Evaluate residual
    * This model has "Virtual resistors". The voltage of the bus divided by its virtual resistance.
    * The components are external to allow for outside components to add inductances to the terms.
    *
    * refernce to equations in class header
    *
    */
   template <class ScalarT, typename IdxT>
-  int MicrogridBusDQ<ScalarT, IdxT>::evaluateResidual()
+  int MicrogridBusDQ<ScalarT, IdxT>::evaluateExternalResidual()
   {
     // bus voltage
     f_[0] = -y_[0] / RN_;

GridKit/Model/PowerElectronics/MicrogridBusDQ/MicrogridBusDQ.hpp

@@ -50,7 +50,8 @@ namespace GridKit
     int initialize();
     int allocate() final;
     int tagDifferentiable();
-    int evaluateResidual();
+    int evaluateInternalResidual() final;
+    int evaluateExternalResidual() final;
     int evaluateJacobian();
     int evaluateIntegrand();
 

GridKit/Model/PowerElectronics/MicrogridLine/MicrogridLine.cpp

@@ -70,23 +70,28 @@ namespace GridKit
    *
    */
   template <class ScalarT, typename IdxT>
-  int MicrogridLine<ScalarT, IdxT>::evaluateResidual()
+  int MicrogridLine<ScalarT, IdxT>::evaluateInternalResidual()
+  {
+    f_[5] = -yp_[5] - (R_ / L_) * y_[5] + y_[0] * y_[6] + (y_[1] - y_[3]) / L_;
+    f_[6] = -yp_[6] - (R_ / L_) * y_[6] - y_[0] * y_[5] + (y_[2] - y_[4]) / L_;
+
+    return 0;
+  }
+
+  template <class ScalarT, typename IdxT>
+  int MicrogridLine<ScalarT, IdxT>::evaluateExternalResidual()
   {
     // ref motor
     f_[0] = 0.0;
 
-    // input
+    // Port 1
     f_[1] = -y_[5];
     f_[2] = -y_[6];
 
-    // output
+    // Port 2
     f_[3] = y_[5];
     f_[4] = y_[6];
 
-    // Internal variables
-    f_[5] = -yp_[5] - (R_ / L_) * y_[5] + y_[0] * y_[6] + (y_[1] - y_[3]) / L_;
-    f_[6] = -yp_[6] - (R_ / L_) * y_[6] - y_[0] * y_[5] + (y_[2] - y_[4]) / L_;
-
     return 0;
   }
 

GridKit/Model/PowerElectronics/MicrogridLine/MicrogridLine.hpp

@@ -50,7 +50,8 @@ namespace GridKit
     int initialize();
     int allocate() final;
     int tagDifferentiable();
-    int evaluateResidual();
+    int evaluateInternalResidual() final;
+    int evaluateExternalResidual() final;
     int evaluateJacobian();
     int evaluateIntegrand();