From 4bff26fd2b3a662ac120776de20c2871d6241f02 Mon Sep 17 00:00:00 2001 From: lukelowry Date: Mon, 24 Nov 2025 14:15:50 -0600 Subject: [PATCH 01/12] Change to logistic function --- GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp | 2 +- GridKit/Model/PhasorDynamics/Exciter/IEEET1/README.md | 5 ----- GridKit/Model/PhasorDynamics/Governor/Tgov1/README.md | 2 +- GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1Impl.hpp | 2 +- 4 files changed, 3 insertions(+), 8 deletions(-) diff --git a/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp b/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp index 1623c17c5..e3b738224 100644 --- a/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp +++ b/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp @@ -253,7 +253,7 @@ namespace GridKit template ScalarT Ieeet1::sigmoid(ScalarT x) { - return ((HALF * mu_ * x) / (ONE + std::abs(mu_ * x))) + HALF; + return ONE / (ONE + std::exp(- mu_ * x)); } /** diff --git a/GridKit/Model/PhasorDynamics/Exciter/IEEET1/README.md b/GridKit/Model/PhasorDynamics/Exciter/IEEET1/README.md index 5f4422f30..666368195 100644 --- a/GridKit/Model/PhasorDynamics/Exciter/IEEET1/README.md +++ b/GridKit/Model/PhasorDynamics/Exciter/IEEET1/README.md @@ -160,11 +160,6 @@ The scale of the sigmoid function ($\alpha$ on the order of $10^3$) should be ch \begin{aligned} \sigma(x) = \dfrac{1}{1+\exp(-\alpha x)} - \approx - \dfrac{1}{2} - \left( - \dfrac{\alpha x}{1+|\alpha x|} + 1 - \right) \end{aligned} ``` diff --git a/GridKit/Model/PhasorDynamics/Governor/Tgov1/README.md b/GridKit/Model/PhasorDynamics/Governor/Tgov1/README.md index 27431fbf0..f27390e25 100644 --- a/GridKit/Model/PhasorDynamics/Governor/Tgov1/README.md +++ b/GridKit/Model/PhasorDynamics/Governor/Tgov1/README.md @@ -104,7 +104,7 @@ The indicator function $\phi$ can be defined in terms of a scaled activation fun The scale of the sigmoid function ($\alpha$ on the order of $10^3$) should be chosen so that for all practical parameters of the TGOV1 model, the sigmoid acts as a step function. This is further approximated by an algebraic form to obtain a practical function during implementation. ```math \begin{aligned} - \sigma(x) =\dfrac{1}{1+\exp(-\alpha x)}\approx \dfrac{1}{2}\left(\dfrac{\alpha x}{1+|\alpha x|}\right) + \sigma(x) =\dfrac{1}{1+\exp(-\alpha x)} \end{aligned} ``` diff --git a/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1Impl.hpp b/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1Impl.hpp index 50c6709ff..323347749 100644 --- a/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1Impl.hpp +++ b/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1Impl.hpp @@ -241,7 +241,7 @@ namespace GridKit template ScalarT Tgov1::sigmoid(ScalarT x) { - return ((HALF * mu_ * x) / (ONE + std::abs(mu_ * x))) + HALF; + return ONE / (ONE + std::exp(- mu_ * x)); } /** From 5dff8e245652d847d10902e1d4a5a6db07dbdcfa Mon Sep 17 00:00:00 2001 From: lukelowry Date: Mon, 24 Nov 2025 20:20:35 +0000 Subject: [PATCH 02/12] Apply pre-commmit fixes --- GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp | 2 +- GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1Impl.hpp | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp b/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp index e3b738224..3d5ce0d13 100644 --- a/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp +++ b/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp @@ -253,7 +253,7 @@ namespace GridKit template ScalarT Ieeet1::sigmoid(ScalarT x) { - return ONE / (ONE + std::exp(- mu_ * x)); + return ONE / (ONE + std::exp(-mu_ * x)); } /** diff --git a/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1Impl.hpp b/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1Impl.hpp index 323347749..ae4d771f9 100644 --- a/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1Impl.hpp +++ b/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1Impl.hpp @@ -241,7 +241,7 @@ namespace GridKit template ScalarT Tgov1::sigmoid(ScalarT x) { - return ONE / (ONE + std::exp(- mu_ * x)); + return ONE / (ONE + std::exp(-mu_ * x)); } /** From f82881c3e226cadb1de222421419730ea8719eaa Mon Sep 17 00:00:00 2001 From: lukelowry Date: Mon, 24 Nov 2025 16:34:26 -0600 Subject: [PATCH 03/12] modify saturation model exciter --- .../Exciter/IEEET1/Ieeet1Impl.hpp | 7 ++-- .../PhasorDynamics/Exciter/IEEET1/README.md | 32 +++---------------- .../PhasorDynamics/Governor/Tgov1/README.md | 4 +-- 3 files changed, 9 insertions(+), 34 deletions(-) diff --git a/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp b/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp index 3d5ce0d13..cfaa60091 100644 --- a/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp +++ b/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp @@ -340,11 +340,8 @@ namespace GridKit f_[6] = -ve + ksat * efdp; f_[7] = -efd + efdp + omega * efdp * Ispdlim_; - // TODO smooth approaximation for Enzyme - // NOTE seems about double PW saturation. - f_[8] = -ksat; - if (efdp > SA_) - f_[8] += SB_ * (efdp - SA_) * (efdp - SA_); + ScalarT efd_sat = (efdp - SA_) * ( this->sigmoid(efdp-SA_) ); + f_[8] = -ksat + SB_ * efd_sat * efd_sat; return 0; } diff --git a/GridKit/Model/PhasorDynamics/Exciter/IEEET1/README.md b/GridKit/Model/PhasorDynamics/Exciter/IEEET1/README.md index 666368195..a5f47a62f 100644 --- a/GridKit/Model/PhasorDynamics/Exciter/IEEET1/README.md +++ b/GridKit/Model/PhasorDynamics/Exciter/IEEET1/README.md @@ -165,25 +165,6 @@ The scale of the sigmoid function ($\alpha$ on the order of $10^3$) should be ch Applying these approximations produces a smooth approximation of the explicit differential equation of $V_R$. The approximation written out below approaches an exact solution as $\alpha\to\infty$. - ```math \begin{aligned} \dot{V}_R @@ -223,23 +204,20 @@ The algebraic equations of the exciter. \end{cases} \\ \end{aligned} ``` - #### Smooth Piecewise Approximation (Algebraic) For the algebraic piecewise functions (non-flags), this implementation is straightforward when the approximation above is used. ```math \begin{aligned} + E_{fd} + &= E_{fd}' + \omega E_{fd}' I_{spdlm} \\ k_{sat} - &=\sigma (E_{fd}' -S_A) \cdot S_B(E_{fd}' -S_A)^2 -\end{aligned} -``` -The approximation written out below approaches an exact solution as $\alpha\to\infty$. -```math -\begin{aligned} - k_{sat}&\approx\dfrac{\alpha S_B}{2} \dfrac{(E_{fd}' -S_A)^3}{1+\alpha|E_{fd}' -S_A|} + &=S_B\left[\sigma (E_{fd}' -S_A) \cdot (E_{fd}' -S_A)\right]^2 \end{aligned} ``` +The approximation approaches an exact solution as $\alpha\to\infty$. + ## Initialization At steady state we assume that $V_R$ is at or within its limits, diff --git a/GridKit/Model/PhasorDynamics/Governor/Tgov1/README.md b/GridKit/Model/PhasorDynamics/Governor/Tgov1/README.md index f27390e25..25fc18f84 100644 --- a/GridKit/Model/PhasorDynamics/Governor/Tgov1/README.md +++ b/GridKit/Model/PhasorDynamics/Governor/Tgov1/README.md @@ -93,13 +93,13 @@ smooth approximation (smooth indicator $\phi$) expressed generically as follows. ``` The indicator function $\phi$ can be defined in terms of a scaled activation function ($\sigma$, sigmoid) and the $P_v$ limits $(P_{vmin}, P_{vmax})$. -```math +$$ \begin{aligned} \phi_L(P_v,f)&= \sigma(P_{vmin}-P_v)\sigma(-f) \\ \phi_U(P_v,f)&= \sigma(P_v-P_{vmax})\sigma(f) \\ \phi(P_v,f) &= \left[1-\phi_L\right]\left[1-\phi_U\right]\\ \end{aligned} -``` +$$ The scale of the sigmoid function ($\alpha$ on the order of $10^3$) should be chosen so that for all practical parameters of the TGOV1 model, the sigmoid acts as a step function. This is further approximated by an algebraic form to obtain a practical function during implementation. ```math From 5329702c48b097e23065ac4226fe3e5865663135 Mon Sep 17 00:00:00 2001 From: lukelowry Date: Mon, 24 Nov 2025 22:34:53 +0000 Subject: [PATCH 04/12] Apply pre-commmit fixes --- GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp b/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp index cfaa60091..ff6613003 100644 --- a/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp +++ b/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp @@ -340,8 +340,8 @@ namespace GridKit f_[6] = -ve + ksat * efdp; f_[7] = -efd + efdp + omega * efdp * Ispdlim_; - ScalarT efd_sat = (efdp - SA_) * ( this->sigmoid(efdp-SA_) ); - f_[8] = -ksat + SB_ * efd_sat * efd_sat; + ScalarT efd_sat = (efdp - SA_) * (this->sigmoid(efdp - SA_)); + f_[8] = -ksat + SB_ * efd_sat * efd_sat; return 0; } From a20dc9972f69d0824149248cd3f17801ce2b84fa Mon Sep 17 00:00:00 2001 From: lukelowry Date: Mon, 24 Nov 2025 16:49:29 -0600 Subject: [PATCH 05/12] Update documentation --- .../PhasorDynamics/Exciter/IEEET1/README.md | 29 ++----------------- .../PhasorDynamics/Governor/Tgov1/README.md | 6 ++-- 2 files changed, 6 insertions(+), 29 deletions(-) diff --git a/GridKit/Model/PhasorDynamics/Exciter/IEEET1/README.md b/GridKit/Model/PhasorDynamics/Exciter/IEEET1/README.md index a5f47a62f..64aa85e3a 100644 --- a/GridKit/Model/PhasorDynamics/Exciter/IEEET1/README.md +++ b/GridKit/Model/PhasorDynamics/Exciter/IEEET1/README.md @@ -155,7 +155,7 @@ The indicator function $\phi$ can be defined in terms of a scaled activation fun \end{aligned} ``` -The scale of the sigmoid function ($\alpha$ on the order of $10^3$) should be chosen so that for all practical parameters of the IEEET1 model, the sigmoid acts as a step function. This is further approximated by an algebraic form to obtain a practical function during implementation. +The scale of the sigmoid function ($\alpha$ on the order of $10^3$) should be chosen so that for all practical parameters of the IEEET1 model, the sigmoid acts as a step function. ```math \begin{aligned} \sigma(x) = @@ -163,29 +163,6 @@ The scale of the sigmoid function ($\alpha$ on the order of $10^3$) should be ch \end{aligned} ``` -Applying these approximations produces a smooth approximation of the explicit differential equation of $V_R$. The approximation written out below approaches an exact solution as $\alpha\to\infty$. - -```math -\begin{aligned} - \dot{V}_R - &\approx f - \left[ - 1 + - \dfrac{\alpha^2}{4} - \dfrac{V_{rmin}\text{-}V_R} - {1+\alpha |V_{rmin}\text{-}V_R|} - \dfrac{f}{1+\alpha |f|} - \right] - \left[ - 1- - \dfrac{\alpha^2}{4} - \dfrac{V_R\text{-}V_{rmax}} - {1+\alpha| V_R\text{-}V_{rmax}|} - \dfrac{f}{1+\alpha |f|} - \right]\\ -\end{aligned} -``` - ### Algebraic Equations The algebraic equations of the exciter. @@ -210,9 +187,9 @@ For the algebraic piecewise functions (non-flags), this implementation is straig ```math \begin{aligned} E_{fd} - &= E_{fd}' + \omega E_{fd}' I_{spdlm} \\ + &=(1 + \omega I_{spdlm})E_{fd}' \\ k_{sat} - &=S_B\left[\sigma (E_{fd}' -S_A) \cdot (E_{fd}' -S_A)\right]^2 + &=S_B\left[(E_{fd}' -S_A) \cdot \sigma (E_{fd}' -S_A)\right]^2 \end{aligned} ``` diff --git a/GridKit/Model/PhasorDynamics/Governor/Tgov1/README.md b/GridKit/Model/PhasorDynamics/Governor/Tgov1/README.md index 25fc18f84..b4d695a64 100644 --- a/GridKit/Model/PhasorDynamics/Governor/Tgov1/README.md +++ b/GridKit/Model/PhasorDynamics/Governor/Tgov1/README.md @@ -93,15 +93,15 @@ smooth approximation (smooth indicator $\phi$) expressed generically as follows. ``` The indicator function $\phi$ can be defined in terms of a scaled activation function ($\sigma$, sigmoid) and the $P_v$ limits $(P_{vmin}, P_{vmax})$. -$$ +```math \begin{aligned} \phi_L(P_v,f)&= \sigma(P_{vmin}-P_v)\sigma(-f) \\ \phi_U(P_v,f)&= \sigma(P_v-P_{vmax})\sigma(f) \\ \phi(P_v,f) &= \left[1-\phi_L\right]\left[1-\phi_U\right]\\ \end{aligned} -$$ +``` -The scale of the sigmoid function ($\alpha$ on the order of $10^3$) should be chosen so that for all practical parameters of the TGOV1 model, the sigmoid acts as a step function. This is further approximated by an algebraic form to obtain a practical function during implementation. +The scale of the sigmoid function ($\alpha$ on the order of $10^3$) should be chosen so that for all practical parameters of the TGOV1 model, the sigmoid acts as a step function. ```math \begin{aligned} \sigma(x) =\dfrac{1}{1+\exp(-\alpha x)} From 5ef7c07c614d472bb5b961f3846d63865f298e72 Mon Sep 17 00:00:00 2001 From: Nicholson Koukpaizan Date: Thu, 11 Dec 2025 15:06:56 -0500 Subject: [PATCH 06/12] Setting mu to 240 for exponential-based logistic function. --- GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1.hpp | 4 ++-- GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1.hpp | 3 ++- 2 files changed, 4 insertions(+), 3 deletions(-) diff --git a/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1.hpp b/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1.hpp index 33fa4d98d..03ec0fd9c 100644 --- a/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1.hpp +++ b/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1.hpp @@ -152,8 +152,8 @@ namespace GridKit // Scale of Sigmoid function // (temporary local implementation) - // This value gave higher precision. - const RealT mu_{400000.0}; + // mu set to 4*60, to provide accuracy and finite derivatives. + const RealT mu_{240.0}; // Activation function (sigmoid approximation) ScalarT sigmoid(ScalarT x); diff --git a/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1.hpp b/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1.hpp index c4412d5fb..4597123d3 100644 --- a/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1.hpp +++ b/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1.hpp @@ -122,7 +122,8 @@ namespace GridKit ScalarT pref_{0}; // Scale of Sigmoid function (temporary local implementation) - const ScalarT mu_{4000.0}; + // mu set to 4*60, to provide accuracy and finite derivatives. + const ScalarT mu_{240.0}; /// Component signal extension ComponentSignals signals_; From 68f213d5ab525d9718c626f280f3d29abcd20c63 Mon Sep 17 00:00:00 2001 From: Nicholson Koukpaizan Date: Fri, 12 Dec 2025 10:51:52 -0500 Subject: [PATCH 07/12] Moved sigmoid function to CommonMath.hpp. --- GridKit/CMakeLists.txt | 1 + GridKit/CommonMath.hpp | 33 +++++++++++++++++++ GridKit/Model/Evaluator.hpp | 1 + .../PhasorDynamics/Exciter/IEEET1/Ieeet1.hpp | 8 +---- .../Exciter/IEEET1/Ieeet1Impl.hpp | 26 ++------------- .../PhasorDynamics/Governor/Tgov1/Tgov1.hpp | 7 ---- .../Governor/Tgov1/Tgov1Impl.hpp | 19 ++--------- 7 files changed, 41 insertions(+), 54 deletions(-) create mode 100644 GridKit/CommonMath.hpp diff --git a/GridKit/CMakeLists.txt b/GridKit/CMakeLists.txt index 2a57de7e7..a01f8a6c4 100644 --- a/GridKit/CMakeLists.txt +++ b/GridKit/CMakeLists.txt @@ -28,6 +28,7 @@ add_subdirectory(Solver) install( FILES Constants.hpp + CommonMath.hpp ScalarTraits.hpp DESTINATION include/GridKit) diff --git a/GridKit/CommonMath.hpp b/GridKit/CommonMath.hpp new file mode 100644 index 000000000..b1b84ee8c --- /dev/null +++ b/GridKit/CommonMath.hpp @@ -0,0 +1,33 @@ +#pragma once + +#include + +#include +#include + +namespace GridKit +{ + namespace Math + { + /** + * @brief Scaled sigmoid activation function + * + * @note The sigmoid constant (mu) value of the constant to balance accuracy + * and finite derivatives. Large values more closely approximate a step + * function, but lead to inf or NaN derivatives. The value of 240 corresponds + * to a time step of 1/4 of a 60Hz cycle. + * + * @tparam ScalarT - scalar data type + * + * @param[in] x + * @return value of the sigmoid function + */ + template + __attribute__((always_inline)) inline ScalarT sigmoid(ScalarT x) + { + using RealT = typename GridKit::ScalarTraits::RealT; + static constexpr RealT MU = 240.0; + return ONE / (ONE + std::exp(-MU * x)); + } + } // namespace Math +} // namespace GridKit diff --git a/GridKit/Model/Evaluator.hpp b/GridKit/Model/Evaluator.hpp index 98789eba1..0200a8bce 100644 --- a/GridKit/Model/Evaluator.hpp +++ b/GridKit/Model/Evaluator.hpp @@ -3,6 +3,7 @@ #include #include +#include #include #include diff --git a/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1.hpp b/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1.hpp index 03ec0fd9c..528113181 100644 --- a/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1.hpp +++ b/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1.hpp @@ -150,13 +150,7 @@ namespace GridKit /// Component signal extension ComponentSignals signals_; - // Scale of Sigmoid function - // (temporary local implementation) - // mu set to 4*60, to provide accuracy and finite derivatives. - const RealT mu_{240.0}; - - // Activation function (sigmoid approximation) - ScalarT sigmoid(ScalarT x); + // Indicator ScalarT indicator(ScalarT x, ScalarT f); // Parameter initialization function diff --git a/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp b/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp index ff6613003..99ec409ed 100644 --- a/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp +++ b/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp @@ -236,26 +236,6 @@ namespace GridKit return 0; } - /** - * @brief Scaled sigmoid activation function - * - * @param[in] x State variable - * @tparam ScalarT Scalar data type - * @tparam IdxT Index data type - * @return Sigmoid approximation evaluated at x - * - * @warning This needs to be abstracted to be used - * across phasor dynamics. Identical pattern is - * being used in TGOV1 model. - * - * Algebraic approximation of transcendental sigmoid. - */ - template - ScalarT Ieeet1::sigmoid(ScalarT x) - { - return ONE / (ONE + std::exp(-mu_ * x)); - } - /** * @brief Net Indicator function for regulator limits * @@ -273,8 +253,8 @@ namespace GridKit ScalarT Ieeet1::indicator(ScalarT x, ScalarT f) { - ScalarT ind_low = (this->sigmoid(Vrmin_ - x)) * (this->sigmoid(-f)); - ScalarT ind_high = (this->sigmoid(x - Vrmax_)) * (this->sigmoid(f)); + ScalarT ind_low = (Math::sigmoid(Vrmin_ - x)) * (Math::sigmoid(-f)); + ScalarT ind_high = (Math::sigmoid(x - Vrmax_)) * (Math::sigmoid(f)); return (ONE - ind_low) * (ONE - ind_high); } @@ -340,7 +320,7 @@ namespace GridKit f_[6] = -ve + ksat * efdp; f_[7] = -efd + efdp + omega * efdp * Ispdlim_; - ScalarT efd_sat = (efdp - SA_) * (this->sigmoid(efdp - SA_)); + ScalarT efd_sat = (efdp - SA_) * (Math::sigmoid(efdp - SA_)); f_[8] = -ksat + SB_ * efd_sat * efd_sat; return 0; diff --git a/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1.hpp b/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1.hpp index 4597123d3..a5d31ec1e 100644 --- a/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1.hpp +++ b/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1.hpp @@ -121,16 +121,9 @@ namespace GridKit // Input States (which can be parameters) ScalarT pref_{0}; - // Scale of Sigmoid function (temporary local implementation) - // mu set to 4*60, to provide accuracy and finite derivatives. - const ScalarT mu_{240.0}; - /// Component signal extension ComponentSignals signals_; - // Activation function (sigmoid approximation) - ScalarT sigmoid(ScalarT x); - // Indicator of Valve limit states ScalarT indicator_low(ScalarT x, ScalarT f); ScalarT indicator_high(ScalarT x, ScalarT f); diff --git a/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1Impl.hpp b/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1Impl.hpp index ae4d771f9..58564668a 100644 --- a/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1Impl.hpp +++ b/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1Impl.hpp @@ -229,28 +229,13 @@ namespace GridKit return 0; } - /** - * @brief Scaled sigmoid activation function - * - * Temporary local implementation of smooth approximation - * of a piecewise differential equation. Ideally this is - * a more abstracted capability with GK. - * - * Algebraic approximation of transcendental sigmoid. - */ - template - ScalarT Tgov1::sigmoid(ScalarT x) - { - return ONE / (ONE + std::exp(-mu_ * x)); - } - /** * @brief Indicator function for lower valve limit violation. */ template ScalarT Tgov1::indicator_low(ScalarT x, ScalarT f) { - return (this->sigmoid(Pvmin_ - x)) * (this->sigmoid(-f)); + return (Math::sigmoid(Pvmin_ - x)) * (Math::sigmoid(-f)); } /** @@ -259,7 +244,7 @@ namespace GridKit template ScalarT Tgov1::indicator_high(ScalarT x, ScalarT f) { - return (this->sigmoid(x - Pvmax_)) * (this->sigmoid(f)); + return (Math::sigmoid(x - Pvmax_)) * (Math::sigmoid(f)); } /** From c58769711d87641431f4ea45f3008a6f1641d973 Mon Sep 17 00:00:00 2001 From: Nicholson Koukpaizan Date: Fri, 12 Dec 2025 11:57:52 -0500 Subject: [PATCH 08/12] Moved indicators to CommonMath.hpp. --- GridKit/CommonMath.hpp | 64 ++++++++++++++++++- .../PhasorDynamics/Exciter/IEEET1/Ieeet1.hpp | 3 - .../Exciter/IEEET1/Ieeet1Impl.hpp | 24 +------ .../PhasorDynamics/Governor/Tgov1/Tgov1.hpp | 6 +- .../Governor/Tgov1/Tgov1Impl.hpp | 29 +-------- 5 files changed, 66 insertions(+), 60 deletions(-) diff --git a/GridKit/CommonMath.hpp b/GridKit/CommonMath.hpp index b1b84ee8c..749dbdce4 100644 --- a/GridKit/CommonMath.hpp +++ b/GridKit/CommonMath.hpp @@ -23,11 +23,73 @@ namespace GridKit * @return value of the sigmoid function */ template - __attribute__((always_inline)) inline ScalarT sigmoid(ScalarT x) + __attribute__((always_inline)) inline ScalarT sigmoid(const ScalarT x) { using RealT = typename GridKit::ScalarTraits::RealT; static constexpr RealT MU = 240.0; return ONE / (ONE + std::exp(-MU * x)); } + + /** + * @brief Low indicator function for regulator limits + * + * @tparam ScalarT - Scalar data type + * @tparam RealT - Real data type (see GridKit::ScalarTraits::RealT) + * + * @param[in] limit_min - Minimum limit + * @param[in] x - State variable + * @param[in] f - Conditional derivative of state variable + * @return Scalar value indicating limit activation + */ + template + __attribute__((always_inline)) inline ScalarT indicator_low( + const RealT limit_min, + const ScalarT x, + const ScalarT f) + { + return sigmoid(limit_min - x) * sigmoid(-f); + } + + /** + * @brief High indicator function for regulator limits + * + * @tparam ScalarT - Scalar data type + * @tparam RealT - Real data type (see GridKit::ScalarTraits::RealT) + * + * @param[in] limit_max - Maximum limit + * @param[in] x - State variable + * @param[in] f - Conditional derivative of state variable + * @return Scalar value indicating limit activation + */ + template + __attribute__((always_inline)) inline ScalarT indicator_high( + const RealT limit_max, + const ScalarT x, + const ScalarT f) + { + return sigmoid(x - limit_max) * sigmoid(f); + } + + /** + * @brief Net Indicator function for regulator limits + * + * @tparam ScalarT - Scalar data type + * @tparam RealT - Real data type (see GridKit::ScalarTraits::RealT) + * + * @param[in] limit_min - Minimum limit + * @param[in] limit_max - Maximum limit + * @param[in] x - State variable + * @param[in] f - Conditional derivative of state variable + * @return Scalar value indicating limit activation + */ + template + __attribute__((always_inline)) inline ScalarT indicator( + const RealT limit_min, + const RealT limit_max, + const ScalarT x, + const ScalarT f) + { + return (ONE - indicator_low(limit_min, x, f)) * (ONE - indicator_high(limit_max, x, f)); + } } // namespace Math } // namespace GridKit diff --git a/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1.hpp b/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1.hpp index 528113181..ff20ca4c2 100644 --- a/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1.hpp +++ b/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1.hpp @@ -150,9 +150,6 @@ namespace GridKit /// Component signal extension ComponentSignals signals_; - // Indicator - ScalarT indicator(ScalarT x, ScalarT f); - // Parameter initialization function void initModelParams(const model_data_type& data); }; diff --git a/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp b/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp index 99ec409ed..17640d4a2 100644 --- a/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp +++ b/GridKit/Model/PhasorDynamics/Exciter/IEEET1/Ieeet1Impl.hpp @@ -236,28 +236,6 @@ namespace GridKit return 0; } - /** - * @brief Net Indicator function for regulator limits - * - * @param[in] x State variable - * @param[in] f Conditional derivative of state variable - * @tparam ScalarT Scalar data type - * @tparam IdxT Index data type - * @return Scalar value indicating limit activation. - * - * @warning This needs to be abstracted to be used - * across phasor dynamics. Identical pattern is - * being used in TGOV1 model. - */ - template - ScalarT Ieeet1::indicator(ScalarT x, ScalarT f) - { - - ScalarT ind_low = (Math::sigmoid(Vrmin_ - x)) * (Math::sigmoid(-f)); - ScalarT ind_high = (Math::sigmoid(x - Vrmax_)) * (Math::sigmoid(f)); - return (ONE - ind_low) * (ONE - ind_high); - } - /** * @brief Residuals of system equations * @@ -306,7 +284,7 @@ namespace GridKit // The 'pre-limit' derivative of Pv ScalarT f = -vr + Ka_ * vtr; - ScalarT vr_ind = this->indicator(vr, f); + ScalarT vr_ind = Math::indicator(Vrmin_, Vrmax_, vr, f); // Internal Differential Equations f_[0] = -vts_dot + (Ec_ - vts) / Tr_; diff --git a/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1.hpp b/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1.hpp index a5d31ec1e..95bb6fda5 100644 --- a/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1.hpp +++ b/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1.hpp @@ -124,11 +124,7 @@ namespace GridKit /// Component signal extension ComponentSignals signals_; - // Indicator of Valve limit states - ScalarT indicator_low(ScalarT x, ScalarT f); - ScalarT indicator_high(ScalarT x, ScalarT f); - ScalarT indicator(ScalarT x, ScalarT f); - + // Parameter initialization function void initializeParameters(const model_data_type& data); /* Local copies of external variables */ diff --git a/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1Impl.hpp b/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1Impl.hpp index 58564668a..f6add0df1 100644 --- a/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1Impl.hpp +++ b/GridKit/Model/PhasorDynamics/Governor/Tgov1/Tgov1Impl.hpp @@ -229,33 +229,6 @@ namespace GridKit return 0; } - /** - * @brief Indicator function for lower valve limit violation. - */ - template - ScalarT Tgov1::indicator_low(ScalarT x, ScalarT f) - { - return (Math::sigmoid(Pvmin_ - x)) * (Math::sigmoid(-f)); - } - - /** - * @brief Indicator function for high valve limit violation. - */ - template - ScalarT Tgov1::indicator_high(ScalarT x, ScalarT f) - { - return (Math::sigmoid(x - Pvmax_)) * (Math::sigmoid(f)); - } - - /** - * @brief Net Indicator function for valve limits. - */ - template - ScalarT Tgov1::indicator(ScalarT x, ScalarT f) - { - return (ONE - this->indicator_low(x, f)) * (ONE - this->indicator_high(x, f)); - } - /** * @brief Internal residuals * @@ -282,7 +255,7 @@ namespace GridKit // The 'pre-limit' derivative of Pv ScalarT func = (-pv + (pref_ - omega) / R_) / T1_; - ScalarT valv_ind = this->indicator(pv, func); + ScalarT valv_ind = Math::indicator(Pvmin_, Pvmax_, pv, func); // Internal Differential Equations f[0] = -ptx_dot + pv - (ptx + T2_ * pv) / T3_; From 150c8a57131dbbf1a1aea42516ca3e3a03927850 Mon Sep 17 00:00:00 2001 From: nkoukpaizan Date: Fri, 12 Dec 2025 16:58:59 +0000 Subject: [PATCH 09/12] Apply pre-commmit fixes --- GridKit/CommonMath.hpp | 16 ++++++++-------- GridKit/Model/Evaluator.hpp | 2 +- 2 files changed, 9 insertions(+), 9 deletions(-) diff --git a/GridKit/CommonMath.hpp b/GridKit/CommonMath.hpp index 749dbdce4..855bc3f3f 100644 --- a/GridKit/CommonMath.hpp +++ b/GridKit/CommonMath.hpp @@ -25,7 +25,7 @@ namespace GridKit template __attribute__((always_inline)) inline ScalarT sigmoid(const ScalarT x) { - using RealT = typename GridKit::ScalarTraits::RealT; + using RealT = typename GridKit::ScalarTraits::RealT; static constexpr RealT MU = 240.0; return ONE / (ONE + std::exp(-MU * x)); } @@ -43,8 +43,8 @@ namespace GridKit */ template __attribute__((always_inline)) inline ScalarT indicator_low( - const RealT limit_min, - const ScalarT x, + const RealT limit_min, + const ScalarT x, const ScalarT f) { return sigmoid(limit_min - x) * sigmoid(-f); @@ -63,8 +63,8 @@ namespace GridKit */ template __attribute__((always_inline)) inline ScalarT indicator_high( - const RealT limit_max, - const ScalarT x, + const RealT limit_max, + const ScalarT x, const ScalarT f) { return sigmoid(x - limit_max) * sigmoid(f); @@ -84,9 +84,9 @@ namespace GridKit */ template __attribute__((always_inline)) inline ScalarT indicator( - const RealT limit_min, - const RealT limit_max, - const ScalarT x, + const RealT limit_min, + const RealT limit_max, + const ScalarT x, const ScalarT f) { return (ONE - indicator_low(limit_min, x, f)) * (ONE - indicator_high(limit_max, x, f)); diff --git a/GridKit/Model/Evaluator.hpp b/GridKit/Model/Evaluator.hpp index 0200a8bce..12b760c6a 100644 --- a/GridKit/Model/Evaluator.hpp +++ b/GridKit/Model/Evaluator.hpp @@ -2,8 +2,8 @@ #include -#include #include +#include #include #include From 81b49f23a6dbae66c6f570c7ff0d2c3d98722d59 Mon Sep 17 00:00:00 2001 From: Nicholson Koukpaizan Date: Fri, 12 Dec 2025 12:03:30 -0500 Subject: [PATCH 10/12] [skip ci] Updated CHANGELOG for CommonMath.hpp. --- CHANGELOG.md | 1 + GridKit/CommonMath.hpp | 5 +++-- 2 files changed, 4 insertions(+), 2 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index d6076d9ff..161b37c3e 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -38,6 +38,7 @@ - Updated Jacobian value storage from `ScalarT` to `RealT`. - Added a header file defining constants to be used throughout the code. - Added `GridKitDocs` target for Doxygen documentation. +- Added a header file defining common math functions (e.g., sigmoid) to be used throughout the code. ## v0.1 diff --git a/GridKit/CommonMath.hpp b/GridKit/CommonMath.hpp index 855bc3f3f..1960a0532 100644 --- a/GridKit/CommonMath.hpp +++ b/GridKit/CommonMath.hpp @@ -14,8 +14,9 @@ namespace GridKit * * @note The sigmoid constant (mu) value of the constant to balance accuracy * and finite derivatives. Large values more closely approximate a step - * function, but lead to inf or NaN derivatives. The value of 240 corresponds - * to a time step of 1/4 of a 60Hz cycle. + * function, but lead to inf or NaN derivatives. + * The value of 240 corresponds to a time step of 1/4 of a 60Hz cycle. + * From experiments, numerical derivatives break down for mu above 700. * * @tparam ScalarT - scalar data type * From 1d6e4ae500828a4b44b3d1383d6ee71cf8e5f2af Mon Sep 17 00:00:00 2001 From: Nicholson Koukpaizan Date: Fri, 12 Dec 2025 12:15:51 -0500 Subject: [PATCH 11/12] Fixed sigmoid documentation. --- GridKit/CommonMath.hpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/GridKit/CommonMath.hpp b/GridKit/CommonMath.hpp index 1960a0532..d80eb4d86 100644 --- a/GridKit/CommonMath.hpp +++ b/GridKit/CommonMath.hpp @@ -12,7 +12,7 @@ namespace GridKit /** * @brief Scaled sigmoid activation function * - * @note The sigmoid constant (mu) value of the constant to balance accuracy + * @note The sigmoid constant (mu) value is chosen to balance accuracy * and finite derivatives. Large values more closely approximate a step * function, but lead to inf or NaN derivatives. * The value of 240 corresponds to a time step of 1/4 of a 60Hz cycle. From 8c92cc24b7e788a52720c239972811472bf5f184 Mon Sep 17 00:00:00 2001 From: nkoukpaizan Date: Fri, 12 Dec 2025 17:16:15 +0000 Subject: [PATCH 12/12] Apply pre-commmit fixes --- GridKit/CommonMath.hpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/GridKit/CommonMath.hpp b/GridKit/CommonMath.hpp index d80eb4d86..7866d8f17 100644 --- a/GridKit/CommonMath.hpp +++ b/GridKit/CommonMath.hpp @@ -14,7 +14,7 @@ namespace GridKit * * @note The sigmoid constant (mu) value is chosen to balance accuracy * and finite derivatives. Large values more closely approximate a step - * function, but lead to inf or NaN derivatives. + * function, but lead to inf or NaN derivatives. * The value of 240 corresponds to a time step of 1/4 of a 60Hz cycle. * From experiments, numerical derivatives break down for mu above 700. *