From 7678a79e7a100e5d39b9834e07d07e36fbc675ec Mon Sep 17 00:00:00 2001 From: lukelowry Date: Thu, 28 May 2026 19:01:57 -0500 Subject: [PATCH 1/2] Add HYGOV phasor dynamics governor --- CHANGELOG.md | 1 + GridKit/CommonMath.hpp | 30 +- GridKit/Model/PhasorDynamics/Component.hpp | 9 + .../Model/PhasorDynamics/ComponentLibrary.hpp | 1 + .../PhasorDynamics/Governor/CMakeLists.txt | 1 + .../Governor/HYGOV/CMakeLists.txt | 54 ++ .../PhasorDynamics/Governor/HYGOV/Hygov.cpp | 27 + .../PhasorDynamics/Governor/HYGOV/Hygov.hpp | 157 +++++ .../Governor/HYGOV/HygovData.hpp | 88 +++ .../HYGOV/HygovDependencyTracking.cpp | 27 + .../Governor/HYGOV/HygovEnzyme.cpp | 90 +++ .../Governor/HYGOV/HygovImpl.hpp | 603 ++++++++++++++++++ .../PhasorDynamics/Governor/HYGOV/README.md | 258 ++++++++ .../Model/PhasorDynamics/Governor/README.md | 1 + GridKit/Model/PhasorDynamics/INPUT_FORMAT.md | 1 + .../Model/PhasorDynamics/SystemModelData.hpp | 17 + .../SystemModelDataJSONParser.hpp | 6 + .../Model/PhasorDynamics/SystemModelImpl.hpp | 36 ++ docs/Figures/PhasorDynamics/HYGOV_diagram.png | Bin 0 -> 72381 bytes .../PhasorDynamics/Governor/HYGOV/README.md | 6 + .../Model/PhasorDynamics/Governor/README.md | 1 + .../Math/SmoothnessIndicatorTests.hpp | 50 +- .../Math/runSmoothnessIndicatorTests.cpp | 3 +- tests/UnitTests/PhasorDynamics/CMakeLists.txt | 10 + .../PhasorDynamics/GovernorHygovTests.hpp | 415 ++++++++++++ .../PhasorDynamics/runGovernorHygovTests.cpp | 17 + 26 files changed, 1892 insertions(+), 17 deletions(-) create mode 100644 GridKit/Model/PhasorDynamics/Governor/HYGOV/CMakeLists.txt create mode 100644 GridKit/Model/PhasorDynamics/Governor/HYGOV/Hygov.cpp create mode 100644 GridKit/Model/PhasorDynamics/Governor/HYGOV/Hygov.hpp create mode 100644 GridKit/Model/PhasorDynamics/Governor/HYGOV/HygovData.hpp create mode 100644 GridKit/Model/PhasorDynamics/Governor/HYGOV/HygovDependencyTracking.cpp create mode 100644 GridKit/Model/PhasorDynamics/Governor/HYGOV/HygovEnzyme.cpp create mode 100644 GridKit/Model/PhasorDynamics/Governor/HYGOV/HygovImpl.hpp create mode 100644 GridKit/Model/PhasorDynamics/Governor/HYGOV/README.md create mode 100755 docs/Figures/PhasorDynamics/HYGOV_diagram.png create mode 100644 docs/GridKit/Model/PhasorDynamics/Governor/HYGOV/README.md create mode 100644 tests/UnitTests/PhasorDynamics/GovernorHygovTests.hpp create mode 100644 tests/UnitTests/PhasorDynamics/runGovernorHygovTests.cpp diff --git a/CHANGELOG.md b/CHANGELOG.md index ba0f1dde6..6761bb82f 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -67,6 +67,7 @@ - Added `ConstantSignalSource` component as first use case for `BusToSignalAdapter`. - Added IDA option to suppress algebraic variables in local error tests. - Removed `COO_Matrix` class. +- Added `HYGOV` governor model implementation for PhasorDynamics. ## v0.1 diff --git a/GridKit/CommonMath.hpp b/GridKit/CommonMath.hpp index e4a7ec87f..8f9f8c5cf 100644 --- a/GridKit/CommonMath.hpp +++ b/GridKit/CommonMath.hpp @@ -162,7 +162,31 @@ namespace GridKit } /** - * @brief Smooth two-sided deadband function + * @brief Smooth Type 1 no-offset two-sided deadband function + * + * Smooth approximation to a deadband that returns zero inside the band and + * passes the input through unchanged outside the band. + * + * @tparam ScalarT - scalar data type + * @tparam RealT - Real data type (see GridKit::ScalarTraits::RealT) + * + * @param[in] x - Input signal + * @param[in] lower - Lower breakpoint + * @param[in] upper - Upper breakpoint + * @return Smooth no-offset deadbanded value + */ + template + __attribute__((always_inline)) inline ScalarT deadband1( + const ScalarT x, + const RealT lower, + const RealT upper) + { + assert(lower <= upper); + return x * (sigmoid(lower - x) + sigmoid(x - upper)); + } + + /** + * @brief Smooth Type 2 offset two-sided deadband function * * Smooth approximation to x - min(max(x, lower), upper), composed from the * smooth ramp function. @@ -173,10 +197,10 @@ namespace GridKit * @param[in] x - Input signal * @param[in] lower - Lower breakpoint * @param[in] upper - Upper breakpoint - * @return Smooth deadbanded value + * @return Smooth offset deadbanded value */ template - __attribute__((always_inline)) inline ScalarT deadband( + __attribute__((always_inline)) inline ScalarT deadband2( const ScalarT x, const RealT lower, const RealT upper) diff --git a/GridKit/Model/PhasorDynamics/Component.hpp b/GridKit/Model/PhasorDynamics/Component.hpp index 20dbd0403..8d2a33117 100644 --- a/GridKit/Model/PhasorDynamics/Component.hpp +++ b/GridKit/Model/PhasorDynamics/Component.hpp @@ -112,6 +112,15 @@ namespace GridKit return abs_tol_; } + int setAbsoluteTolerance(RealT rel_tol) override + { + for (IdxT i = 0; i < size_; ++i) + { + abs_tol_[static_cast(i)] = rel_tol; + } + return 0; + } + std::vector& getResidual() override { return f_; diff --git a/GridKit/Model/PhasorDynamics/ComponentLibrary.hpp b/GridKit/Model/PhasorDynamics/ComponentLibrary.hpp index d8fe3b2dc..64b81f1b2 100644 --- a/GridKit/Model/PhasorDynamics/ComponentLibrary.hpp +++ b/GridKit/Model/PhasorDynamics/ComponentLibrary.hpp @@ -7,6 +7,7 @@ #include #include #include +#include #include #include #include diff --git a/GridKit/Model/PhasorDynamics/Governor/CMakeLists.txt b/GridKit/Model/PhasorDynamics/Governor/CMakeLists.txt index 7c7269784..fbe71c740 100644 --- a/GridKit/Model/PhasorDynamics/Governor/CMakeLists.txt +++ b/GridKit/Model/PhasorDynamics/Governor/CMakeLists.txt @@ -4,3 +4,4 @@ # ]] add_subdirectory(Tgov1) +add_subdirectory(HYGOV) diff --git a/GridKit/Model/PhasorDynamics/Governor/HYGOV/CMakeLists.txt b/GridKit/Model/PhasorDynamics/Governor/HYGOV/CMakeLists.txt new file mode 100644 index 000000000..1719101b0 --- /dev/null +++ b/GridKit/Model/PhasorDynamics/Governor/HYGOV/CMakeLists.txt @@ -0,0 +1,54 @@ +# [[ +# Author(s): +# - Luke Lowery +# ]] + +set(_install_headers Hygov.hpp HygovData.hpp) + +if(GRIDKIT_ENABLE_ENZYME) + gridkit_add_library( + phasor_dynamics_governor_hygov + SOURCES HygovEnzyme.cpp + HEADERS ${_install_headers} + INCLUDE_DIRECTORIES PRIVATE ${GRIDKIT_THIRD_PARTY_DIR}/magic-enum/include + LINK_LIBRARIES + PUBLIC + GridKit::phasor_dynamics_core + PUBLIC + GridKit::phasor_dynamics_signal + PRIVATE + ClangEnzymeFlags + COMPILE_OPTIONS + PRIVATE + -mllvm + -enzyme-auto-sparsity=1 + -fno-math-errno) +else() + gridkit_add_library( + phasor_dynamics_governor_hygov + SOURCES Hygov.cpp + HEADERS ${_install_headers} + INCLUDE_DIRECTORIES PRIVATE ${GRIDKIT_THIRD_PARTY_DIR}/magic-enum/include + LINK_LIBRARIES + PUBLIC + GridKit::phasor_dynamics_core + PUBLIC + GridKit::phasor_dynamics_signal) +endif() + +gridkit_add_library( + phasor_dynamics_governor_hygov_dependency_tracking + SOURCES HygovDependencyTracking.cpp + INCLUDE_DIRECTORIES PRIVATE ${GRIDKIT_THIRD_PARTY_DIR}/magic-enum/include + LINK_LIBRARIES + PUBLIC + GridKit::phasor_dynamics_core + PUBLIC + GridKit::phasor_dynamics_signal_dependency_tracking) + +target_link_libraries( + phasor_dynamics_components + INTERFACE GridKit::phasor_dynamics_governor_hygov) +target_link_libraries( + phasor_dynamics_components_dependency_tracking + INTERFACE GridKit::phasor_dynamics_governor_hygov_dependency_tracking) diff --git a/GridKit/Model/PhasorDynamics/Governor/HYGOV/Hygov.cpp b/GridKit/Model/PhasorDynamics/Governor/HYGOV/Hygov.cpp new file mode 100644 index 000000000..bdb7e15c0 --- /dev/null +++ b/GridKit/Model/PhasorDynamics/Governor/HYGOV/Hygov.cpp @@ -0,0 +1,27 @@ +/** + * @file Hygov.cpp + * @author Luke Lowery (lukel@tamu.edu) + * @brief Non-Enzyme instantiation for the HYGOV governor model. + */ + +#include "HygovImpl.hpp" + +namespace GridKit +{ + namespace PhasorDynamics + { + namespace Governor + { + template + int Hygov::evaluateJacobian() + { + Log::misc() << "Evaluate Jacobian for Hygov..." << std::endl; + Log::misc() << "Jacobian evaluation is not implemented!" << std::endl; + return 0; + } + + template class Hygov; + template class Hygov; + } // namespace Governor + } // namespace PhasorDynamics +} // namespace GridKit diff --git a/GridKit/Model/PhasorDynamics/Governor/HYGOV/Hygov.hpp b/GridKit/Model/PhasorDynamics/Governor/HYGOV/Hygov.hpp new file mode 100644 index 000000000..bd75dc940 --- /dev/null +++ b/GridKit/Model/PhasorDynamics/Governor/HYGOV/Hygov.hpp @@ -0,0 +1,157 @@ +/** + * @file Hygov.hpp + * @author Luke Lowery (lukel@tamu.edu) + * @brief Declaration of the HYGOV governor model. + */ + +#pragma once + +#include +#include +#include +#include + +#include +#include +#include +#include + +namespace GridKit +{ + namespace PhasorDynamics + { + template + class SignalNode; + + namespace Governor + { + /// Internal variables of a `Hygov`. + enum class HygovInternalVariables : size_t + { + XN, ///< Speed lead-lag denominator state + XF, ///< Governor error filter output + C, ///< Desired-gate position + G, ///< Gate position + Q, ///< Turbine flow + OMEGADB, ///< Deadbanded speed deviation + EF, ///< Governor error into the filter + FC, ///< Desired-gate derivative target + RC, ///< Rate-limited desired-gate derivative target + PGV, ///< Nonlinear gate-to-power curve output + H, ///< Turbine head + PMECH, ///< Mechanical-power output + MAXIMUM, + }; + + /// External variables of a `Hygov`. + enum class HygovExternalVariables : size_t + { + OMEGA, ///< Machine speed deviation + PREF, ///< Active-power/load reference + PAUX, ///< Auxiliary power input + MAXIMUM, + }; + + template + class Hygov : public Component + { + using Component::alpha_; + using Component::abs_tol_; + using Component::f_; + using Component::gridkit_component_id_; + using Component::J_cols_buffer_; + using Component::J_rows_buffer_; + using Component::J_vals_buffer_; + using Component::nnz_; + using Component::residual_indices_; + using Component::size_; + using Component::tag_; + using Component::va_system_base_; + using Component::variable_indices_; + using Component::wb_; + using Component::y_; + using Component::yp_; + + public: + using ScalarT = scalar_type; + using IdxT = index_type; + using RealT = typename Component::RealT; + using SignalT = SignalNode; + using ModelDataT = HygovData; + using MonitorT = Model::VariableMonitor; + + Hygov(); + Hygov(const ModelDataT& data); + ~Hygov(); + + int setGridKitComponentID(IdxT) override final; + int allocate() override final; + int verify() const override final; + int initialize() override final; + int tagDifferentiable() override final; + int setAbsoluteTolerance(RealT) override final; + int evaluateResidual() override final; + int evaluateJacobian() override final; + + auto getSignals() + -> ComponentSignals& + { + return signals_; + } + + const Model::VariableMonitorBase* getMonitor() const override; + + __attribute__((always_inline)) inline int evaluateInternalResidual( + const ScalarT*, const ScalarT*, const ScalarT*, const ScalarT*, ScalarT*); + + private: + void initModelParams(const ModelDataT& data); + void initializeMonitor(); + void setDerivedParams(); + + ScalarT gatePower(ScalarT gate) const; + RealT invertGatePower(RealT pgv) const; + ScalarT toComponentBase(ScalarT value) const; + ScalarT toSystemBase(ScalarT value) const; + + static constexpr RealT TIME_CONSTANT_MINIMUM = static_cast(1.0e-3); + + RealT Trate_{0}; + RealT Rperm_{0}; + RealT Rtemp_{0}; + RealT Tr_{0}; + RealT Tf_{0}; + RealT Tg_{0}; + RealT Velm_{0}; + RealT Gmax_{0}; + RealT Gmin_{0}; + RealT Tw_{0}; + RealT At_{0}; + RealT Dturb_{0}; + RealT Qnl_{0}; + RealT Tn_{0}; + RealT Tnp_{0}; + RealT leadlag_gain_{0}; + RealT db1_{0}; + RealT db2_{0}; + RealT Hdam_{1}; + std::array Gv_{}; + std::array Pgv_{}; + + RealT va_component_base_{0}; + + ScalarT pref_set_{0}; + ScalarT paux_set_{0}; + + ComponentSignals signals_; + std::unique_ptr monitor_; + + std::vector ws_; + std::vector ws_indices_; + }; + } // namespace Governor + } // namespace PhasorDynamics +} // namespace GridKit diff --git a/GridKit/Model/PhasorDynamics/Governor/HYGOV/HygovData.hpp b/GridKit/Model/PhasorDynamics/Governor/HYGOV/HygovData.hpp new file mode 100644 index 000000000..5c44a03a4 --- /dev/null +++ b/GridKit/Model/PhasorDynamics/Governor/HYGOV/HygovData.hpp @@ -0,0 +1,88 @@ +/** + * @file HygovData.hpp + * @author Luke Lowery (lukel@tamu.edu) + * @brief Modeling data for the HYGOV governor model. + */ + +#pragma once + +#include + +namespace GridKit +{ + namespace PhasorDynamics + { + namespace Governor + { + /// Parameter keys for the HYGOV governor model. + enum class HygovParameters + { + Trate, ///< Turbine-rating power base + Rperm, ///< Permanent droop + Rtemp, ///< Temporary droop + Tr, ///< Temporary-droop reset time constant + Tf, ///< Governor error filter time constant + Tg, ///< Gate servo time constant + Velm, ///< Maximum desired-gate velocity magnitude + Gmax, ///< Maximum desired-gate position + Gmin, ///< Minimum desired-gate position + Tw, ///< Water inertia time constant + At, ///< Turbine gain + Dturb, ///< Turbine damping coefficient + Qnl, ///< No-load flow at nominal head + Tn, ///< Speed lead-lag numerator time constant + Tnp, ///< Speed lead-lag denominator time constant + db1, ///< Type 1 speed deadband threshold + db2, ///< Unsupported mechanical backlash deadband + Hdam, ///< Head available at dam + Gv0, ///< Gate point 0 + Gv1, ///< Gate point 1 + Gv2, ///< Gate point 2 + Gv3, ///< Gate point 3 + Gv4, ///< Gate point 4 + Gv5, ///< Gate point 5 + Pgv0, ///< Power point 0 + Pgv1, ///< Power point 1 + Pgv2, ///< Power point 2 + Pgv3, ///< Power point 3 + Pgv4, ///< Power point 4 + Pgv5 ///< Power point 5 + }; + + /// Ports for the HYGOV governor model. + enum class HygovPorts + { + bus, ///< Optional terminal bus ID for case-format compatibility + speed, ///< Machine speed-deviation signal ID + pmech, ///< Mechanical-power output signal ID + pref, ///< Optional active-power/load reference signal ID + paux ///< Optional auxiliary power input signal ID + }; + + /// Variables available through the monitor interface. + enum class HygovMonitorableVariables + { + pmech, ///< Mechanical power output + filter, ///< Governor error filter output + desiredgate, ///< Desired-gate position + gate, ///< Gate position + flow, ///< Turbine flow + head ///< Turbine head + }; + + template + struct HygovData : public ComponentData + { + HygovData() = default; + + using Parameters = HygovParameters; + using Ports = HygovPorts; + using MonitorableVariables = HygovMonitorableVariables; + }; + } // namespace Governor + } // namespace PhasorDynamics +} // namespace GridKit diff --git a/GridKit/Model/PhasorDynamics/Governor/HYGOV/HygovDependencyTracking.cpp b/GridKit/Model/PhasorDynamics/Governor/HYGOV/HygovDependencyTracking.cpp new file mode 100644 index 000000000..760bac957 --- /dev/null +++ b/GridKit/Model/PhasorDynamics/Governor/HYGOV/HygovDependencyTracking.cpp @@ -0,0 +1,27 @@ +/** + * @file HygovDependencyTracking.cpp + * @author Luke Lowery (lukel@tamu.edu) + * @brief Dependency-tracking instantiations for the HYGOV governor model. + */ + +#include "HygovImpl.hpp" + +namespace GridKit +{ + namespace PhasorDynamics + { + namespace Governor + { + template + int Hygov::evaluateJacobian() + { + Log::misc() << "Evaluate Jacobian for Hygov..." << std::endl; + Log::misc() << "Jacobian evaluation is not implemented!" << std::endl; + return 0; + } + + template class Hygov; + template class Hygov; + } // namespace Governor + } // namespace PhasorDynamics +} // namespace GridKit diff --git a/GridKit/Model/PhasorDynamics/Governor/HYGOV/HygovEnzyme.cpp b/GridKit/Model/PhasorDynamics/Governor/HYGOV/HygovEnzyme.cpp new file mode 100644 index 000000000..cc6755968 --- /dev/null +++ b/GridKit/Model/PhasorDynamics/Governor/HYGOV/HygovEnzyme.cpp @@ -0,0 +1,90 @@ +/** + * @file HygovEnzyme.cpp + * @author Luke Lowery (lukel@tamu.edu) + * @brief Enzyme sparse Jacobian for the HYGOV governor model. + */ + +#include + +#include "HygovImpl.hpp" + +namespace GridKit +{ + namespace PhasorDynamics + { + namespace Governor + { + template + int Hygov::evaluateJacobian() + { + Log::misc() << "Evaluate Jacobian for Hygov..." << std::endl; + Log::misc() << "Jacobian evaluation is experimental!" << std::endl; + + if (J_rows_buffer_ == nullptr) + { + auto size = static_cast(size_); + auto signal_size = static_cast(ws_.size()); + auto buffer_size = 2 * size * size + size * signal_size; + J_rows_buffer_ = new IdxT[buffer_size]; + J_cols_buffer_ = new IdxT[buffer_size]; + J_vals_buffer_ = new RealT[buffer_size]; + } + + nnz_ = 0; + + GridKit::Enzyme::Sparse::DfDy, + GridKit::Enzyme::Sparse::MemberFunctions::InternalResidualWithSignal>::eval(this, + f_.size(), + y_.size(), + (this->getResidualIndices()).data(), + (this->getVariableIndices()).data(), + y_.data(), + yp_.data(), + wb_.data(), + ws_.data(), + J_rows_buffer_, + J_cols_buffer_, + J_vals_buffer_, + nnz_); + + GridKit::Enzyme::Sparse::DfDyp, + GridKit::Enzyme::Sparse::MemberFunctions::InternalResidualWithSignal>::eval(this, + f_.size(), + y_.size(), + (this->getResidualIndices()).data(), + (this->getVariableIndices()).data(), + y_.data(), + yp_.data(), + wb_.data(), + ws_.data(), + alpha_, + J_rows_buffer_, + J_cols_buffer_, + J_vals_buffer_, + nnz_); + + GridKit::Enzyme::Sparse::DfDws, + GridKit::Enzyme::Sparse::MemberFunctions::InternalResidualWithSignal>::eval(this, + f_.size(), + ws_.size(), + (this->getResidualIndices()).data(), + ws_indices_.data(), + y_.data(), + yp_.data(), + wb_.data(), + ws_.data(), + J_rows_buffer_, + J_cols_buffer_, + J_vals_buffer_, + nnz_); + + this->constructCoo(); + + return 0; + } + + template class Hygov; + template class Hygov; + } // namespace Governor + } // namespace PhasorDynamics +} // namespace GridKit diff --git a/GridKit/Model/PhasorDynamics/Governor/HYGOV/HygovImpl.hpp b/GridKit/Model/PhasorDynamics/Governor/HYGOV/HygovImpl.hpp new file mode 100644 index 000000000..cfb48f7cf --- /dev/null +++ b/GridKit/Model/PhasorDynamics/Governor/HYGOV/HygovImpl.hpp @@ -0,0 +1,603 @@ +/** + * @file HygovImpl.hpp + * @author Luke Lowery (lukel@tamu.edu) + * @brief Definition of the HYGOV governor model. + */ + +#pragma once + +#include +#include +#include + +#include +#include +#include +#include +#include + +namespace GridKit +{ + namespace PhasorDynamics + { + namespace Governor + { + using Log = ::GridKit::Utilities::Logger; + + template + Hygov::Hygov() + { + size_ = static_cast(HygovInternalVariables::MAXIMUM); + setDerivedParams(); + } + + template + Hygov::Hygov(const ModelDataT& data) + : monitor_(std::make_unique(data)) + { + initModelParams(data); + initializeMonitor(); + + size_ = static_cast(HygovInternalVariables::MAXIMUM); + setDerivedParams(); + } + + template + Hygov::~Hygov() + { + } + + template + void Hygov::setDerivedParams() + { + va_component_base_ = Trate_ * static_cast(1.0e6); + } + + template + scalar_type Hygov::toComponentBase(scalar_type value) const + { + return value * va_system_base_ / va_component_base_; + } + + template + scalar_type Hygov::toSystemBase(scalar_type value) const + { + return value / toComponentBase(static_cast(ONE)); + } + + template + void Hygov::initModelParams(const ModelDataT& data) + { + using Params = typename ModelDataT::Parameters; + + Trate_ = ZERO; + Rperm_ = static_cast(0.04); + Rtemp_ = static_cast(0.3); + Tr_ = static_cast(5.0); + Tf_ = static_cast(0.05); + Tg_ = static_cast(0.5); + Velm_ = static_cast(0.2); + Gmax_ = ONE; + Gmin_ = ZERO; + Tw_ = ONE; + At_ = static_cast(1.2); + Dturb_ = static_cast(0.5); + Qnl_ = static_cast(0.05); + Tn_ = ZERO; + Tnp_ = ZERO; + db1_ = ZERO; + db2_ = ZERO; + Hdam_ = ONE; + Gv_.fill(ZERO); + Pgv_.fill(ZERO); + + if (data.parameters.contains(Params::Trate)) + { + Trate_ = std::get(data.parameters.at(Params::Trate)); + } + if (data.parameters.contains(Params::Rperm)) + { + Rperm_ = std::get(data.parameters.at(Params::Rperm)); + } + if (data.parameters.contains(Params::Rtemp)) + { + Rtemp_ = std::get(data.parameters.at(Params::Rtemp)); + } + if (data.parameters.contains(Params::Tr)) + { + Tr_ = std::get(data.parameters.at(Params::Tr)); + } + if (data.parameters.contains(Params::Tf)) + { + Tf_ = std::get(data.parameters.at(Params::Tf)); + } + if (data.parameters.contains(Params::Tg)) + { + Tg_ = std::get(data.parameters.at(Params::Tg)); + } + if (data.parameters.contains(Params::Velm)) + { + Velm_ = std::get(data.parameters.at(Params::Velm)); + } + if (data.parameters.contains(Params::Gmax)) + { + Gmax_ = std::get(data.parameters.at(Params::Gmax)); + } + if (data.parameters.contains(Params::Gmin)) + { + Gmin_ = std::get(data.parameters.at(Params::Gmin)); + } + if (data.parameters.contains(Params::Tw)) + { + Tw_ = std::get(data.parameters.at(Params::Tw)); + } + if (data.parameters.contains(Params::At)) + { + At_ = std::get(data.parameters.at(Params::At)); + } + if (data.parameters.contains(Params::Dturb)) + { + Dturb_ = std::get(data.parameters.at(Params::Dturb)); + } + if (data.parameters.contains(Params::Qnl)) + { + Qnl_ = std::get(data.parameters.at(Params::Qnl)); + } + if (data.parameters.contains(Params::Tn)) + { + Tn_ = std::get(data.parameters.at(Params::Tn)); + } + if (data.parameters.contains(Params::Tnp)) + { + Tnp_ = std::get(data.parameters.at(Params::Tnp)); + } + if (data.parameters.contains(Params::db1)) + { + db1_ = std::get(data.parameters.at(Params::db1)); + } + if (data.parameters.contains(Params::db2)) + { + db2_ = std::get(data.parameters.at(Params::db2)); + } + if (data.parameters.contains(Params::Hdam)) + { + Hdam_ = std::get(data.parameters.at(Params::Hdam)); + } + if (data.parameters.contains(Params::Gv0)) + { + Gv_[0] = std::get(data.parameters.at(Params::Gv0)); + } + if (data.parameters.contains(Params::Gv1)) + { + Gv_[1] = std::get(data.parameters.at(Params::Gv1)); + } + if (data.parameters.contains(Params::Gv2)) + { + Gv_[2] = std::get(data.parameters.at(Params::Gv2)); + } + if (data.parameters.contains(Params::Gv3)) + { + Gv_[3] = std::get(data.parameters.at(Params::Gv3)); + } + if (data.parameters.contains(Params::Gv4)) + { + Gv_[4] = std::get(data.parameters.at(Params::Gv4)); + } + if (data.parameters.contains(Params::Gv5)) + { + Gv_[5] = std::get(data.parameters.at(Params::Gv5)); + } + if (data.parameters.contains(Params::Pgv0)) + { + Pgv_[0] = std::get(data.parameters.at(Params::Pgv0)); + } + if (data.parameters.contains(Params::Pgv1)) + { + Pgv_[1] = std::get(data.parameters.at(Params::Pgv1)); + } + if (data.parameters.contains(Params::Pgv2)) + { + Pgv_[2] = std::get(data.parameters.at(Params::Pgv2)); + } + if (data.parameters.contains(Params::Pgv3)) + { + Pgv_[3] = std::get(data.parameters.at(Params::Pgv3)); + } + if (data.parameters.contains(Params::Pgv4)) + { + Pgv_[4] = std::get(data.parameters.at(Params::Pgv4)); + } + if (data.parameters.contains(Params::Pgv5)) + { + Pgv_[5] = std::get(data.parameters.at(Params::Pgv5)); + } + + const bool source_default_curve = + std::all_of(Gv_.begin(), Gv_.end(), [](RealT value) + { return value == ZERO; }) + && std::all_of(Pgv_.begin(), Pgv_.end(), [](RealT value) + { return value == ZERO; }); + if (source_default_curve) + { + Gv_ = {ZERO, + static_cast(0.2), + static_cast(0.4), + static_cast(0.6), + static_cast(0.8), + ONE}; + Pgv_ = Gv_; + } + + Tr_ = std::max(Tr_, TIME_CONSTANT_MINIMUM); + Tf_ = std::max(Tf_, TIME_CONSTANT_MINIMUM); + Tg_ = std::max(Tg_, TIME_CONSTANT_MINIMUM); + Tw_ = std::max(Tw_, TIME_CONSTANT_MINIMUM); + Tnp_ = std::max(Tnp_, TIME_CONSTANT_MINIMUM); + + leadlag_gain_ = Tn_ / Tnp_; + } + + template + const Model::VariableMonitorBase* Hygov::getMonitor() const + { + return monitor_.get(); + } + + template + void Hygov::initializeMonitor() + { + using Variable = typename ModelDataT::MonitorableVariables; + auto index = [](HygovInternalVariables variable) + { + return static_cast(variable); + }; + + monitor_->set(Variable::pmech, [this, index] + { return y_[index(HygovInternalVariables::PMECH)]; }); + monitor_->set(Variable::filter, [this, index] + { return y_[index(HygovInternalVariables::XF)]; }); + monitor_->set(Variable::desiredgate, [this, index] + { return y_[index(HygovInternalVariables::C)]; }); + monitor_->set(Variable::gate, [this, index] + { return y_[index(HygovInternalVariables::G)]; }); + monitor_->set(Variable::flow, [this, index] + { return y_[index(HygovInternalVariables::Q)]; }); + monitor_->set(Variable::head, [this, index] + { return y_[index(HygovInternalVariables::H)]; }); + } + + template + int Hygov::setGridKitComponentID(IdxT component_id) + { + gridkit_component_id_ = component_id; + return 0; + } + + template + int Hygov::allocate() + { + size_ = static_cast(HygovInternalVariables::MAXIMUM); + auto size = static_cast(size_); + + f_.resize(size); + y_.resize(size); + yp_.resize(size); + tag_.resize(size); + abs_tol_.resize(size); + variable_indices_.resize(size); + residual_indices_.resize(size); + + auto signal_size = static_cast(HygovExternalVariables::MAXIMUM); + ws_.resize(signal_size); + ws_indices_.resize(signal_size); + std::fill(ws_.begin(), ws_.end(), ScalarT{0}); + std::fill(ws_indices_.begin(), ws_indices_.end(), INVALID_INDEX); + + for (IdxT j = 0; j < size_; ++j) + { + this->setVariableIndex(j, j); + this->setResidualIndex(j, j); + } + + if (signals_.template isAssigned()) + { + signals_.template getSignalNode()->set( + &y_[static_cast(HygovInternalVariables::PMECH)], + &(this->getVariableIndex(static_cast(HygovInternalVariables::PMECH)))); + } + + return 0; + } + + template + int Hygov::verify() const + { + int ret = 0; + + auto check = [&](bool condition, const char* message) + { + if (!condition) + { + Log::error() << "Hygov: " << message << '\n'; + ret += 1; + } + }; + + check(Trate_ > ZERO, "Trate must be positive"); + check(Rtemp_ != ZERO, "Rtemp must be nonzero"); + check(Velm_ >= ZERO, "Velm must be non-negative"); + check(Gmin_ <= Gmax_, "Gmin must be less than or equal to Gmax"); + check(db1_ >= ZERO, "db1 must be non-negative"); + check(Hdam_ > ZERO, "Hdam must be positive"); + + for (size_t i = 1; i < Gv_.size(); ++i) + { + check(Gv_[i - 1] < Gv_[i], "Gv points must be strictly increasing"); + check(Pgv_[i - 1] <= Pgv_[i], "Pgv points must be non-decreasing"); + } + + if (signals_.template isAttached() + && !signals_.template isLinked()) + { + Log::error() << "Hygov: omega signal attached with no linked source\n"; + ret += 1; + } + + if (signals_.template isAttached() + && !signals_.template isLinked()) + { + Log::error() << "Hygov: pref signal attached with no linked source\n"; + ret += 1; + } + + if (signals_.template isAttached() + && !signals_.template isLinked()) + { + Log::error() << "Hygov: paux signal attached with no linked source\n"; + ret += 1; + } + + return ret; + } + + template + scalar_type Hygov::gatePower(scalar_type gate) const + { + return ScalarT{Pgv_[0]} + + Math::linseg(gate, Gv_[0], Gv_[1], Pgv_[1] - Pgv_[0]) + + Math::linseg(gate, Gv_[1], Gv_[2], Pgv_[2] - Pgv_[1]) + + Math::linseg(gate, Gv_[2], Gv_[3], Pgv_[3] - Pgv_[2]) + + Math::linseg(gate, Gv_[3], Gv_[4], Pgv_[4] - Pgv_[3]) + + Math::linseg(gate, Gv_[4], Gv_[5], Pgv_[5] - Pgv_[4]); + } + + template + typename Hygov::RealT + Hygov::invertGatePower( + typename Hygov::RealT pgv) const + { + static constexpr RealT tol = static_cast(1.0e-10); + + if (std::abs(pgv - Pgv_[0]) <= tol) + { + return Gv_[0]; + } + + for (size_t i = 0; i < 5; ++i) + { + if (Pgv_[i + 1] <= Pgv_[i]) + { + continue; + } + + if (Pgv_[i] - tol <= pgv && pgv <= Pgv_[i + 1] + tol) + { + const RealT fraction = (pgv - Pgv_[i]) / (Pgv_[i + 1] - Pgv_[i]); + return Gv_[i] + fraction * (Gv_[i + 1] - Gv_[i]); + } + } + + return std::numeric_limits::quiet_NaN(); + } + + template + int Hygov::initialize() + { + if (verify() > 0) + { + Log::error() << "Hygov: cannot initialize with invalid configuration\n"; + return 1; + } + + const auto XN = static_cast(HygovInternalVariables::XN); + const auto XF = static_cast(HygovInternalVariables::XF); + const auto C = static_cast(HygovInternalVariables::C); + const auto G = static_cast(HygovInternalVariables::G); + const auto Q = static_cast(HygovInternalVariables::Q); + const auto OMEGADB = static_cast(HygovInternalVariables::OMEGADB); + const auto EF = static_cast(HygovInternalVariables::EF); + const auto FC = static_cast(HygovInternalVariables::FC); + const auto RC = static_cast(HygovInternalVariables::RC); + const auto PGV = static_cast(HygovInternalVariables::PGV); + const auto H = static_cast(HygovInternalVariables::H); + const auto PMECH = static_cast(HygovInternalVariables::PMECH); + + ScalarT omega0{ZERO}; + if (signals_.template isAttached()) + { + omega0 = signals_.template readExternalVariable(); + } + + paux_set_ = ScalarT{ZERO}; + if (signals_.template isAttached()) + { + paux_set_ = signals_.template readExternalVariable(); + } + + const ScalarT pmech0 = toComponentBase(y_[PMECH]); + y_[H] = Hdam_; + y_[Q] = Qnl_ + pmech0 / (At_ * y_[H]); + y_[PGV] = y_[Q] / std::sqrt(y_[H]); + + const RealT gate0 = invertGatePower(static_cast(y_[PGV])); + if (std::isnan(gate0)) + { + Log::error() << "Hygov: initial Pgv is outside the invertible gate curve\n"; + return 1; + } + + y_[G] = gate0; + y_[C] = y_[G]; + + if (y_[C] < Gmin_ || y_[C] > Gmax_) + { + Log::error() << "Hygov: initialized gate is outside Gmin/Gmax\n"; + return 1; + } + + y_[OMEGADB] = Math::deadband1(omega0, -db1_, db1_); + y_[XN] = y_[OMEGADB]; + y_[XF] = ZERO; + y_[EF] = ZERO; + y_[FC] = ZERO; + y_[RC] = ZERO; + y_[PMECH] = toSystemBase(pmech0); + + const ScalarT yomega = y_[XN] + leadlag_gain_ * (y_[OMEGADB] - y_[XN]); + pref_set_ = y_[EF] - paux_set_ + yomega + Rperm_ * y_[C]; + if (signals_.template isAttached()) + { + const ScalarT pref0 = + signals_.template readExternalVariable(); + const RealT pref_err = static_cast(pref0 - pref_set_); + if (std::abs(pref_err) > static_cast(1.0e-10)) + { + Log::error() << "Hygov: pref initial condition is not steady state\n"; + return 1; + } + pref_set_ = pref0; + } + + std::fill(yp_.begin(), yp_.end(), ZERO); + return 0; + } + + template + int Hygov::tagDifferentiable() + { + std::fill(tag_.begin(), tag_.end(), false); + tag_[static_cast(HygovInternalVariables::XN)] = true; + tag_[static_cast(HygovInternalVariables::XF)] = true; + tag_[static_cast(HygovInternalVariables::C)] = true; + tag_[static_cast(HygovInternalVariables::G)] = true; + tag_[static_cast(HygovInternalVariables::Q)] = true; + return 0; + } + + template + int Hygov::setAbsoluteTolerance(RealT rel_tol) + { + std::fill(abs_tol_.begin(), abs_tol_.end(), rel_tol); + return 0; + } + + template + __attribute__((always_inline)) inline int + Hygov::evaluateInternalResidual( + const ScalarT* y, + const ScalarT* yp, + [[maybe_unused]] const ScalarT* wb, + const ScalarT* ws, + ScalarT* f) + { + const auto XN = static_cast(HygovInternalVariables::XN); + const auto XF = static_cast(HygovInternalVariables::XF); + const auto C = static_cast(HygovInternalVariables::C); + const auto G = static_cast(HygovInternalVariables::G); + const auto Q = static_cast(HygovInternalVariables::Q); + const auto OMEGADB = static_cast(HygovInternalVariables::OMEGADB); + const auto EF = static_cast(HygovInternalVariables::EF); + const auto FC = static_cast(HygovInternalVariables::FC); + const auto RC = static_cast(HygovInternalVariables::RC); + const auto PGV = static_cast(HygovInternalVariables::PGV); + const auto H = static_cast(HygovInternalVariables::H); + const auto PMECH = static_cast(HygovInternalVariables::PMECH); + + const auto OMEGA = static_cast(HygovExternalVariables::OMEGA); + const auto PREF = static_cast(HygovExternalVariables::PREF); + const auto PAUX = static_cast(HygovExternalVariables::PAUX); + + const ScalarT xn = y[XN]; + const ScalarT xf = y[XF]; + const ScalarT c = y[C]; + const ScalarT g = y[G]; + const ScalarT q = y[Q]; + const ScalarT omegadb = y[OMEGADB]; + const ScalarT ef = y[EF]; + const ScalarT fc = y[FC]; + const ScalarT rc = y[RC]; + const ScalarT pgv = y[PGV]; + const ScalarT head = y[H]; + const ScalarT pmech = y[PMECH]; + + const ScalarT omega = ws[OMEGA]; + const ScalarT pref = ws[PREF]; + const ScalarT paux = ws[PAUX]; + + const ScalarT yomega = xn + leadlag_gain_ * (omegadb - xn); + + f[XN] = -yp[XN] + (omegadb - xn) / Tnp_; + f[XF] = -yp[XF] + (ef - xf) / Tf_; + f[C] = -yp[C] + Math::antiwindup(c, rc, Gmin_, Gmax_); + f[G] = -yp[G] + (c - g) / Tg_; + f[Q] = -yp[Q] + (Hdam_ - head) / Tw_; + f[OMEGADB] = -omegadb + Math::deadband1(omega, -db1_, db1_); + f[EF] = -ef + pref + paux - yomega - Rperm_ * c; + f[FC] = -fc + (xf / Tr_ + (ef - xf) / Tf_) / Rtemp_; + f[RC] = -rc + Math::clamp(fc, -Velm_, Velm_); + f[PGV] = -pgv + gatePower(g); + f[H] = -q * q + head * pgv * pgv; + f[PMECH] = -toComponentBase(pmech) + At_ * head * (q - Qnl_) - Dturb_ * omega * g; + + return 0; + } + + template + int Hygov::evaluateResidual() + { + const auto OMEGA = static_cast(HygovExternalVariables::OMEGA); + const auto PREF = static_cast(HygovExternalVariables::PREF); + const auto PAUX = static_cast(HygovExternalVariables::PAUX); + + ws_[OMEGA] = ZERO; + ws_[PREF] = pref_set_; + ws_[PAUX] = paux_set_; + std::fill(ws_indices_.begin(), ws_indices_.end(), INVALID_INDEX); + + if (signals_.template isAttached()) + { + ws_[OMEGA] = signals_.template readExternalVariable(); + ws_indices_[OMEGA] = + signals_.template readExternalVariableIndex(); + } + + if (signals_.template isAttached()) + { + ws_[PREF] = signals_.template readExternalVariable(); + ws_indices_[PREF] = + signals_.template readExternalVariableIndex(); + } + + if (signals_.template isAttached()) + { + ws_[PAUX] = signals_.template readExternalVariable(); + ws_indices_[PAUX] = + signals_.template readExternalVariableIndex(); + } + + evaluateInternalResidual(y_.data(), yp_.data(), wb_.data(), ws_.data(), f_.data()); + return 0; + } + } // namespace Governor + } // namespace PhasorDynamics +} // namespace GridKit diff --git a/GridKit/Model/PhasorDynamics/Governor/HYGOV/README.md b/GridKit/Model/PhasorDynamics/Governor/HYGOV/README.md new file mode 100644 index 000000000..07df81da7 --- /dev/null +++ b/GridKit/Model/PhasorDynamics/Governor/HYGOV/README.md @@ -0,0 +1,258 @@ +# **Hydro Turbine-Governor Model (HYGOV)** + +HYGOV is a hydro turbine-governor model with temporary droop, a gate servo, and +a nonlinear single-penstock turbine. It reads machine speed deviation and +outputs mechanical power. + +Notes: +- PowerWorld uses the connected machine base when `Trate = 0`; GridKit requires + `Trate` to be set explicitly. +- HYGOVD `dbL`/`dbH`, `db2` backlash, and Kaplan blade-servo fields are not + modeled. + +## Block Diagram + +Standard HYGOV block diagram. + +![](../../../../../docs/Figures/PhasorDynamics/HYGOV_diagram.png) + +Figure 1: Governor HYGOV model. Figure courtesy of [PowerWorld](https://www.powerworld.com/WebHelp/) + +## Model Parameters + +Symbol | Units | JSON | Description | Typical Value | Note +--------------------------------|----------|----------|----------------------------------------------|---------------|------ +$T_{\mathrm{rate}}$ | [MW] | `Trate` | Component-base turbine rating | 100.0 | Required positive value +$R_{\mathrm{perm}}$ | [p.u.] | `Rperm` | Permanent droop | 0.04 | Source diagram label: `R` +$R_{\mathrm{temp}}$ | [p.u.] | `Rtemp` | Temporary droop | 0.3 | Source diagram label: `r` +$T_r$ | [sec] | `Tr` | Temporary-droop reset time constant | 5.0 | +$T_f$ | [sec] | `Tf` | Governor error filter time constant | 0.05 | State 1 +$T_g$ | [sec] | `Tg` | Gate servo time constant | 0.5 | State 3 +$V_{\mathrm{elm}}$ | [p.u./s] | `Velm` | Maximum desired-gate velocity magnitude | 0.2 | Symmetric rate limit on State 2 +$G^{\max}$ | [p.u.] | `Gmax` | Maximum desired-gate position | 1.0 | +$G^{\min}$ | [p.u.] | `Gmin` | Minimum desired-gate position | 0.0 | +$T_w$ | [sec] | `Tw` | Water inertia time constant | 1.0 | State 4 +$A_t$ | [p.u.] | `At` | Turbine gain | 1.2 | +$D_{\mathrm{turb}}$ | [p.u.] | `Dturb` | Turbine damping coefficient | 0.5 | Multiplied by speed deviation and gate +$q_{\mathrm{NL}}$ | [p.u.] | `Qnl` | No-load flow at nominal head | 0.05 | +$T_n$ | [sec] | `Tn` | Speed lead-lag numerator time constant | 0.0 | +$T_{np}$ | [sec] | `Tnp` | Speed lead-lag denominator time constant | 0.0 | +$D_{\omega}$ | [p.u.] | `db1` | Type 1 speed deadband threshold | 0.0 | Uses CommonMath `deadband1` +$H_{\mathrm{dam}}$ | [p.u.] | `Hdam` | Head available at dam | 1.0 | +$G_V^{(k)}$ | [p.u.] | `Gv0`-`Gv5` | Gate point $k$ of the gain curve | 0.0 | $k=0,\ldots,5$ +$P_{\mathrm{GV}}^{(k)}$ | [p.u.] | `Pgv0`-`Pgv5` | Power point $k$ of the gain curve | 0.0 | $k=0,\ldots,5$ + +All-zero `Gv` and `Pgv` source points select the identity curve. + +### Parameter Validation + +Invalid HYGOV parameter sets are rejected by the following checks. The displayed +equations use effective time constants with $\epsilon_T=10^{-3}$. + +```math +\begin{aligned} + T &\leftarrow \max\!\left(T, \epsilon_T\right) + \quad T\in\{T_r,T_f,T_g,T_w,T_{np}\} \\ + T_{\mathrm{rate}}, H_{\mathrm{dam}} + &> 0 \\ + R_{\mathrm{temp}} + &\ne 0 \\ + V_{\mathrm{elm}}, D_{\omega} + &\ge 0 \\ + G^{\min} &\le G^{\max} \\ + G_V^{(k)} &< G_V^{(k+1)} + \quad k\in\{0,\ldots,4\} \\ + P_{\mathrm{GV}}^{(k)} &\le P_{\mathrm{GV}}^{(k+1)} + \quad k\in\{0,\ldots,4\} +\end{aligned} +``` + +Initialization also requires $N_{\mathrm{GV}}^{-1}$ to be single-valued at the +initial operating point. + +### Model Derived Parameters + +The speed lead-lag coefficient and nonlinear gate-to-power curve are: + +```math +\begin{aligned} + k_n &= + \dfrac{T_n}{T_{np}} \\ + N_{\mathrm{GV}}(x) + &= + P_{\mathrm{GV}}^{(0)} + + \sum_{k\in\{0,\ldots,4\}} + \text{linseg}\!\left( + x;\, + G_V^{(k)},\, + G_V^{(k+1)},\, + P_{\mathrm{GV}}^{(k+1)} - P_{\mathrm{GV}}^{(k)} + \right) +\end{aligned} +``` + +CommonMath defines the [linear segment](../../../../CommonMath.md#derived-functions) +helper used by $N_{\mathrm{GV}}$. + +## Model Variables + +### Internal Variables + +#### Differential + +Symbol | Units | Description | Note +------------------------|--------|-------------------------------------|------ +$x_n$ | [p.u.] | Speed lead-lag denominator state | Not circled in Fig. 1; realizes the `Tn`/`Tnp` block +$x_f$ | [p.u.] | Governor error filter output | State 1 in Fig. 1 +$c$ | [p.u.] | Desired-gate position | State 2 in Fig. 1 +$g$ | [p.u.] | Gate position | State 3 in Fig. 1 +$q$ | [p.u.] | Turbine flow | State 4 in Fig. 1 + +#### Algebraic + +Symbol | Units | Description | Note +----------------------------------|----------|-------------------------------------|------ +$\omega_{\mathrm{db}}$ | [p.u.] | Type 1 deadbanded speed deviation | Defined by CommonMath `deadband1` +$e_f$ | [p.u.] | Governor error into the filter | Reference path less conditioned speed and permanent-droop feedback +$f_c$ | [p.u./s] | Desired-gate derivative target | Before rate and position limits +$r_c$ | [p.u./s] | Rate-limited desired-gate derivative target | Limited by $\pm V_{\mathrm{elm}}$ +$P_{\mathrm{GV}}$ | [p.u.] | Nonlinear gate-to-power curve output | $N_{\mathrm{GV}}(g)$ +$H$ | [p.u.] | Turbine head | Implicit water-column head +$P_{\text{m}}$ | [p.u.] | Mechanical power to generator | Read by a machine model + +### External Variables + +#### Differential +None. + +#### Algebraic + +Symbol | Units | Description | Note +--------------------------------|--------|--------------------------------|------ +$\omega$ | [p.u.] | Machine speed deviation | Read from a machine model +$P_{\mathrm{ref}}$ | [p.u.] | Active-power/load reference | HYGOV component base; external setpoint or constant parameter; source label: `Pref` +$P_{\mathrm{aux}}$ | [p.u.] | Auxiliary power input | HYGOV component base; optional, defaults to zero; source label: `Paux` + +## Model Equations + +### Differential Equations + +```math +\begin{aligned} + 0 &= + -\dot{x}_n + + \dfrac{1}{T_{np}} + \left(\omega_{\mathrm{db}} - x_n\right) \\ + 0 &= + -\dot{x}_f + + \dfrac{1}{T_f} + \left(e_f - x_f\right) \\ + 0 &= + -\dot{c} + + \text{antiwindup} + \left(c, r_c;\, G^{\min}, G^{\max}\right) \\ + 0 &= + -\dot{g} + + \dfrac{1}{T_g} + \left(c - g\right) \\ + 0 &= + -\dot{q} + + \dfrac{1}{T_w} + \left(H_{\mathrm{dam}} - H\right) +\end{aligned} +``` + +CommonMath defines the [Anti-Windup](../../../../CommonMath.md#antiwindup) +helper. + +### Algebraic Equations + +```math +\begin{aligned} + 0 &= + -\omega_{\mathrm{db}} + + \text{deadband1} + \left(\omega;\, -D_{\omega}, D_{\omega}\right) \\ + 0 &= + -e_f + + P_{\mathrm{ref}} + + P_{\mathrm{aux}} + - x_n + - k_n\left(\omega_{\mathrm{db}} - x_n\right) + - R_{\mathrm{perm}}c \\ + 0 &= + -f_c + + \dfrac{1}{R_{\mathrm{temp}}} + \left[ + \dfrac{x_f}{T_r} + + \dfrac{e_f - x_f}{T_f} + \right] \\ + 0 &= + -r_c + + \text{clamp} + \left(f_c;\, -V_{\mathrm{elm}}, V_{\mathrm{elm}}\right) \\ + 0 &= + -P_{\mathrm{GV}} + + N_{\mathrm{GV}}(g) \\ + 0 &= + -q^2 + + H P_{\mathrm{GV}}^2 \\ + 0 &= + -\dfrac{S_{\mathrm{sys}}}{T_{\mathrm{rate}}}P_{\text{m}} + + \left[ + A_t H\left(q - q_{\mathrm{NL}}\right) + - D_{\mathrm{turb}}\omega g + \right] +\end{aligned} +``` + +CommonMath defines helper targets and smooth approximations for +[deadband1 and clamp](../../../../CommonMath.md#derived-functions). + +## Initialization + +Initialization is performed by evaluating the steady-state residuals in +dependency order. Let subscript $0$ denote initial values and set all internal +derivatives to zero: + +```math +\begin{aligned} + \omega_0 &= 0 \\ + P_{\mathrm{aux},0} &= 0 \\ + H_0 &= H_{\mathrm{dam}} \\ + q_0 + &= q_{\mathrm{NL}} + + \dfrac{S_{\mathrm{sys}}P_{\text{m},0}} + {T_{\mathrm{rate}}A_tH_0} \\ + P_{\mathrm{GV},0} &= \dfrac{q_0}{\sqrt{H_0}} \\ + g_0 &= N_{\mathrm{GV}}^{-1}\!\left(P_{\mathrm{GV},0}\right) \\ + c_0 &= g_0 \\ + \omega_{\mathrm{db},0} + &= \text{deadband1}\!\left(\omega_0, -D_{\omega}, D_{\omega}\right) \\ + x_{n,0} &= \omega_{\mathrm{db},0} \\ + x_{f,0} &= 0 \\ + e_{f,0} &= x_{f,0} \\ + f_{c,0} &= 0 \\ + r_{c,0} &= 0 \\ + P_{\mathrm{ref},0} + &= e_{f,0} + - P_{\mathrm{aux},0} + + x_{n,0} + + k_n\left(\omega_{\mathrm{db},0} - x_{n,0}\right) + + R_{\mathrm{perm}}c_0 +\end{aligned} +``` + +This closed-form start requires $G^{\min} \le c_0 \le G^{\max}$ and a +single-valued inverse $N_{\mathrm{GV}}^{-1}(P_{\mathrm{GV},0})$. + +## Monitorable Output + +Output | Units | Description | Note +---------------|--------|-------------------------------------|------ +`pmech` | [p.u.] | Mechanical-power output | $P_{\text{m}}$ +`filter` | [p.u.] | Governor error filter output | State 1 +`desiredgate` | [p.u.] | Desired-gate position | State 2 +`gate` | [p.u.] | Gate position | State 3 +`flow` | [p.u.] | Turbine flow | State 4 +`head` | [p.u.] | Turbine head | Implicit water-column head diff --git a/GridKit/Model/PhasorDynamics/Governor/README.md b/GridKit/Model/PhasorDynamics/Governor/README.md index fc46e5eb3..6d3983adb 100644 --- a/GridKit/Model/PhasorDynamics/Governor/README.md +++ b/GridKit/Model/PhasorDynamics/Governor/README.md @@ -9,5 +9,6 @@ A governor models the control system that regulates the output power of a machin There are a few standard Governor models - Turbine Governor (See [TGOV1](Tgov1/README.md)) +- Hydro Turbine Governor (See [HYGOV](HYGOV/README.md)) - IEEE Type G1 Turbine Governor (See [IEEEG1](IEEEG1/README.md)) - General Governor (See [GGOV1](GGOV1/README.md)) diff --git a/GridKit/Model/PhasorDynamics/INPUT_FORMAT.md b/GridKit/Model/PhasorDynamics/INPUT_FORMAT.md index 5bff41ad4..8abd23e7a 100644 --- a/GridKit/Model/PhasorDynamics/INPUT_FORMAT.md +++ b/GridKit/Model/PhasorDynamics/INPUT_FORMAT.md @@ -145,6 +145,7 @@ are specified: `Gensal` | 5th order salient-pole machine model | `bus`, `pmech`\*, `speed`\*, `efd`\* | `p0`, `q0`, `H`, `D`, `Ra`, `Tdop`, `Tdopp`, `Tqopp`, `Xd`, `Xdp`, `Xdpp`, `Xq`, `Xl`, `S10`, `S12`, `mva` | `ir`, `ii`, `p`, `q`, `delta`, `omega`, `speed`, `Eqp`, `psidp`, `psiqpp`, `psidpp`, `vd`, `vq`, `te`, `id`, `iq` `GenClassical` | the classical machine model | `bus`, `pmech`\*, `speed`\*, `efd`\* | `p0`, `q0`, `H`, `D`, `Ra`, `Xdp`, `mva` | `ir`, `ii`, `p`, `q`, `delta`, `omega` `Tgov1 ` | the TGOV1 governor model | `pmech`, `speed` | `Trate`, `R`, `T1`, `T2`, `T3`, `Pvmax`, `Pvmin`, `Dt` | `none` + `Hygov` | the HYGOV hydro turbine-governor model | `pmech`, `speed`\*, `pref`\*, `paux`\* | `Trate`, `Rperm`, `Rtemp`, `Tr`, `Tf`, `Tg`, `Velm`, `Gmax`, `Gmin`, `Tw`, `At`, `Dturb`, `Qnl`, `Tn`, `Tnp`, `db1`, `db2`, `Hdam`, `Gv0`-`Gv5`, `Pgv0`-`Pgv5` | `pmech`, `filter`, `desiredgate`, `gate`, `flow`, `head` `Ieeet1` | the IEEET1 exciter model | `bus`, `speed`, `efd`, `vs`\* | `Tr`, `Ka`, `Ta`, `Ke`, `Te`, `Kf`, `Tf`, `Vrmin`, `Vrmax`, `E1`, `E2`, `Se1`, `Se2`, `Ispdlim` | `efd`, `ksat` `SexsPti` | the SEXS-PTI simplified exciter model | `bus`, `efd`, `vs`\* | `Ta`, `Tb`, `Te`, `K`, `Efdmax`, `Efdmin` | `efd` `Ieeest` | the IEEEST stabilizer model | `input`, `output` | `A1`, `A2`, `A3`, `A4`, `A5`, `A6`, `T1`, `T2`, `T3`, `T4`, `T5`, `T6`, `Ks`, `Lsmin`, `Lsmax`, `Vcl`, `Vcu`, `Tdelay` | `vss` diff --git a/GridKit/Model/PhasorDynamics/SystemModelData.hpp b/GridKit/Model/PhasorDynamics/SystemModelData.hpp index 8a6e73056..32c1a5c5c 100644 --- a/GridKit/Model/PhasorDynamics/SystemModelData.hpp +++ b/GridKit/Model/PhasorDynamics/SystemModelData.hpp @@ -12,6 +12,7 @@ #include #include #include +#include #include #include #include @@ -42,6 +43,7 @@ namespace GridKit using BusToSignalAdapterDataT = BusToSignalAdapterData; using BusFaultDataT = BusFaultData; using Tgov1DataT = Governor::Tgov1Data; + using HygovDataT = Governor::HygovData; using Ieeet1DataT = Exciter::Ieeet1Data; using SexsPtiDataT = Exciter::SexsPtiData; using IeeestDataT = Stabilizer::IeeestData; @@ -106,6 +108,21 @@ namespace GridKit std::vector stabilizer; ///< Stabilizers within the model std::vector constant_source; ///< Constant signal sources within the model std::vector signal; ///< Signal nodes + std::vector bus; ///< Buses within the model + std::vector adapter; ///< bus-to-signal adapters within the model + std::vector branch; ///< Branches within the model + std::vector bus_fault; ///< Bus faults within the model + std::vector genrou; ///< GENROU instances within the model + std::vector gensal; ///< GENSAL instances within the model + std::vector genclassical; ///< Classical generator instances within the model + std::vector loadz; ///< LoadZ instances within the model + std::vector loadzip; ///< LoadZIP instances within the model + std::vector gov; ///< Governors within the model + std::vector hygov; ///< HYGOV governors within the model + std::vector exciter; ///< Exciters within the model + std::vector sexspti; ///< SEXS-PTI exciters within the model + std::vector stabilizer; ///< Stabilizers within the model + std::vector signal; ///< Signal nodes /// Monitor sink specs std::vector monitor_sink; diff --git a/GridKit/Model/PhasorDynamics/SystemModelDataJSONParser.hpp b/GridKit/Model/PhasorDynamics/SystemModelDataJSONParser.hpp index a3682feca..3227aa95c 100644 --- a/GridKit/Model/PhasorDynamics/SystemModelDataJSONParser.hpp +++ b/GridKit/Model/PhasorDynamics/SystemModelDataJSONParser.hpp @@ -142,6 +142,12 @@ namespace GridKit raw_component.get_to(gov); sm.gov.push_back(gov); } + else if (kind == "Hygov") + { + typename SystemModelData::HygovDataT gov; + raw_component.get_to(gov); + sm.hygov.push_back(gov); + } else if (kind == "Ieeet1") { typename SystemModelData::Ieeet1DataT exciter; diff --git a/GridKit/Model/PhasorDynamics/SystemModelImpl.hpp b/GridKit/Model/PhasorDynamics/SystemModelImpl.hpp index 287b328dc..9367531fb 100644 --- a/GridKit/Model/PhasorDynamics/SystemModelImpl.hpp +++ b/GridKit/Model/PhasorDynamics/SystemModelImpl.hpp @@ -259,6 +259,42 @@ namespace GridKit addComponent(gov); } + // Add HYGOV governors + for (const auto& govdata : data.hygov) + { + auto* gov = new Hygov(govdata); + + if (govdata.ports.contains(HygovPorts::speed)) + { + IdxT speed = govdata.ports.at(HygovPorts::speed); + constexpr auto OMEGA = HygovExternalVariables::OMEGA; + gov->getSignals().template attachSignalNode(getSignal(speed)); + } + + if (govdata.ports.contains(HygovPorts::pmech)) + { + IdxT pmech = govdata.ports.at(HygovPorts::pmech); + constexpr auto PMECH = HygovInternalVariables::PMECH; + gov->getSignals().template assignSignalNode(getSignal(pmech)); + } + + if (govdata.ports.contains(HygovPorts::pref)) + { + IdxT pref = govdata.ports.at(HygovPorts::pref); + constexpr auto PREF = HygovExternalVariables::PREF; + gov->getSignals().template attachSignalNode(getSignal(pref)); + } + + if (govdata.ports.contains(HygovPorts::paux)) + { + IdxT paux = govdata.ports.at(HygovPorts::paux); + constexpr auto PAUX = HygovExternalVariables::PAUX; + gov->getSignals().template attachSignalNode(getSignal(paux)); + } + + addComponent(gov); + } + for (const auto& excitedata : data.exciter) { IdxT bus_index = 0; diff --git a/docs/Figures/PhasorDynamics/HYGOV_diagram.png b/docs/Figures/PhasorDynamics/HYGOV_diagram.png new file mode 100755 index 0000000000000000000000000000000000000000..fc7d7916d2887f6758718be182a29573e7599993 GIT binary patch literal 72381 zcmbrmWmHws8ZNwPqyzy$I;2}Vq(Qor+DJ=xhe%4dq;v|1G@DLAy1N?zkr3(p=JuR3 zzA^64dl@LOS$nM+Z#?hwypwP>l~>ppBp46~1Y2HCS_1+>CWAl_Swxb%NMrrux){P*?vMXt5~y=>p>kMzH*$&qmWd-0b0|Lw((@)#pd zWNEC)N0N62(?t~G_<3qZ@YGM}Y2Y6xxTs`d(5g&^+5EgOZ&rSD1<*{j>Y@7O3!fFL7I6-fIX zZ-t%=ubNGS2DuXf13iz?kMZAXz$611Pk#9q!xp@pUou-OY!Blw5zr^@?5uaU*=}hg zBOuU7=i^QDnv3H_RpM|pyfY%pYYwK#SwhFRri3C_ZJT8L|Fg<_YRWFm;c&tf!I!zE z43^lP*pbTm$6IUx^gw_A9BTsG#jsVPvJQu4Os5v!w3hfVXz1%g!D24$!F2K` zGAMzI9#Zd=7(Qxeb&u3!=8`#K@~Fmiz*2(quaAZ=H)nkaWL%|glj3)&AQ%z|64c_k zxH-5+nLi5E7^&t6d~kF9C3D67Hy)6l1X79THq6=Iz^Kr1krA+?lR}U~sH>H~wR{dM z{*kIU+L!JVk*BZig^!09pPAgeEGI20!^48RR3&V3eJ)#Jpg&%v=i%D^`^8(QGjS|Z zL4thQAJhfj9Cda}v9yhdRE)Se5Gq7-JX{onKu0bkBLW0;2qJQrnmDo%s!<383=KjR zfhdhE{t*!u5dxuSGC5@8F>bovUp` zh?0-1P*b7rIIdN+GUDcqyYx}LojJ#fskkE)`*4Bwg4@IR=DtSY89_j~L@?*Kk?$O( zb(@X|V|-07B^7I`?%u!jh1R1LMIhlV?IAk*>4ox^-hWRXMswUyqW=03gY3=RyB0)xJ$k>_l+pICX0|N8zs0M`DXp*QwVP@RtVeGTcMj=Re zUZfV!aCAxQP40TiE{!^O=3p08GY(4;8Tej|%uGLa6n(`f#w1QcYtg|bK51==sunZu_=;{wCN#L^)Cd*WP$K7vGgeH=Qgs zpAho$fNEZY4ub1ppnc|IzH*K~L_MHm4c#?jFRU*H4;>rJzZhNe=2GCT(?e^+QC-`v zqY>fLQNw_jnH61{v>}@d4Qd}zd0HHwlFG-}#;05{(}^aE4ivjsP%R0%KW-l3gAyzazTrk{7=HO|O3H>5{aF<3rLeB_VP! zZa~x{Q_w>LyDkG<-$pHqdCN6*=gVAsp z{~HZY7b~-oz2776=mwMvF)uNmN&2;w@Y(gR4SeXOB@&PUHD6nn8Frlz2Vaq>d>I6E z3=GIeadZqkTVhO0)vo|nFbuozhN)9wG(ue@5W4g zK#R8%h{Uv$JV>5?PIKZW$Z*1G)}z$&V$**v#3o`lFQgu@>cEms$T#3Tas-sy-&WAx4 z^rk}`kyCM^#Jft8Y!3FW=^Bdexh-+B4`(dz$T51&Dn!Ed-XnG1Xsy-o&zkmr)_MDI zYNw{Z`b4uEGG4@*gYj;W-hx61{ zlyo(p^=~bW^B>V*o&J_dAoR6_aXDn4J=Cmb#_Ar8bOm}Cb``w$p)JNGo_w_e;{3bN zXE@SZPb8@Gnzlq#DEoWNOVdM%`f<=1gNlB1E)~z;rrvr5cjW-gdfFtDCR=A@rH5M+yaC4sdtAyaJ%Sx zWUcg>3mIIxV_evGq%{ip{Ov->5UnQl-4X8dlDlD`Rr zxJ;X)cLEd=|IjBG&)|RvX%R7HHezDF8hA|HW`jaG*2}cEy zUaaXj;dYDT<#9pQ$5V1tH4MA!F~!k`IFZnwFsr*Wqxn^ceauJME{?K{tTPUG z8Y*-MGK#eN$C4f_ZUS^QWGXdkNo1p^2rmB5Flrz*{pjcc=vZpvoTR`-BnS|3K_v5+ z{HV@}B+DH|&yAdZ1SJ9++_yM~O6IJxwchid2jcIxYE7Jv#B-Y>yg#6|4>8!E&#E_0W7(@-6JP1&+U>D_}cGZq)=R(M3um@ zg2KXzMd9|eFO8Qsk5$u;m0-r5&SB0tWgdOzV)xwLvMslyGd z|DCg!8-_!b_nuJbVvKTge6}7HrHXXf*dX^8ljV9SqiH@D2x?#qqF_z*hfFTr1+f|# z$jgm#qoPBoFW!|&ej3o}MF^BxW$$HCRW9SDNQQJR&Pnu~a3r(z(iY=s+XL*@SbAzej_B4VWtVg5 zQ)`?H&0{cMf$Bfnp=aa;L&%F!!A#udLf8abI*3FAbC#^v?u{*^@WnFwE^tcbT{F2u zMML9^0>}Ep&zbq2s8_{9{Pm97@bbYw?$ysRkBQn5AUenQ+H+M_^}p4nr=WD$1H9iB8(VxezoN*+!Et_AHlYFaGW<>;o9yZ)0)G;F?qNOk9>>>Ur zL(WurvT|XOoxr-_DPYKZ^ol3!o3b94)pv5W9$CGIBnv^1FeBEb(2T*Nf)RU_`GU-# z4KGu+P|Dlt#2|1eMWV!f%xEU~T;C=e+qm!LBb{nekVHa(UF`E&eHUnEMJp+n0$2h5 zEO;#bcJWjA+Xbj}-cSP54{Pc8RcM2tScV4f%y#WA_$dUZZpJH>XMKYJ3j&Xdz)}DC zMl9Q6JZx{n6;n^Jks%(;*1TPq^2qnwk372bc`UHnawJ0N`)sdy!@|NyBym5KO$p;Q z4*xHGofAvgHWFJjkE!e~eh6!sH$8u9EBPn5=ENMng(ohUm9AZg(m!Nr#OfiKs}C5Z zDdkdcY@a9F>rna0<^LbarKP_Jx#u*VAHdvskh_Zsg4V3I&W{#6T(@tM}*g7S!m0w~HXpG{URQEut2- zT{XiFF%bz zUtaD~kcLky&KeJ#-i*JjBXzf3h<@FO0)KVfV>Jx^oAnz* z_vfrVQO|{T2^l<%0#2yTG6w-GfBT(()hDPZ+n1zNtL%G{wgJ--D5g)mhlXis+VkHq$7g0@ z2X*>yZM`GLLcPMRphJ0FVAS=ia6c~l{lcep$EN;C&agg@Z~DzoXs{|#9&r|pw zr1%K#8nJ)~Ak25rqd9wCZ-3Q@a$+pKWpRGA8tb+{v))7MkmvzS?+Qz$#4?zfnK|zB zXC?>w$9oGf7SzWStPWqyO@vz|OQ$6Zt+T99qcb%SAhRqEHF|t53m+bON|q+U%n&jL z*gCfMF7$c@^$r3-qVQP^rSwHF(~3-D(WVg+g06s)cnhD{@!u1Gg@6H~dl5D2dU*-* zIjju%OJUlC=SZ0%;b=I?F9v=BgPh5W*Qk_$KEa2SH`V!-3E{xauFh$*Cs(_|!1Lyy zep?^xt!NNJsFoatz_H|9&86u3l%F_4JbY_w?C^`rzkOqkdd7{an8N;}SUrDh`QP=t zz|J)*%l2}{Epuk~x1~VL_}=8U?|64MYr54Jf5HxeQucRXF)$J=0I2F&SXpb!Rq%+3 z8UJkw6S%W_>GvE^)6{_F>onqYbi)8!h52oL45O+}{BPXK&a9yy`OEj`wRNwgrHxK^ z$K6d!u6R9^Lqo>)&*W58LbOVQp=HS)`0{$wSqxwW=ASuh^oK$qr`|7Lym*l_ZN04@oj1E5pOFz&Tgy}M<;!;?&M6FUOP)07 zs%oKzM9G^4*j2=IskTeuW14pSWs~3+F|_t<0G6^!zRnswhv(allhyYRf>>UUwdJp2Dur6 z%?i{gIi3Gd%D&Q+LNBM}N|zkPQqDZtK5$gj!Hk^nS~ueBAPQ(c!az@9oBNXK5%Fk~EdrkpQwp zNB>gffi;W-JkyuZ^LWdti#WgvQ=7!(IkPE?pNVN-H)V z5^TuS03Yovn39sJ>MJOeB+q1n>+f)^B@jxOdFv%F$GoCP%89Gw*v>yC1v@lK}+ ztGS~k?AHTC%08)rryaj6HZJ_aBhba!Jp}dN@~>alltNS zh4sdc-#W@1OI5>Z^|>_ymVjDmj}cO1gaq0N1M~C`tTZkU6TgIjJ(QOt^UM`qnk}v$ep5OcJGbT3t&hR&W)4w3qm?A+XzfDvr3X;ut zx>Bi@Zpvpppj9x-^Fdb+vf?56HafAn>M->(IFEgz&Fcqm*YR2b0Y!6I8lH_ntz=et z=-Q)0-{E0ATgz|~9F4C=iX*v6TDFeC!aG5w)7HNi>mOe%xy0(9sHwjIN-JUwRmt-@tnfN47J^Uc*)S;M8y>ewjCaY`m8 z>s{A--OY}48-d+3ZxE|19)^k?vJ=mHnN(Zm)`n1DfxB(~L3>3@|5dooZkYz;pT~eq zg9*W>W3Tbc2JpeHL}gjYAk4HH#1<4hFDNLWVq$tEl!>EI_yq2KopB*=85_S_^5nI& zEZP4>_OK9QOd9Z3#n68yb(lN5p9_-UxaTank2*nmt>6Bc_pYmulGnx$dK6Nl*My=k zUzSb3{dN<19h1=C-{0VUTnYch{44V|F3M-mgeqF7E-Xl*QNiy-Yo;9aZ{Ar$Ny;fi zhS)k$b=q*)i(ZO@{7=yVHlO*P*w%j=CxAb5%*?_yUI*fkME7h*%(Y zApV@mJ%DO#^xrL_J}F2*@!|`rAByxqvyecuDsQ>HD~YUOMY9M%r_nYuGi&6DHUY*8 zk8WGTD!{F2sW@?)J&t-J2wAq4x1`=SHY8Ykz-_gn^*B{Aprkm}v%BY_gVOIde0wlu zusJ+qV!!6s5+y>;d)EwDs@)qUkvfNaebMgRQq;#?jWrfVoWwdGKlgQkdZfMADCj*M z1A|n>tX`4ga|VXNP#yV6rd#UtT=tV3DQD;Ev#n(#&IU49$HeW9->A-^!BBG9RbKU# zw~;I!_1ma8&tj@iFz4F+#M}>;06eD59sT;SVSjdEXk2Z!9Z@`o%<9S|;5Hpcs`geR z5TIZqkZ1rH35_wcpk?Nh00s~>>~!#pntyGoAK(3MWDeCRmP~M>IGsa!iv*RduYzlC z%|Jso$c`bZODc%b zW=}M*M7+uC1J|69Y8t9>vpi96%mR5l`^UoXCT`0gQnHG8^*M!0 zXQ>{GPyOH<%q_o2>C>E8{i}Jl)N0SVM*}M^Mn^$GfoJ;1)P__MX9i32=Hk*nx|?~F z(fB+PA|oQy3=BqS^0z*H{D_?L;webe{#K7;Bz)~-26@RKYp%XoXAdJI@574s;G5NL z6lAI~jB1r2YP5YY0*XgV7!*V(0G<%tlv3ipysk~gG%qm^qK&Dh_(fwTSo_5Y+bgDs zJE|bEABS8BlB{Wc8!%G6R9P>ls#?!J3tgR=i!8c+L1tTaf6tQgWluc04M!0}Sp5UY z+&DxwFP6M_korc2W2RUA=#HG_21epz^ejd);-jf(JODDh#R;8RvsOjHyW7BMH}(W6J$7;kW=>BMHe zIIEbp7%q!cv+HgaW_Qo7=?xp2(?PfgdH>U=PrqX(fsZpWF^yTnX85cY^rU_;FjjT_ zV0r>^WIk#b+}l5i6fVaQFDH)w=FDWe{;RW#6J${053|xKkt+(lQdZeK>B9|<>w&fe zc#TL%I9u8Px>Ol;A|j)qanG%pPw-&q9CnnxVJ;%-Dzg&;HGr(UJ0Hlm8Rk9$-W3Ho z(NySiN=niIXuBpSu|-5g#DRH$nhz{H^d_8yR#YN9Q|}>;zgf9xB(v;+gNdnX|1U3d z5(dar>O=>}YS!+uvZ&QHHD!ueEJ+FgQUKTrjTKZ-)X(gp@F2u!@ZlE(XsM}LAQKTM zP=9A6v*sqPFCL@~RJKkbvjt38wN=J^_*Gh7GVFgH#9nMa(B`y}tk-{geegf|DH!UKO z1D|hOHrZAG&ey5OrWV_~+Zh`hyBp%izX!#QuS|prW+2tki4EB1XV0FgmAsjm^3hh| z0~j|hPVFH>->S)9Hos~mSi9s6OvU}o`z1hH#sCNpn$Pf0cO=d3*<_$lf@%S5i+!0N zH46T(Ch@pKbA?e0QqSYwixE<}&kU8wqhw}9muVmy+pgX=F9zYjL8R$Er8`j4?Pp4?oN7s@d-<}rdoN+^?l zQ*p|gcLP81Y=tRNWl2AC@L>TnHurVWi2KcpBNP6pH{6m!{yC^wl;w5x>yhHS&%f_XjpP6Q^MF#?S>34lcM$2oV!+2vIge*`{W z?Er~RK}O)q_dSRb6>(?fje)gwX(BeCl(Fii$;S zL=*=6y_)pZR+L5wywKS56}Yr7Y2Uz2 zyAWnBa7kh{fWqp*k&&r$6M1^-3P5H33R^w4Q1P^N3cF>oTf$)u$K^40sMRq@KEy7w z()(Jo{mrh_UelV|(?U*1)1@iPv-b)`lhe=`jfz=e+X|Z`5~hy0w7ksNhp=>ke56+{8Z^R$k z>-=82m+_xo^p|r836Z22c$est;x#_580f3u?CQMPa~0Z8d<^i>MKE1qUR>R0?1^9a z(eqmYC?%;&*!**=X7tX_K^%GXNI`_lDyj=pLECn<8sRm=N(4p9X*apYlJV?pvxkFg zoo?0lY+xk&r;Jt!TH$f^vxH=7P_2|W_VfKT7Y;K;5EA>Fm+0x}P+oGVZ=4Yo@p)6m zA3k>gcdK9JnlCa{s{Wy%B(t>b=1J^>X>FG++N%Pg#~S8c1!xb{lpA%V z2+~``@Ccu}zDBwgJaYS)hkx<>Be?Mo=$K+!WSK?@nf-)_oz*vyH7XQ)z*M~S#*@77 zGV1H;!|P8yUn6Esy)h3ue9E4BL;vNa04km_HIC|1D7H%$p#*GRB|q*(&M)R#mS~T{ zIR{gwB!zdouBk6{qk)-EZDR_DyD|vw%t_cd&Og1tvt|Vp@9*EXZ@G_|S;77e`QGD! zyVGlZxF}q9QcCl!GCE8r7RwXkD|;YhPWq~4j>LLnYSTYJW6jh8eBwJ6+L6V3L}q2Q zV{M^9YielOOu2qjLmaRlkQj6@d*6Ubq5d z=!#vn=DKP0Gr38|?(?W#xr6)juVxM|w~aEL<8@LHYEyakB(6B4h`GL3{Nm>+xR*4Q z3wEKa>y>$hl}s2mhUM`tNr8^uX3BGGoeBOHFkrvdd7KHE6=3r@roAejCJUReKy1LP zb;D*2Y-VMra6yqHR(^@F>U+a>f{pfVlC>vxH&!^|daY@2!Spul9)6$8!$H*?jgsT; z+FtVeVJtEMT>cf$e4y(tDd86CeV?)h{0kt>5Yzb`de%CFwwDe1r)M*=vtxI6Eq|}J zQ?|apNcP}o+2a+vJ??U6SO>h+_ucXQgC3!LtN9 zP@2FrTJ{yp&TG57rT1IBUV%JPru&FFX`ty-F!$cxo}8SVqo;%kNG+=zFJF@X`QxVj z>zuDtyQ1ggBYXhV*pI%>Zv2ScRyH_lyE#nYpO!D0T%Whd{!wc+`Fo*|i8mVNabzHJ zIwmFHv&A@)`4S1zV13Lf-Ag|Om38UxM*s_sqn%4E;7HXI15?esK#7T|0z8ARuI`L@ z5YNHj_-bZGV(`|0&SD8ZKuBuLNqqwYU)@^bVy@uV@L5#Idk=%V5G7uD040D zi{CgLipDJ=Ia+FsxjqNPGh>_D`Xn`~3U&6&a6T=+xVW!LzL(*ujl-i0o8PEA3=dG6 zCx8`)QfqUpw1$^fTy6CygoTH91|Xq1&YPw5H>YcUKnOiOJ@pn-BhSPt0>seom~{S6 zKlx^hm>AfwLCyx$+xo;`j6i?j7)1=Cs6n0fNsnMq{hjbq*YWk=?->{vHoL>|J%MQP zr;j$auUStCz%GCSf8D<1qg&&&nPE*#O6pHxhV3if zw(%-Ug=OVOP6x2rjbuIhc>a2UXxX5tsVO1D=6q1|e~M9_aE*+Lih@cE^-wSX1b+Od zm(uj_&tbsZsv#riH8c?CG3I?}u7h)QWl4;-GOu2l0@_aH#xAJhQ$vFQcuw2MYBob- z^+N;88k~9PPj7)r#TvR<{T5{sbo)gB$V(gw5ym5DMs{`_V5zdMuCAHwdNdw@lJXeb zfmN7vZV#Vy{%KM`qQzfsR{tok- z?vCZ2P_s$qkBb*gN?38L>X$PigwAv%AfWYEwBJWWlkyLwXZ{+?MUsso4!FO&8EkF| zq(Yw&-X38myE`D(Gt|<29g|>Q^9=3+hJJT+M_(s%B!W1trx=T1Ng-Iz0^1+%Ni{V! zPuN|5O+6Y+W^?fzI_A+nr}RJ}gVuccoltl3_ivTMswjY+1Ox;*O-&>aTwGlLa{Xo+ zQBex;&QxxTm@|P{n2k120jRjSC)bZs`t%Dn<@I-N?}z@3XZO95C-PV^=#mf`BW9iM zNk~Z0)BhP*k>KkQ9v+TUI$&XDCUGDMg=v)3He;4d?q#TaE;&L{D0!2YlLL{Fk(sNt zN@!{RA>Ww=^`X^FmBK}?aU zlxU+tzFnhAr5ZqLQ+;U9%cDt5OkC>-Ksuo|2nawsnl4l%T5=v_OBMYfbbnUix88|H z*$v{DX76&NQur}V=o zA&D3pQ|ovVl~AotCI~dS!}bSXFYfqPwLP1Uir&;pmaC{bnoob#bobyZco4RBDQpLh z7f6-!;pSYWS4&gVc&b?aayLKL@oGCom<$O5%EKPoSdl=uz&SoX4!gE{I4T~Ve*mJw zJ3E)88YaTR!a@O$!?&4=+45Ey8hA>6@)u{%Us-+5V_YDFx}ZtD)Z8yj6wFi!otk=+ zQcd=}?r&;R@^(83(AF?2Gi+gYi%)^dnf%zwy}ZRDR`i?d)lR%sQ8L#}%hKgySvk6KILK6;KG zXvFq)&#G>!3Hn=dCM&g(%<|mr<{!p;FaS!{dy=B;vhlYx(W8tS|R7;tla;6X=h+YZ& zvxB-F#``m6Ixy>`Wt~?bewdh;IPMju)TJe?e{p>3*w@!r>8n30Qy~nHi|g(f!o&4! zI|+fKKkF_K7w6X+2jk=8896v!IXZG8A|j5}+38*33J;%UO<8PohhN^`96fV%lmscS zs+Xda@GX2Qh)75W3%@qm+W#igi}lP95&uO4RU8Nvyg>FkSjZ%Fu|@auJKBN|#?*AY z#zGa9EDh)?iyu_#3xI*R9VkmN(TMy-yZKumDBky zkGmd|x{26gSX2zVKm2(!F6So_A~;Ewy*%QZ;_R{dfj^m5m&SqK`h~?!&cUDVY%VZA zFTDVqhn*%684%E3#qafGgpshqu-`Ci22ofj;vy?$i!OWjo&brh!cDuS zQ4m9*PkYnxRR_H4fBCXiJs_(GmsPXK4j93Ho5Ad-V(cIOQeBOk#QP3u`vQAkL5ZXr=0T#apg>Xi^Xno z1)D0bzvf;w zz|hNb*at{Y>i_W(i$d7qrPDOr-1tD1jBYbM8cs|$i3J!J=F7L3z^5(G;1~UJJ)ZrS z-DJKpeZp{9GQ$A__WZ0qpv}{H=?p-oE7PJC7G8wf(DL$5H;Hx7U&^wTR`u$fw18s% z^k|h*LPFw%y%c0X=*rZ$G@A#ja?=#ReD->`e2R!x$eslchKWx4%;(37LTyOmGv;mW zlk-Zans{gbbT8Jfl!}~?coOh4F>lg^ZY#Xrd<>dy0u6Iz`X;IsV%GWy7+a9@A%$mX=#-q7QqP!!+HW$6q4tl z5RhQ#E5A4RMS~JjU=&)*A;&;4RSn$_&*T8LxlWU7iQ6x2OH0PUpdbVYsB@#gnr^ZW zw>|(e@Ht@V(=)Aq?RZK`N|s${VSM5{g?S1mDFg$UIHP+cTc$PZEf%&^*ci*O5`utc z?Q}Te7?=IKbo!T2ij*}Im^YjASbtxiqfG#PP^a2{OiQOxvD=>_jgs4quD696)6Xp9 z7};2&8=N&5*@E^IDuGD}x9S*x=(l{D2|99wR4daWwMy;cKeV#fDW>v;9Ne?r#n^F^ z+|o7t9E4AzLLK_007wj*X*T5&Zf?A7_t!RCKYz}=&p3PTu;Pb>DwPmxTz6y#Uzkzb ztJyZCk~bTR4(&hwj*7(s3R6dnsL*A!wzd|!!zF#VHQbxz1H_?8SkoLhKi1^a3w*K38{BFf!t!WXmT$uJv|Zu z2~}rtOOa8dOa9qwAf14OjuZ4UMrB5C^Im0&06=a~;2^WQ9O}`klJOTnu_v6AOm13% z{F+2c@ z3&xIGYz)-K+-iaiMR|Fk+q-y6(wM1Uet++Kye6^aHmiST4~JEuP0|?7f|yT0%?J8$ zQYJ$ROC(}eZc+?=8)sa)dwUl@C}XDpG7@Aozz^-ZHbIZLSleHD`A5sqEJwSWB_P&G z9|E=?u6To{4?uo}KJ%Cm-PxV{8T%%A&qKe-RRWBZ+&Ak<{QJfu=GbDqA#?7v9B&&b zd5_gX;kdeS>r54kR;g^D`zB4Yw0w4bHOuExP;L}xP*G8{nhtUuO1AH{RjZz#$w6Qhm8e_5vleS9);do{(TuhD;+51WVxLW%xY?t={us&nnzA)cC1QmF zDt13uUs5vA4Fkm%<2k7ICG1KV9o{3XAHi>okAJ`ivOe5q4o)yl!Y8Vgk_m>?XR0Cs zU5mQ456^*cF<)556ihU24THsm)vi8c9W?w-59oY00IjcY;tfv~+EV~PT!8iV^ejxf zSRR1N@!;G8!0!uKjd%gbAI{D%fW>vy>g{<9B)B=t0y(Ah9&pN_(E}m5-E8mjYxPt^ z$$nrQv>u2;(Cqn>QFlQmT&VWXToIPUu)bvyzob`NsLI2RaB&wXTW3<5_nn{P z4~|tBiZ$2GsI;{Y@(V5w1vjb8ow3h2IG*a9*m)~Ts9ObnGYIPDgVK0a^sdR*oNjmeazr^A$1%5|3WXrlB7%kbHhDF{7EGB* zp0KZAR}ssD%>sr7q~0;TUy`_jCZtAk%F1pHstS|uw|e<#XgomEDfG?j7i0h!VI}kU zwHC4mX_R~%sSXT0g%Nnd)tPi(g!IUFhGy>(*LUCLF}Ws7hEBYHj^ebn*Qt^1iZGmK zF~Yx9 z0dGM1$GDXd#kzfH2Z{{ntKs{Fn9!2#+5!R@Xos4q^{@!A+w=UjD*z?3602<$PxQ-* zKPo+{$)P^=ZaupKQoJ>JP=O2Wt87|n>BFhxQV|;d>y1c`Ng=1DVdE(M; zhwCeF0}Bg1q{m=A097>6pU&~=Fi(C^UtS{RT8?mhpLRJulx~ummbOdCVJev*RSw2V zo*Eaoc3dD&VnLeWuRJaV08q%x8mQEJe@&N;xBL0|fW|PmtGsZ>o2^@1f`Ym-QCOM4kW+~g|;+ma_USB+dWYH+898|K3GliX&HI(K{3sGYu7>5 zqNHlQIQBMZr25U*gOkA%(}IHXHO&u+3WaWFn0q6>djh%B*80~XDYKCHNoyFO zW+c76o7|YQ+4fxY29o0wVoD|~`0>_m4X!fhcI=eDId<6y>!#k52v-`2XX$PIN(BNn0O#1Cacv=nd-^ z<#h8|cEaa{Jyru;T1}TJ5r1bGUcOw0R@QTUel&6aoG`*}$Reell7ScAbdlM$R>xm` zcYDGQ`VWAz1VUr=+1~|xsCD-5$)d*u7%1nz=V?L58^0-Nng|O)ChZdg`K|a3Ga$s9 zc~thwU2UemD3dPvJ=}S&1)wTDdo%Ti#-Me~eXHw9RN+A)S2=0OwH$mUXVpA8iTc`#s5Y<+S!ny@ql5SqCee_?v! zgWh=%)ZDT7YDI+<^|1>JdqNUBjc)t|bTkWK#5|U`=0m9gi}uafU?luTUK~IP*pdPV z9=uO?1Gr3k*7l0Bdn%2)wxs$Zc}_Pz`GmUOI&YhNt5QVWXFCght}v*;?}kDXo(9MW z8U_Y5Z$tF)r$V%5*# zt_iRurae)eb5*8X#@~@e0MY3Qs=4E}pvN0<^kig|v%J~gEa&m1|3U7V-q~<|eI`;x z^|zeFtgJ&$VTrGQsO}7}wLv3CO6Sy`YhQoA)mg8?iO8!PMqr~ro*3>~48%qaLXel2 zzg4qC)ga-tO?G`<@&~ep5wL4LWiHbuOCi>x;N=T?loux_b zFE=q(*$cfQSllnZDUTtI(rFMl;ZmYK+|R4=LGA_V-?)moW{A5z1{+ z4!SWTlQN@y?Q)>ko(f4_F|PV?;ayzzs5Ja|ZCtOHX=MP^0QAOjE7=Mt*8HrFXnUm) z#Nez6@WsWm6AYlof{D&3k)o=E7ncAvJ=UZ$AR~Sv=(DkNeEx6LM>S+Vm=`wIjF~FE zn8GI{6K|*v6Qs95e`{PX=bZQ!iN!qj6LUe_D~_N}h&sXcFz{{P00guXbX3x?x7yn7 zK_7K-r0>J$IYGbZD7b$;7*Me>P+mcou-&qeyw&i?2z6`of|R(OFqG3dJG+{@;eh_pEbkhqFs$1+&0_9@Em(_biB7j!{KX zyj-yV#tmWzptk^jk>n9S1OjH2%`E$pQ!lyuL7UEvI4&}9%XjAH(8aG`zpAOJ8S#&N z1`McI!mM7aH!r=)~>0cZyJGTy^06 z3hzGq3|elB5oy5MY!4eaoc)?k<#)(xc?McC7+E%#6#GA1EL&XgBH<`cF7?i(rICP2 zBoAcD-MK0=72nn?fxSpjVFwtg9*Et4r8}1g+(m3DLTG!vZAVg2Xn*>02tsIbm%kLt zR=7XX^v@XoAlf`elbLV2C+5!ICmGnX))?2uplte01n(~nllT+B!H>h8e$(IQ2aA7p zvi&vyqjI6+-^SIlEwE>UfMx;jAXJ_b&st%4=4dHej!fozb-Ee~bex`MCx%=&gyN;QrzU7DFD_&Q7Z#nwY#yv~0S6*wWYO$|W?x4*j{_Pl7{8NC@T$o*}uibq23y-Z1NSA&7fRzKOpttPH-&0X=pW&i%5}bCQ2g3W6bB3(< zcI@z34}BKI*2bX-@o;0=l0s07GaebTO655dd@JZl+|=8xx9E3e%msQ|d%-ASb&Ft3 z|2CrAU+VAi8_<{nKXGwy_GjYpc!Ti1tyc-?~sD(WaYr%1Yu-!bFzdJfeSWlQz9idu-5?{w12ODn6>^@E6 z{NlfDHAC7L^2CE*7#1e-Za9rMEIRr}Upy_m<-Y8bNVCKQm`VN86};)UmqU!TuI{1f z)~4C%6BrIG5VuYO^90_m58Ac8#Ldj;fSwqcZr}9ipZDYA!4zVuNuO&p+a*iZ?+VNp zM{c=M^a=kv2NVoy8_*i=GnB?_QwdN7yfs`GoM()iGZ$p7o^p=s$rAGV?rZD%w#?`w z4`EOz`VqX_9h`ba13lyXOO8P;k)XXF%HxM)-?%B7F#v@3=Ixjpi0<)PvM^-_+$4v+ z;fr6Sl2yqb$`maIzGq6_h~FY_U^QpyHseZ1OZidIqoAaO3;|$rusN*n6`7zja{a2W zaF+KkCU|x{v3K?#G|_NbuTCBAb=MsNe(KnOD5J0OMbu(gD=xV;G5NdCA*>NRXQ_$# z$6m`Oe2VU;*xI16`b5Os!F{9^@ z8{(xyS=aTg}R$fkT>#)-Zna zeBh5rULO;eLEq!&ZD8hQJ>brc_qSef<3VU(WzHL4B>KD35u#~0Ywa)u@*po+pEU^8 zCtqZnx0ChaQc}<$sX}hdUsN)gDJQ<39V}|ENWc3*f+O8#Gf^Oqj)~dDknKZrl;>R; z6664nN`z0JhJ@oY6+Q7#^~E zz%~I1&nq}M*!cGH&o1cDfR#C{w)uJu;UOb1Cq*EjQT8$sZn=|0%dWs-ju(EsD(Fsr z4_t4YbrF;RF&P;d(1WanH_S=5y>!;R((2Obe0J~yILF71AA^Ib0k-Yz>;#9uB>^)H z$oL#JM(_!*YC1Z^kx@{n1q6se5=6(wmH}mI$YaVFYN-(ThMYJ^NlQ;Lb9DtQ4X6W5 z;=$UN5=}1r6sU@SnuXDUjI2}r?oplXJk1e}HLq5M0k#D0+XVs)nKHgbHXyj)U7scR zVLBs#h8%ouvNR%Yvc6e32mu-yI_QjdK!t$Pk(!rRb46(jzreo>uEXk!9sp}I2PZ{Z z5Rq`&fFJ#wUO)tWxg^jstrZ#3kFS~g>_FK#X@taz5ZZJN>t1%;$Zy-M+9L0utf|EjIFM) zmbAvp#8F!lv#Q%|@Fy^?=m0>n19jW7mZ317SLnoq#y5l31$$ax^VU`^?#CbC5fFw8 z?+pS03cLYyc>2VRTT1KdJ{!pvi*?xInJcVL+b) z`t3UCuIRFlhuj)9A>gM(?EzvpYJoVOJSo<4;OpQ|4($PN7N z+(2Q!5zP;r44YITuu&sq23`DbjQB+-YeV%dITW8@qnL$K9IN0{dBwY zH4NTozk=EIo2aO$Am%mmX6k!B1v?M+2L1p;pFb!oUOVJNSuaUJmo(SO$%zqP zI$W3j28Ce&$NL^n`^piVZqr{2N zxVdhYul^lc{m&VGlk?`k+g;&(s*eIY(%(t5`6Hk`#+SjHo}Ol8WtFzKX9vvyVL%h^ zu;kgez4{B!AqF4LG~pr(=xK$wPD7}G7L+OE1}Br?vig2`GY-%dU#Y3#UZ3rU{(o$} z2RN2*{60)FvXv-=XlM{+g^a9r84=2CT2aZ~86l&RCPkvGY$7u|iAa%6WUsP$&&&7s z`@jG9I9^A`caZ0K?)$#3&p1Eld0xA;-gNqK8Wuq{CwhDM+^vz8jTns)eG_C3%X4+^ zrih1|WNrzyWLX6QBsZQN?{XL~rWJFXym;cdG#mD`GDZeiF1{eXZ#aGEFMFi>vg3(J zxy=#hZ>+-*FDx|F11+ub=f`4%oe~`#&C1TcHvMN1s^Kn#(Nf3;{n2bwLA0lO?b%ezas!`^Y6L@@lWplI`S?dLF2*&qIL3*zFNCh+GSQWNKH_C_r}rO zeqQafXSZV2k>v@*2ZPP83A+1$fXARo#*|<5S?$%W$0E`RW`gN9W%4q!oQJM zXGDLi7MP$I7yx?l`|}gQF>`SxR&bXvGreZ@?LX-tKpLHdOzyFD*9(^vy{Bf@ zy!^u~Wb6L@4axT^pKYG{_tW(dSvP-2<#>M1?pfb$Iy}GakY9U{Ah&zCl=E}e^}Ze= zcLCC!Ugc|Sm6d%K{>Rn-|AK!EW35`V&gDCZLEkDkme!JEEBHv{LICzyW!|9nxpV$C z`wDp&IQE9x=UzB>Zc}*(UkKoKz1I1EBZ2>=!juPJ{4(^DBo-dHY)JX{>vdY%`kkkq z7iC$szQ6TH(*?``lpc)WeuGLzB!H-FKu73$myn!hXah9+#FPEh5azI`1@qj+mc)-k zLs^MD7@-Q|;FEzg+!W~LUZSETfC08D;1aY3wmtWxy0coXOFz!EGhiO|Uy{=sYrs2m zV6%VU#g&_zI}2oEd^B&Y8@*9E!ji_>v)fm0tGJ!6Hhy0VH9K^pMLWJhi!Rcy!h89n zEKAB~hhNv;Q&Jc@n!h8;#d8;C^iD^~Ghx*6`SS@CRn=EUg>E8N=xP+vAw#;{*r`ap zr0~&@MTs@h0Fb@$6H{{oV}rPFReZ+`gDL&g?(4V;2@CrO1yt75VEkHYaQ)+@dqIZC zIOvj*LW_Gs6oQE!cxARi!KN==FKnNdGDEq@nWDe^1t*+=FOh{r)5>%vp}lzX$EKR|$_AE16|@hqA27!!HyeUIV^ z1nzJeQ{!!p6OX7(OIvt0hr0nx8-QC^%bZf9y zz1u%J*F-%!j|m0OWyz`t*8aJ$F_WTKl-O8bA`b=PIdx%Y7gg}N2Pv#9Zote51)a)? z6N~`9eehi#r~&x|1O@RJ@lV@YiRPEM!RU|85}vBD^700|EzQ3|GCm{3?jIbyO3>gm za$(bYbD|g%dHe|~`6Jk)L z6jR@L@T~WFtq$!5%)uNt%z1pj+1GjmB(67aXe1;gkSV>48L=eEwF%x-tu17P`U@&#EebA-9A+qSI*-APt{wv)>8 zXBstL|J^%!EL$)pA{4fS)!nC8!{`ahV?_+wWo7f)Sx^SNNLMVWIUsT>me_ne*b#J- z>oLE;N;Y~Z2u4}6+in{qIjPzG%)PlYLthwOZ7QsHvDWoev(mI8Vf_@ ze2zz{A4?Au%w;SJi-bKV#7+G2-9!*R_`6y zh+VxJFVPd@zOw9qP&)qkKwT&3i^`Ukd)F!{Dk)KtJWQG!qTxkq%vszqC;ryW8G%2#@T2NXf%qDsZzksiVod(4?W|J1*tz3!LSc;@1q@g{0 zhwt_74fK?{=)y_G%nZ_A7y0meQu4iK+MQXrB`p<^Jy4U00*fL!3hTZc*Q(DFKXEutkJ+@4m0_yz)xrkQC5k%)l*V z*(7S=_Mu9=pgrTr4;;#Q$kT z%&v`K)8c)-{$qL#&ksgJS%_^42Awh3o0Hb&YyLd1d~_`K-Pw!JTxe4s3fG*Al(t<1 z@X5RW>L)Z;p)dq^P<02h!ox{AF(CYq8U|q8SCFw#wO2*?{u07EdiClR1zA)7hf9Sw z;RP)LX$;ASTgv(l7A#^y?`fspzm^Je;NYTsmRLX(*tb=roa+`bVyRK`gX5KP<)ORS zDAjrYV|_(juqcjSzB6ryJ8(aF%blR@H6vck#DwJ9fbo!-f2H$Yy|hpI*xv=IqOO3ami)91fWlO=J>t9 zJ11qXUAx8{GxPq5wUWxo1AiiI&&?!pz&C~1-4NdoJU9W^%*nZ5YWt*l_3DJ3+SfiCr@TZ?mQ@SMWCa>oQfG0KB={@M7jFHD*}U(js77a>Zm!%9*R+CA#9N? zvbP^9Jk&pvZQaJo#}|SXh^_Yv))<=(-q^;>Cm)30eE5X9_wQ>DFFA=bPYuJu!!I?T zQn?qw%ITeX{MNbgQHxF1Rr_J3M4p$;#hw)Mtw}${&D`EcHy@I6j%Zk65(uP^=J6gltz^5sjKCZqVbZ$TAMMrn;b_Tl!(Zob3T+b$ghsiZQm zvVtz!W&Z=dgm>>su=R+ZUJ#fuk>fx~W+N?uls?sZ>T;A?m2ValH$ zOn5A98;wyMdhGvC>z)JI@JPADC(umb7D-Pe6KAaM{l0+^)co|H(wEEt5AUjrCAeGf zZ;$9F8e}tj)3TJ~HsP*%mXsU@j3@fFld_i0Ce%OLe{%B19bey)i~BLtq&og+uM%oj zQpEzD--I&%zCuvtrBOQRBRAzu-u{Hi!N~qCuOC>UFhXysug4k*N$wnjJuMFkW7b!VKlI_e3FS%kdw5<g~HE80qQ03|+%al0BSS)_wi*Awmj)vWzZ#J|Mg*9$M7{L$$#s`SH;qr$>{G_P08>Qa$b&;~J*f6YUnre< z*R`1C>;2QR9I*w>9co7()TVM6{yO}Aa>WKfTb17fJ`Evukm^ypEqg-lvJSAfD)^UY z&~}jY3ZcQIAFJbYaEV2RF6JHNE`4AX{X2W0wePXG$@9$R5;z;B)3B1fB4QJIv)hzk zNuT4uR>_u@osN#>+9#fNI#o-uIjY&(ilLdIJZ@am1jOU#U0@j4bHcoq{FAC~oHn9d z_oQt)DLAmP%K|w|u}VO=)rvb+hH-=+=JucZx{As2*(U0F z>l!*zlQm@+X;JfQS4sAe(iy*>|kkVxYIJfm>fjwi+CRR-Ne_{ZQS1dF;Z%IMbn#-P zUl>yqi!@q3GU7qcIMR|?3ZcP7^_~=>Mnmxmd}h+`%%1IZYe*N3p_DO*&JMj{>a5({ zf#2dzn_+xJNZyNEf!LQ}&Yk)7U9{O=_AO=;3%pc@O9Ko%`@{@df8{O}RCv_sUXu7T z=QI9glF!HgS7Kr!BT5SZ0w`W#&6NDTy7|v7veOXFL_SVyI!a~UK<$l6GTT~PDvDVe za1z0M{x%ztGJxQMh~??0sG5$Vy8=Sx9~h{NViiKx38dTZj+JF8xA{Li(Ir=r2Vn@8 z6n-dJU>zG%zvy1HHGd0x*NEs)ZFbml3-`G~$7h+fFIs!m`6_zV-$>5QaUV`;hR%*( zOsxI`5rfC&V=>HIF+zL z7TJ`u)89_l7q~j#xpN1tAb^b%^3N-GDY1TqNeYIJTCW{XW&`Ci5gthR8F`Q+ZSU8L{ZJpFX7rYxIJ=6cGqmP#9lGtA)YOoTPcU zPV)J6NjeOWc&yv7gW}Q6f_cl1f-Sh&s3$MC(FjOQ{PZY|JmhRM*?W|fnP^g#5&HT{ z<|>~ZX;Wa^HWm3GIEz&oi#l9%OGuA33qKaCzVH3n-QnN*XL_jv zK+-{@?{>k+28wxzXlt(-`&Jr_wB>O!aEg3F&ttxH^1Ge#TaelDAoB@1|x70 zqzQb0obCWP2dFe2k>6kwpP2ZGw0j`9(*fPgnQk~v=rkgBN3AyiT}Dy*3Q;j0Wds7n zT@oaTXV*{5NBAg{7Me+GKN?Trbji;TdSmRvO_($oOA}a}hv#!5%_Jy8 ztmkk4yad?`6-oDd&%LKW2z|gJC@?9jsdL!WRDdwGI2%!R66LJnzzzuKA8lA9^hp>B z+0SWdY9_Ae?hXgy%q*&ALqGPbp=-8fQ3ia7 zG#pcA%e?q9`f;Z@6R>rwq*fb1u~m65P&f0zki?gH%iygK>=MtX z^|p2}K9UiSq_4g|!=`fOmHa*ee!*^R17dqIaFReo(oJa4eGYg1O}(#!1&{5?W}ew% z{pD>6x`goXa11N>IE3|qq|u#tCdmX@#N?~SV?kSFb=udSd-VIS-NMXRwxP#SUgu#k z6*aY%OpPIdw#Sb5}g4H-`b!(lt|DJ~Sk7VYVi! zR4*^R|4~aU*SZF?*UicQVP0E62E{1AcqY9yORC~p`1CSm{NI~{w`9O|Y+6`$6#Qeq z)BSo?9t5ALl+os%X+K+YfBSJc&2O;?gy@D2BaM)-=l3^`U)sv0g}6G?5=Kpx{c* z{mlidfufc-+PgRRCV4j-W$t08CzUNSgB+$x0pK9ea98#q)8ajnbGqr}mZrZE&e2!+ zSIvoHxcw0f=@88{+q-VAq$}B4{uXaB#u;Uc&ZSolG|6F;{;WOpCq5DG`4+j2v zd7k-C7t=G|wqA<3;EN!_4|oU@Q!1qEclyL_${Q1(PvD(qJ~^?2sEe8norRP0MnC9> z*!F*adL)W7i<0rrO@hf|^hZp+Xl}-yIH3ORN9By8l`m-*Fhe9M1XU@WSCM(D{HNvv zJeIIDMl1a=A_9=sa?z6mS8OV_*gJ^RpDG|Y=~(_$aDHwTmMn6wsO$VKN8NN2U#g?_ z%YF7=tv2u7a|T@%`t1Q>Eq2QxM(Jvs&^LxZE;Cm)*_({*eNc>r6Q%w%+iKy%b|fJX zpuos)`7M-3#q`;R6lA?8DLgRrj3V&|2JZ)}Z;ba1(?E56MPgfN@do+jS+ji7QD97f zHdAOUYY^iS%-<|7?6`@+w*%^OF+kz@`1Q#y!y3;0Jlp>3OcLyP!XK0u?>|}BmQwoS zao1NH!za7X030P-YpTCa9q4sS#x-ziGB;r$ihB!KA^l1nsEY+2Iti*88euL%8V`={ z@+IQ2+?#Rv&nhk8aV2>0;AZfVaeuA9Rz0!`NJ_7bvTD!Yj40%ml(dxxeOO$a)WD2l z?+__#j@WSz(o0sNd1M5Qju5tWjegC`?D9X;;cRI@tNt!VdAW3G4Kr$jRoD~x71SIW z%={*|s$dQ$u%^_PsT2iiDst^!!~9KseSPo;JZbvu*l1P8jp#Q=ve?5chhJQc)m-Nt zD;gvCrDdrMpM@C=L5rIbc}fD8_74bC3&-YuWh&m_Jb4z|4rP-n6l~^ATW_dUD>nSO zC2iZi7Ik3UMkoJ~)a{_uB-gREG9EbHIQQ~V&sj=Z4fYABNP18hdg>>YTV=QXem}sV zahzZpE-rG*6Cv(oV8FW8v*yxzOibId`uh7@T@{J73bd1%3YBN9pQ64(-BA_z5N50sxMzSDE-Yd z_1pWVP=uGEatp+CiT$2hO#ZX@dJJ6YimX9$#m%fe+`@O@z=4sL^z?MqSj}B)Sidzi zd;+#6hP*Cgn51LS8e6yKkE3O4?8LEHO>Jljk<(EM)GXvrJ0$50Y35Bo?NbQRNwAxD z>>IfL_2<+NH=}D7%-Dwl6dJLbRe6Kdil2$iTlS_*{SmvAyMv%&3z&;@Ii{<-0ueRE z{3GXTy?|qcR@(==4dm@XfQPHxlW3B0C)NQdg7j)b84Oj1vak~7-qSz@Vy$rbXvH$V z=NO@!6~#H(lPq=j-o1?#_qeajEDzn1nG;RxlGsep4g7)WL=kKjzowe0!_OuKMfM#c zzudony`!V!?EFv)>K^7f(P)~mQLdV%@Y|<)*C9Q8c=$x~ZDOX^orcSEN8YA6fZ1@t z$sPa9z_SeIz3{2)z z6QR5*+3a}_M{{E1;>NhSOhB(fIUA*-!M!7aoeD3Q!!Y8Tu%+?R9A=eqfqv%7QzK0s0|!y|x)>6)U+3k>_z1taRK|ZU;AFKf z(23tq3BFy8&Otm%w!=M>$CGMtbGzAOCV-px4*Nq@%Boa!SV-BLV>C;=GQ% zYy=aQ4~jq4(AN0DVPU_9O7Dx;8)mWOph-Z#p?5(0wFh|U8OhSo@FS68vSb)W^VFCd z^j221v-8>0J)LN+Nzc>?hH3?A$qQ~WQ*^jYF%c~A;X`(mT#D@@S}t2V`7ectZ?oZF zS}MK#*62ETS*4i{MP1%3^6!&hEQGi(QV|u_)>CJnsaPI(_%gBK4n@ino^C^>;;qt8 zo_WVceGJJFVUynwcF1*J95IRztPpFI+uHhh&kqJNpdBAc_6h=6Z@qlVU z=DkGMNN{DWY8jZ(F?M?f`bPNeG_q$w}VMyGH z{5D+9)z7!+|A@VBoS67s=uX4Km0UKphCDk>2e9=B*! zKcjw!Bh|CDMn{Z3^*ZZ^$a&1+NQYHLydRa zpB0IJ-f!eldjKi|spv>{;-0a$W5mmKckCq@oy(VZXliK@+$}EyY#_vS=tw@K`He`& zUOLMbnW%Lfqb4N#atv3vZ;8#9vM8KmU~306_GH2L&_RYD-?xRfIM^=;v<`33HMjhc zOxHU&9at`Uj=l#?#iaE2HTv7&&0N7OH(4n-i!8d8?(hQ^3YB(6dKn&hmrbNw1dKv_ zZvFM|G^FvoRnb3{UU2e5;GU(Z_=87PPf2i{{8$!lHS#l^NNNd(uIsPbqh{QU=kVp9XPfsT_5DJHGs5OFB* z4SaCucVmFMc+^&}OUAhdC<@VZ4GY_sel;{V>^Rwf zY?TFr#c%|oH=zLB;mNpk)FFts_obu)3{&iP4SD@a`ANk*)SOZQ<*XlJ-xI*n-g)Dn zU^}Yd(~Z@c0E~#+2Nd{3xrdUj#t+pko5mIpPMGltfBZOyZjOeU0VtPF;yKEDYTx3+ zOJ9qxAkJY(1V*#P2wkh?xBd7y2t>M4%T^O#7 zgf#v~41!O+9>^UHyxH*s(MtMT$6k*KFB!p=XI|e)9+Ys%vp7|vf7cIF`hmhg_=vd5 zYhTwtu05@I@+8xj7YBVBx(^GnQ%+K|OPYkXI9ncKG zfGWVwFlI4D4~a=ib>0w8TcGB))LM*f*)hp*@lUrdz4(!Dbgw6pe8kTdho6+stl4}m zXgczkrw>YlJVj549~bni(}*AORq^k`bk}51jCCH^xid~l>sI`=r-wR8qkRJUafpW| zMN(j!!XL0>xFdJtaXUVfmE=NJ_Vs%KO`O5@m$^0Ifi@GaEMT#$9$LW`M)f~X;*Qgi zK_Jjn=e>}*(MyMx0j|47d!#8i-)J8xrkR_36ilyt%)T<#G$HG??4L#l5W|s3xn93_ zRfwH20P>X5+zPR4r}uPT*c*!@l>EEQw!^VUMAl0#;pZA=0AtdQliPrGv2k(fkMX$% zu3xKtR>yS4h&$}>)D51tLhiFzl&s4uBkTG3MMdd~x4xpIWAu7?h93bK3FXXMQSRu> z4QoNm;r#@-EEtMKR7OYtTj`dNe^FxpBcr){KPqBZWy5v@xeJIbqVOY@w-GUGLK|M( zaLl{af1dY8sMm(9J$h$Kk9tzI&WNE78eKYt8fQ;8yeATkTC24#`mdu zIyKL~(YC;au|huBB;3{IY}mj~zscgAZHdj=E~^$r$3K12@11FHFsZZgxQNv#ZKd-# ztk^YoK$NFFuo zkA{|zdJssFT8VuV^&jA*{9`nGrMw59YMMA_^GN!zn5^4-7X`=nm{#*+M1g`e zGmSUidtx0lDmqCWe2%AZIa${`f)^KE&b$6?GIZT4 z(rJ3TwCcIp)sEB++A-x_pWn(Yezc>y;}Lk=fu5dvl2JW%d`VZVzvJ8cT~k-sd;3lF z$}1ipm|C~E*cKtQONK?M-tziB$5({Rb+cBY5;n_-=Bsp3m0kdCW!NOZg9z6LKre!V zA9m*VbyxRZ3ogt^J51b_-_bNIy8ju?8KF*Oe-4pdl3lUI`(1`dl8*x6S(Z26V0=a| zK>k=EZFp<})!srr^(Wh%eh#t+@8_*df7%z2gNCDUagrgqVBVaOTauA}-$z!?wUn11 zu=Stn6Wk`w7-4j9(?=&GPf)0(b1dNf8A84rH+3_0vF5af&3+b^S_kSXvvNUkIi8sl zk-i3FPyK84EJAagHj6xr*y{0ZkA4!F^n#KFd!ZuilQVjHyH?Ti#z9xKvHwvT2=Le2 zcIv1e8xiuIERNq(3X)}$6Wi{MUT8eX$PHkh#N$`MU|zhqhp0tjE(h*W(nd>cZ+JTp z-X6F%dIySBHAO4_83DR_p|&-ZX-!r7CgXMgjwsQ+xMUXqELz`AhFp7fA+8;ueo<0( z`(V=^EySdE?~e6GO4EbhH2jmh&9ujME}Y8tP3N}Y(GI5cD*s{Sg)j^00MF*_U1-dvxc0?0nTI`n zYR+VM#7-dI8teb|SOR-*>IJ4xx38Uq6RoJ3()^a6*1>7~xU}t^U$d_p5lRDL!F~J} zT*yxcK4GP75UcNTQpEnS=-&yQx9s0BPYo`3AYCucM%aLkmXWT zjACw!4nX0)$Bqy&7oaa<75?p{R8>^V5X!l{r&KCQJ{44`ln9Ck?u?;q_!?ANrm2!HH8}j?LyV1!5&{DyRbc0Td#wE)H0+%K^H;6lqjlaw zl;Y5Rcjv0bK78fDZwskug`7a_QuRZcCaZXhV?H~NQnU1FiTLYNg{V`Ib_eVeRy@7w z7KS5Atj3Z?pRJ{D+wcYR$@<#u(dD0t8=GjY?)x#nuwzp~La4g{t^4&0A+x*_1Gd-r z7BhD-;M6y2=8-nf?vkk=>8{*vOI@kulIg*D1cpZo9N-*({j8v z(zga!@kIP!JnLXzq$D20ufRh4Z@7qWNXM+pImNU-|KL~G*&PM50}4M1Eoh@n!-e*} zy;HR*Ag@WQKG(uNabA*YyK28}O%AQi)eiZ^4~x#z9~S9W zbRn;(9XaB!d*ze33J`c(?Ccnzoq*U^Z@0jtGVm$VSaL|h(>y!*onoWD!*swGOG!a! z0}0s^&2{tl>tIWqp`$;a@L}_H*AMI~3vP$xmrSVSFDx&Vjres;RM|T%&eWWBnEj*6 zx_P5FqLPEtKXQ*BRoG2@rv$7sYs|1urLbZQ)suCiSP{w z3i{NfYs93%Zj0W9<9PontvYF~($()bMO!!SN=r)%4G;GLH!zNR$KAznbt-1?7mDz< zwvaAHW5o$Ss_xyLUUlc`;BM`+hKhD)FC`_tOIu(5lndPiHY54I?%9&7iJny#RZgCi zc$ctSL+irGMBq)`jhvi4#0mmUq(6#7PQW@UZf^4Jt}_?p{{1kU#k365c$1mOWDuMW3JCwb`$Doq>5rvMwZXioNlBo|>b0FC)v=#ICNcW3v+v0n;TTPw8K z)2%HcB)!S8g6Z5`(~9Vg0g8C;ob#c@#zmyP$54k0tS|0Jm?ppxi+oIcv{mDMoShMa z&g*taYx)%yB0VL9To8M~wu`w?n75Igy%S_T!bt_TkguOou>~wXaESNi&vpl%%e}Fk zh!&Ln^ACgAs3meYz9}r&8G$mNVOSXrG<3HZI17ylun^${ln(Vssul*o>1@IAYrSXT z4}UBgmGkWIBC|-&rCIeV#Y`hvQIj@=hR66V&CpxfCPwMTlkBQsIFNBU?wC5JuH&;k zgtf>^K0?9LKLbQQ4#J+pHsz=D`fH!*x^P<#tWlEUM8QbDHTp|$dD>bqK!A9Bud)_V z-`ap8>U2>?^*AjLE(^UeKunI6#i9|keu#lYAl0J^zO z)MXWbNJNc`Ci579di(lJpvVBx{~0)psPF4mKfFMv=;hCw z95=u;%NA2f$D8TIA1fd*kPY-Abk&abfrD*%QiPdN?@@{onCHsUn7acOu^7tr?2gRK zMa8MG`y$+RxF1r%zl!KWpz||>N#A2hD=t(46KE?s)9d7wF=`m9IForXu6IDB-jAusM0O_d% zM(!o8tG^2~6?1LP2FN5CX$bkIx3>23++-hs^PO)_M+Jk7UYYlio{HIS>T~#uG6+m3 zK+giV>{z=huPLx3Jc1rV_sy~GdG9b`N(jb~wIBneo_c~wIU*I|l@eM4!VVF=1aj{4 zEjGsCs8_#gKZfqWvL$03!Nx`6fbo4#GME6J784xaqK&%eyHbsnk263bgH;OLvoUew z{r{vV`T^=GqgZByU@|-FQ1~%)TJnG9X&45D@E-aK0XnSRsIXcN+t|Ff`^yMy1{;JZ zI)I#lQCx&s?-?B(Uku}$^qfmAC7VB7)e&N+Ter@3ev>*q0|QPaAgsS%j{Pr~A_KaP zV>^%$xDb1QeoZLO(8Jf9V2BzemOA+nP*|D^{+c~!b#;4?qcNt#({>0rX@+q z5{wotEvodNmTfY2gLgn>{&6pTsSPunVT_)kLprHHuf?uP=Dj4iVZURtyy(QYj+vRT zDp#p^2)_TNBAi9jQ_sb0`q4 ztnUF(6XyanND8ALXAH?XPAMP%aS`%C?}O8(`2TY_=Z3IIiE|C1Y@9iNJ_yfw3t>G2 zzP5$xOKxE}J?HOt16#~PUt#%CcpzqU6KMg7MSzYhE(BxhZGwT90|^JjP}1K+edr{x991QcH#Q?NKv6!^%!cd9OeSPPWCU ze^%`!kiUs^6i+b}%AdhObHooEc{4rK%z|moLBh5Gz|n&wfFoNzJr*;2f2EH1e`~Q(BDn|xaZU`VK9RvAcBUB^ zlXNUB71%hyOuXS4B>+VYhRVddkSNSRw}Bi8pOt6=f@McZ5R7q;0|lB;SqOLtj9s5a zXgLw0N2ibNT#H+A4bULl^S;pRx;FY$v;GG||C{yf^mvU62Fr@&`PuI+TBH0==B#aI z8_2n3At~yGsN$yQGun9aOTh)Z^&gr z1_GRc$_S&PL#yYdQBodeuoPgY2a3Sz=}FkQl?Rv%RAlNaf>fM7@eJlp+2M1 ztDMK{>pv643PJ#A;xba#&Jt*>aWK6NfshEZN!mbJcRdR+C?PSC@SVUaU38LQCZ9?` z;I8@s9{g2efJJDgF>3d2I1jqO5@xLug$e{wy~d?2;dp;V`PZB2C{)-ZWxWL#?0YPW zcjbA7QF?kF8MZsObG z#}t7em6Qm6O&ALw_}PDqePeHz=R~%IMNt=4jsxf&6cqThP&?jdNtLl5<#Ah_P)8Y& z_2uanvOsI4I)V_Uhv+$LW4BW=qo$Sy(1iHygmko<{t}{T*oqkn11Mt)%3~<~K*+?f zfC>p6UKamN2=2dG%i|G9%}qr`#kQU+z@zB*$60EMva)Rq!#@^xa^y>K8@L-W_y|w{ zmx*rf)9w7g79;8~dqilB_xMl)V?51Yi!3hEv5{m<)N+ecQ7e#@K4eL~i7dWafY!M$ zdjqkkⅈ6DHDP>{8{t#{hD)Qk%kXT45-U3^JH60a>d!ki z?5?iX(AWNzXrO1i?pWY)4NdK!$%gg*mk+dLjCd?JGnPR+Ue?HIRGM>_r6f+s59NSxUxaWUwPwe)k4W!URggmeHNK zMbfA~5KsItdLZVv7#Tv^@c#YRlL|YWHrDGVaq|DXqM$AKQ&2NYpGuWJ2(>)YR0!1yYz|mjF|^ zfS)6;06HO<3KfeS1xEQv2H9jFL}p!ZD!D4vW(vhLfyp3N14)Reb-Ib~DWR}~@xz1B zJgjXoVRcxXs1Um~a3*hT*R{8d8Vat0_-nj1Fhla?(l-%swGQ3k9iB*mZ+<) zNxWO>SfO9wYY(Z_D-_tj8@RqWNA%TlpMV<~79DIRY;2+}7wOkh(@I$sA=_ZWMnk-J zVeR9I*!n@+xwOfnl!Q852<^&3#|ohrbt$aOv5=5MSY@QkBhDXal-iPu2Jm2sER<*s zb0*5U#jIMyiAxgpZdrBs9%Q%o&NJ3z@orLZ%m3yTB)rjvu!`DUrcNT^b)o3WT=0O# zQ(*EEl#t};52+sP*~uL@9-uhf->E1cqql~csOrG0fc($En>m?=N2L7)XQ=(6&?Djb zO+HsnDccRb0NNMg?QDN|tB_3vl@=la#bY^XJ0;e&y}x-%XN;CvxW6Z=OvcN;a_gqf zuN-u2s5T&y-^m;v^RKiUljs8NXnp9ZF$ zpiG)#bWBpO4VKV%{d&}Gw1}`J#6cFR&Y<<=jtvABAhD8E$rmraZtJnk)Gk^*1`VwU zF()G-WYwaHq_Iiv|F!m*WCgqVS?<_iA#)-0Kty^4T8?3R(}}hIthS3kOomsSxz*a` zDeZqndDU;wv`HG|P^s8?n2nY0+?%s>UL}($Sgw>x@7sF|k#>oJ zkkTP%Y|3SPBltxi6p@$&142YcEwdhn3~oh3HIjDWb&CDyue)TqkhbP`UDS2X*X!3( zBJM$joMR_z+kGc<_H~|@;3%kq5nX6%K4BjA%yripvdsmeA2mRujtb`$=IU@lt1W0L zQ5B$|w7vV$rc9!_tADJp-1*vp_pT3k&xUh;El|FK)k?XsFLu_fZTEn@%Z7wUdrLU^ z3cc+Ooz{1_PbjcQ3gcLB)-O>9efT)`eR4N+ThI5!jafl>VU17!7e%(0iq`)$T$T!F zCcYd$-5n!kR5E6v6CJ(SePO0_;kMFD@sPtp?9Np4s(zt8B6m8=|7<^b7McUlXJOG! zhYDK_#Za&EVBg+OqI~%LOo|1ktB`t+=zYOv$rHnPW5+p#0Er9%|3P)gI-5YGswuoaW^jWf8vN; z0->W3wK=JC37L}Y6cTLe(Izg!RmlNnbU3q3^inlDP;`hyx-m8;N?k(JAaW_9U;wAn z!qp808piU>2M-=3?$$(PWD#}#X0r{f;^7T|BWWr`t$#_!mF*>+>f&9YfcW^F$6JT* zLBUUhB=>iy8Ac;IDk^tyJaMY8_sZ2r2E!odNn{|hpqum*iy^qga{k2l}=kdP= zgCjY;`z`efIgCEL(_-JMnWog(Gd}|6aGl`wLeA{A6eYx)+txolKDK5X&lo+! zp~7Q7^4*00jkO2g?>^icCE8h(Y&Gj#!}p^=i~)s`pBe4UtDJb&gB8ps6;XOZ^6t!1 z!5_J$+T2)@bLBU;nRbhO_Z{~te1%GEadD~+qKlkAZ=S94Dnp&@#@e^k$xSgrF+1F1 zVQvl%kT=pl(JA3l08^+!Q8)5N5fwVIsU{~OptELV`~-pnf!u~}JjV7%6Iq63XW2%$ zA=4wRoM`e4?Gwq z@@~+(qG3lVzhms5(4UfjD*sFt3x2C-jQ!YrH2arlMWp-2?6$b#t|fP~<}3xa&23D! zhn!xS^&h!4!h5eX)EMm=|9%+`R8IYM951z+kHAyLHRnpT_6euOgK{h7azT-z>jFF| zkOjK4oV7%i4s2wpJ5rXO;^uGR*KsGU!4&mfiC{~y<$f|;vTC)zQ^N3u0?k(9`P}s+ z%)w?3;YEV;&j+u@(vO*r)(7;^X5L(PuzM%z|jF#A?%3LX93r zkzJMWJaH@lo&$lSGXJUD*=4p(Q7Phvbb4OWeTHN}Dk@S7TIY1pfN%QWkQX3D=H-uv zYwGGM0LiDCH@)M??QtHxCBHnc6_Oq|MX8~pZMlZVFz)dA5JMpm8J6^;DktqWifBLV zIjnk6r9NnKQuATNqiENPQv*Dk?+(7(opBB0jg&1xO4$0BmQvIrtlLmL1|U_bsM%E% z4JDTr4c)gnI_7q@wLlOcpoWfSkk{e{Q7#KoB?hEvK77RfD&sY4^pC|cAClLJppaxg zlTzQ@k6$2m6~bQe#p!+ov>?`f@OQlnbc`@$3zz5jpd3W+ivqPvI{u^$s z-Y~RP(B=|JJNR`#%Jm$x)sNO_{{GEVMmSuh8?p*A$zqAg(@UU@7t8UwD|h!CTHTx-<^pUQrD}qj20(bxtrDn3lDtF-GdjH zn6mHauIq3s?Q+c{$2Cd24YQhVaPB_bM?K)aOi6ot#V3aLAAkRV((tl{Alc?z$2C6v zCWa$G7ujG8SwX;W3}LfVrXd_Bw~ycfi1(|8gJduOTZ9~toiKg`kanloo0c$>NjN0_ zqkD}Ccz}?S0t+5L7`gEg2LJ*@5m1P5mI=*u%zY&g2*@hC_~Q)k?cpC=-PZ{dRFahJ zWQ3FVl!HJY@Cs5*XTSFlC;@DTdSIMW3gX&0B+lZ04!1T%hk5Z+2-xhG;?B#qQn}Pe zZTUr?3i^^a7i7;%nTq>%xX7N^A?Rk1$!MYL5;*QJDf(SoP*iAbbB8pSR7f7e1pP%m zdH+j+^5u8U)LFTBgRSdRb8Ls#PIF7I#Y>jVM7akw+&4uPuKaWA#xJKC?}1;L+~fuI z6!xhrBKx#t_y~xyP8<>`tF`|XSQ&pV?k}yJp7iX@9`6COA9|*Nkq1gh% z28V~ZoZz=eK1N{Z?clxoCOQneAIL3GQ@rvy?+A|Vq>`x{{J-rs_u=Re#*c<$(_D95 z)YAdz=O)PGG!?vJConAYM;_jvEcX|qM||^v{>rJSwQDI~y)qeSjO8v7XnC0y%{ZhZ zQuh0=h3MjMS(|Hxl)z@UqCMa+U>*Q&t5ExIqvU#4}~F-)&m(tWpLgb zH~%%hBxoo=m|wo z)|3y|zqvp=J$!S!lVLkw>beO2ZQ|_ZPtAd)1TOyKdC0yIb0#&!uk|9ij9p!9xPio= zk@8k>P|5a>3o*J;d)N|&?{cW1{)2R8Co}sdI8jnl*QZr`2^sE{+4Xm`uoC?_5w(Pb-TVHXJ^Sj- z8?GHaWkGiz2?(FA%5a}S?NN?RZ-Xgz5SIi!@-0Z{t6xcKfoplTktqXllfcTX$$J|9 z)6)0w*3N3BQc<(-?qX9fPK#Y@nId(@CJXc4(z`quNF0Xi4bdB7T>a3Khpa{c^Co@q z&ab&QuOyyJT|*8H!NjAlt@ebX%Cp-r$MdTkkiPYY!{_VxILEV5#R64!RSb$Yv&fB8 zpbNKEQCI(u<^^}c-NIZp#NrGU*EkA*)dNLGFJF;-hKs(3`5`_AXsRexhq3(N0yJd4 z8Ha~aV`TqKn=$;q-i9w@;H(COhQvJfm&;Tw;(oaO2FH(pVc{-dI{GUYZLm=~=i$|@ zBrvIi_oVPZgGD-NVw{UIf=lD)O7}zB`qNVWoMX`j1p@mt5mEs|fXGgaexTX>KhHN4 zR9^}T6y_KZ6oY?_BVz|wEXL}3w!a-*m{@5_uM=gAFyS_SY``S=jDb_{VP=B(wt(Wz z%^?D>6N;|d+#PpRTzjJ7tlNE~6XJZ{T=y$ygN}2t&$z!YaoeQk72|uYH zZho}+iOCx`_Y)s$zQrXt<+aV+?pSb==UsZ78}}~3VkCJg`DQyU%C8{#<8OF8&WkMD zGZ%zD7ztca{My)@I<9-(Wf^quAJL4RWulqzR~L+}V?yvnJxg z?alJzo-Cci>2>88j+GkU@}At+^V#<3B!0c}r;3=X5cdEwe!_ZeiA*4iP@F4RZt@5bbCENl z1e`+Z02)IrrK&kq?l`%m35@h1j|;Lw^?+}u-((NPdw+L$Rs#l+ZH-hS{n zDnP36^0idVF?DHs@wc5mCZCQLJD^EPOEKG&;^uSOeg4}I%*e-gx=-$e-Or(Ez27fR zK8q}qDwsIA*h#C&|$Fgc7oe+9OI&`li3J{<=QbeAw8pvS+BIkVK)G>BdGH83Yk*Nr1 znKSb*)Vvs$Osp;VmzmP#%WX}3zdEIvA(DJscZ}Ri!5~V%`U;QIUt%{XX`6rZe z0zUgvKdJjb)Er~M-=ynZ&4tSnQRCl&TJXULSquQm>g$hWd>aR7992xS8V#bE)nq!! zWa+u%?fSA&95-@#@LfNKefg|oF8#8pq#uTMoUW7=1BK=*JD+C3&8O#!?ARyHf^ll{ zBkA5yB(+R8Cf_tB~+uOikAS2CPf98`9wcYDq!Snl%9>?$GrZwsf4t z*9FpkV<1`0>lvkkEvsgvvmN5rEd2Y+ra@_uNPqvuIf$*);$V9`V2ItkA73aQx;oN5! zgsn!6#U$3g{4d+=XWAu#{tUNba zXl!hMP4-eTQ;@U#PMkVblD2mko+TuolEu;)uqU|)i2F@)q}IHFKQIRqXn-6btFf&S z&Nc=VOgwMLY2C+Z|98PO`lT{S$KzS!L>@03a@T|MyrN%W#c>g*m$RTxz&4rsMS=u5 z_lkHQMGKX{W^*FgKA+0{-xYUXY%h%{a?}aF=EWZ!=NEHnqfWadK_YWv(5 zdsFnV4tTz)I2?RQSUiNdF_0m?Etm_P<=odlVMy&y)*JMm^;b?xBrZ}EjKyfa2X{-*dCK}aunj*%c}aP7Bxr0@L4eAMNjJc zW>WKQ<20en)LzHBisYv=)_Y|5>-Bk&(k!cwkJWxXAJBZI4+O>V?>pVCQFS1F5RGLH zg41|zohviuREJ3Aj~hg*s&S|GnN$xrF69VJggZ~x|IX(&Ylv1ff9l29u*?t96DFGk z7)u}RJKdTEuT3xhtGgkQ77%!9@jCW8N2#d(Xap`l*TdVqK3bK;l?mK2Zfq2|$p2b} zVz5lelfiuVZ}!*6fzhi8Z!~}El*@>PQK@dh~dIkB4hYV+AmaNj> zXuhS`b-K&7!_Q=s=*jfA&&eXVDzZDG#!2IFq3)Xru{bJbMEW91S{{#?{Qg5Z1JQxw zS=k?%uF~RxH+)KoM6ao1I(uldtT$c1DbvYM7BxP9xI@2e z$6e%7UfJ8jgkB!B&wZ|D7LGh2Deg-{4zpc;f3xd1g<2o@>$HX38IlHI{bmNCkn~Vu zIugdCeMZTJR@;u)(|;4kU6!k(>K?Gp)IYgYvRV#3O=#e4Nwq6=c*aRddsZZ(G+H%K z^x(m}m`rcP>gFj}8oAy>A}*E|*-QulVan#gIS8SGjLM76?n|5n)3N)yU|3a#9pC{h zv!`v$O-iI-)-->iZLB=OXvv%ihus&BgpRV`r-;+uU7Te_!?k+D4I2y&zPBr#eP)j= zd@jdk>yDi}OI9m+b)Qfh@$2h zcQrCt!}OGgw!+p}Pv^#oZM1$o8=gpw#0HBAmxvCV%VQfv=bqAEOtUNfC7MFTED%v+ zQtU-daLDf|Mj`h3JzGecpFlnoU{ofJmO++f@u8Vt-h9jBeuNq{t@<^TkP>omme<6L z`^wB}I>QI2pt#-PYXAJ^mxDQ-UJT>NVVmb%+%D#w8hG6hkQ4nxayqu(5z319R@`r{ zm1J452bNuI3=n11*vsK5Evk{EGv_0xQZOcnTdMgM!?>5)yI<)i?Wib9+uJ)e_lDo| z+EW+-9QDBY(HY_!*!K(E?RTuo4U@znZ*2AFx6hOnjG396)1wxt@_wALaEXg0H~v{x z$J7aSRdSl5MoKfu*hmtlFnE~70@M9p=vm2R_O<+pRM@m6b zgMrG$t}s`1dc4MHc|CUgH%|2}#9d{!i^bhP#qDpx>~x&MJ@6^khsaIRvAU-sc8wL6 zaT)~y^0}_##k6#lbMJsk*1yz zfe|E2{#@<#a$+GCg0M^KqxbNNGs9RR2gXDiA`!M^U;_)1;uYj;V;na+FF;oWx{p6D zm<0ycq*)|kmTaw#rBVA;+ev7zDu``4{tgv|;nJXPPtbucMG$iN-Mzcc)%K?5t20;U zTSeh|8OpffwNz_d7mB{*p*5S8#s4stUOjkQQ>5)eBM-|52*kwVNz@e|_2NIWhH1U_ zK%=?RXxzJ|-B7f6@f$OzixuNaIQE|Nz!3#q?@O6Wcf~#T0;3B0WWGQS;S%|UVQ#Q50-*?NsAy*bLSHMz}*7-qX@UN56jC_{W&eCe< z$mV~Z%(e2UG(1REP8vi)UnO-HAq$h(LKLjTm%wd0hHo{!_U$pci}Pi><70wLFqiNn zOOR;Z;G2f~e)ph8t{ZA2JO%qs(lLF)IOrNC_$q>jImMxD4=JU4a29XV&U#iR89(n{YozGy`@Aa$cZ7YHIDF&a zYhko-h)bV6NIYj^zA62gxo!%V`mV{m##OWM;*%FQtm=nN82l)4%HA$!xpBo4Cr@@g z)H+qHp}xh3C~2cy{~=JV7VZE!5zM->3zvimysQ@Exw?;j^~?LYFZa8b6Jtp-mn6kA zj4bZK{fGkGc~Q~^pvKZ#snp9v14hER;sR2P0?MfZYUBw$u*<0SnbM4-(;Q{Fy=OsN8#bkEJK&#+%`1+v^k%K zEcw+69a-mG##J^rqTH*hIywwQm-MZqRU;)qfBpZEcF=?PG#v=yAq}5Je}5IC9o`BdeO-$IqXCQkS1(<8c~zMfiKR zn9|`5cy3f~t+GI$JFgxRZfGFw0wz-h=o4VTSeI|dnBXm zjBZ;-_6j93QYp8Yoh^jQ%FLF1CzR|hO0su$Rzi{;!gF3ezwh^Xp5yudaeR)?(NWy) z`+8s3d7bC$e7#<0+$@tI86lV{Kt8nX4TG1nB!fn-v4(yiB!;+9J8Vnie_iP3If_NW z_)U-6Wto>gfj_Bi)U4+qTY{<+rcEd+3J3BrKDbF6(J*Ms6R1O%ka2opchQ}>?(Zi7 zlqUvxUM0-K8~au|>_7TNnf2doEz-5PKK8RgPNDr{&;7Do{$FK?T#53;4t* z&Z$H;hOWzIc_*~FzuK}8XzOWfX-y}jIK_BXir`wcrB^vV@lNAYp*PoZ&!gMW@*SWb zhEeevE4pubSFqI_m5g}b!PoB7+FaWz6xPaVVGpRQ^~nrrDaZXQblG6L2BCQb;B);L z1s-&*BRf2p)qYg0f=s*V@L=+g;Rw$H$;(}&bAUsRLKlEG=O8a<-va#~FN|uzL)2|3 z=;|BcMDrkf=%dsQMLs+RrYi zf7|W3g0}hfyqo#&c=P9K<19+UfF2}Q{^)8=z4xRtZIErzxZQ0r);Z7DegDwKbjxN+ z&`F}yPG!mNW!s#QQM{Vaiq^I!+z5!Vc(9OmZ~}#$jaCs-h+!)IapZ3n*o6(>-NA=E z5c3Oq?_%Hm2;V^-f$Cues9AvIfUF!*pb&$U_r2YC)%_80pALt#423qZUlKzTQov;x z9}i_Gv_UUtjc(O?xbzB{ql^m{0V*z`q5$$tR-Mm0apr*@IF(=rUHKFGc?FE}6fLh* zPf^o33=yA@(FmpwTPYo8g^OpiYh4P4{z2l{!tBYrAfZ|teN~<->H_bJd3Q}dXX!5W zDiYGkZE7Mpg93S`mPVviu=m{?MX<*vd^*}kwElu*)paD;-*opMafT`|u)Dym=Wnlnh)=S5F@S+N~ z;-g?Q%7>0XRJ2bJlMeFeM~uV|Q(Vri4DSIrz}seZd;TS{=rXfSY};JicYhyW_U*Uh z?BA|vo(2b0k&mIT+ps?X6(kCuxu5zyV3AZSvp~$}!8OBPenqwv=hq}NGWPo=R#dw*m}#yh^JMZ0lm3L_gg!^pj+n|bKsBC+s-;_4!kx4nnA5dI+R877#w&2E&cemBPGpB9+=>=(jyiCFa+Gr(; zCo6kU4Uq_cj+~3`SnX5q{Q_G`iq1TXCCab6Y`JjvjSJiUSr?&S_Xc-yALE znArfT;$nFWTwP>&8jtVSiNv>UdTH~zZ)o&2Pk63`QTqOE+T|4qfZfA}6Z?AYLL>3M zd`8a@FmU}z`5TC3|H^smDw`+chD18BF0ZNTR+=0b_45bfSbCut71L6mcp|9vIPs<( z?QyS3m_n!$#s8FP_eV#}TfYF3XVPq8hZlz3}WeGYPf#mCx7 zo9{G#-)oz5)TS&0&@jbq0%eB;c>GmA*&x`+LqI2ifc|pUT4&EYKa-FOcbv&4@~R(E z#Rh9UHWAlzw`*OuzJG2XfP4UeNoM|0W>VRFiE5hzfE>m)em*7-)$cHPMNxzS4t(?4d z20EzRv1)BzRlt-#^9Uh@NjZJ(cN@%AW~c8y4Xv7LnIScRC&1dE5%oPVx!&rX#sN#v zdDtIawh-y*-LVDYF3OA!L@U<2G3P9Ga!W5<*rN;TaZK*N`gE!4%bepQ#@spKUMVPo zHSvQXDvN05CVcX1id10_@9Evandk|pD7u*pr3~B&xES9Cu|bZiVS7L49A${FRMk!~ z!qpOOv9X2-7(q-_bamuvXHbkBqYsp7R#X4Zx1IP?=puW|dp#-zqicEVFzHM`A>aW$ z8~|0BG~#JaX^{gDReJYkYUS_uSoUZSX0;C#(_++i3&`2*o6=RcpE*RviJT(*m4VDd zXVj~nq5AvCC!df)3=iGbkJAM28*W#cOG$aYztvtd1hRf?R*Bu27uOJ0thpV)vfCZ} z4u*<|bSMi4^Dzdoh6-mxFH>`*^5-6?WQP{s346&_P9v*9)G9 zqVZqy9(eE&V=s^XHvKU4&rOc%JJiTR;S!5^SvG)&kzHvA8ClTPDm^NfX?JFgzF)g9 zVvbNKV-@Ci3A8z5E+#0q)RYT_%r_6qoahd_=ZWGUTG5pHss)_LVgXDAgRh)?WBv;y z{;!WSW$k@}s<2nz=bL%3=)NY(2=Svi6_yi;d&UNB@{Dp0%tyvw`L$(C76 zZbU@{Q66-+Hn@P=-GL&z-J}Q1gcO9Irv!1g>sFl3APEs$;=0XfknJcucI=n7g{V7C z{rzWha&i)LYC${E0vYn4{!s#o4@Pu&)q3roDOLJkBBjYR(CO{*(&0S^FxttdAFoy! z{H;qPVgL;In+VcNdcDKIyeL}1A(z{UYg7yW2NCP>V!fb;a1e-yD}xkJ@;fj>)P!?g zM#O;k9hwgBv1~cc?hM|e5p3>M$c6Ck7uY!HiCq<95w4-VZ`Z_Hrr$kjN!Y)yHc7#_$ zWpZcqit*4V!-pxc^p+^KoPeHz;5R{kShq+)GU}l7JW>7c)vr+zT|D=~G5GPL!5J=>i( z=72|WyAi!x>gA}ec|hA%m3=pq9|My8;uz_Oft~?WT+p(h28jhTS|q?F7Cv-`So*T# ztfPKHy+UGJr4MScR2bq5UzY&q?-ZLML-pSOAmGRz1NAh}k3;{K6x4b$GP^wy(o2J3-~LM!3siI0BHPZe=v!|o#(38P^&3qL6y9X~iMk4)lY zRji`onx5?ZceX5L%p1%Tzoqy)8=7!?O+GyMuf>>%fhvNNjw1T3s+F?RW%K{}YXS7v zliTLB$KjDEuIM)YBmeKWPQq^`N9_d#*W7eYtO@t4`lC}A@&hI>=Ezm)0lv5x9g+z1 zQS2~BKVotx|9fQselRxDXUEG?GX*}v|M?e@qFJ`u58^HyWC{R-!Lzk0eWRTt{<+K+ zHL}6DDO(2j3C<^DOJj@7)@9bt^(nYzO@k}hNC@D01e!x0m$-HV(nQC~l@SW2UaRY-ptZNJm6bP(h{)& z0Fq1r4x9I>FuLQ6yHFrNF1S?nw!59MmWUCRgx|`>&pd*B#YsqTP7F-+uT~n{;l7E7 z6T*8aq&!jAzqFZ?62l+-;of@A@)vWb_;#M+sD>GtFr z6&!yjq+kM8J5I$ciMtBv#0*rjCWGrYp%!J$Hk&Xe{}P(-SQ&cp`R6UsA$zUQ`+xty z{XMhY>s$UTqzb^mG~=PUcma#&F+mDWm(8Ts`) zGRyfq?Zf4ufK6`}SfaG9rsB78M1fI)r|_h@lAKickACo{8w}N>FRnae7*r*|Y*@4xMBG5zBSI9_1U42Om1?E4nUT#n3GAot!efyWrT#}R z_JbK_Hy+nm)`jfeO5M-l!RQsTCcl+EK2_Ex!N8S0TqBrJ`&C#djyR)u&X#x^{qhV! zGCe1zP)WLzed|(0a(<{gAj8;#KYjkPOYa-Gl?O+iAVsIs0759RX&=s84LyP7b*voI zjT1F3$0pha_pqav%cz7=ei`_Xg@t6~;Yb^eA0$2yvQn?MwZ(+}-|sIdtH-h(X09D= zaK#ju(I&z!WOOuq7A}cqKf)d=mH3l%@=-@WVp-RE?fRL8r;7dU|F5l$>WQNk7$8-E z+=>&;lD)TO_+MbW^ZMl~uvy`4qu*PqCizf44kjMDWbei_c$ixv9?) zGs5HY^+^a#R1vz@L{9Uz$3ci7sH^Pl^gy3Bs8y#G8W1G@PsTmrw);yrd=8u^VT_YVEF?`secXm|kh9;+?<*ww zs+Z$YD^nn|3;oz0OLp1jPj@dHV!VxB>-fR5_lLm5XBcEo3xT}V0}d$Q^fg8L+@)_b z;ROb)Nc>rMq?!bwA9kWy`?*W%JO_q)Emi(?w@A_3B#(P9w|0to6Zn^sb962?2!>)= zJ}fMjK5ar#vB-PbnzpCSTI()p@B-B5J)h4}2)8&evNA8ajUgs5YSTX`c^dF9yt*(c zhN;5!RD|@1%7H0H{FxJD?n3I*&&0RIx;50Tr~cryy-&)5Lub;DJ`fK%Q~->+c{fqB zhs4nYBSe;odkgN-g-c2~SS5h7utu1#bB;u;Z*J!AW~B}+ZRQ#CFMaZMFu@Q|RTk!J zKhH2J9f+pe_E zL(s^B77zS=$=h`e?sG3*&WL$9yqO7?8WI?W#4U{ zG-lo(B23CW^tj{XysYV_JT&3#w(+dKQp%%kbdrYYovVp+cU7ynn{6jCyq&R#Q{>Uu zM85ekKxW4WUMRJFDSZUlLg7+Q8FOUr%^WQFd|uI9QcN8q&v@Fz+Uj8;{R~bGZq(ZF z96G8SOPTokJ#~(R8<+{@9@%*Cpl+CD?w* zU&o)LfYbz`deGu|2R@#FUBiE@GO31#+B{2}hggg8f3gK{lMZuZ7Y^$|HoW@Ppzc%? zY=FL&a}=3j04!467IQuqbpB+{vtS1DsCN57RNWUqTozYlh3!qoW-cpumciep(x z9?WA>Jm9B9`bGxiENepb8|3Un{BON|*>8V1~kX~=aN#F6` zn6C9Gh7G%k3uedX94oQ)Kh}3<)#d9G{zw&Nt6Eh15yrT*Cd$~rec{jA0B{@Cwd=jv z@`t`ghJ%NFUyBdYV3zq|NwMi++RdE9LrDD|HkfNT@5<2Akco9WU;f^m>SgfD_q&AG z)SJVG1LDHIhuEcOEzF;lJAFvc2^jg&2T}_e3!wdlb-?YdXb55%+`xh3z{8Ca4%IiUT+a5zvc)X-f~rEp z^~4;yO^*uO(Qv$y<*%$Q8k3WTExJ{Af8Lk<*+$)~C}F)9R1nX9Rz&h!_T60xZ%Ii})ha7ev*reYQ+%`ZQ9%<-W9dvqb;7^$b2lDfZSg0b=S- z<^8r%29khLi*&G)mW)WqaV!n}J|RCc|Mg0#XBvBA|F+8|t#df8fVcYG-qHzr5t7s| zzlwCsR^LASt$+A?zB<^Faeu44{u6Y~Uez4O(V8inRV>PGs7enXHn{w4&OJB*3hyJD zp2ml_w|L9Eau>IsI|^uukYtO!^cq}wOdE7zKrN*T6NMk_de)da$Xv+d3D6;+!X6Yb z^INo}kVs%LvKZdpFT-T#U^!!LO379wAlrvzR{tME*rbJ}U9K*!VZA2KfUB51e1=s= zDN@H}q^OttLRuGBM|`)%oQ!os$hLYqdGJ}dCUvqH;K+J`gjVhZ+e7Sjt(M|ABte)V zt3HY=A^ViZ!A`OsC+hZU51Y|r(atese=^3Kl5=`fZPD2 z2|q-O{myxBQEG#KlH(OQs(Tvl(}Q7m-$R-;wzk+TNLEq!;+^bBMHZ7!qu@4HZS|6LJ#h^uKGK45P;;h}#_$&rb8%`?8pN{%hUM(l0;! zc@q$T1$)WU#D$WyhyaDmSG|QgGF&Hmv`Z3WXENCwB6wzs#j)K5bz)>Khp^7U<-he#e?o$+xuk6mv<0 z+h1%rV>kU?@|dF&xn4_37ENx|zH<5W*RI1b)3LjUghNxWPM>`Idxz#vhrQ2f)ct+< z+AwP1`T6{0%c^1{@js1?j~Q@&0b{&nZJw;+_{Da-o*ZV0fD<;Uc3yLA8s~qN8BKb- zt0_HdjLHzRVZ}o{xvf{PSoQMCeJ<*k+y$|O_PUlaq3b&rw;vRFk6(fn@bzAG_qd*u z86ito@ACtTrX~HI*sjZkrxW8^mlI1(u4cqTykDO3xeZnzRIJE3kj~DwsXUy7@xxTN!zipu~^Q zD!1$YNC~uRU;Z+qpy(9AEMm7^t02<1OI{U0;C`p4gX)#p?#J9moyz=!FN9l#Ft#Is zY5!)qkQ4cap=I39ZPVRJ#ef`HDcbp~e9J<>?`*~d7cN1|&T&?g~@JuwbEF za0>OGxbk{!0lU@adz6HDTV|ulcdAp+uE*WeZU)*G9va_dbacVLKjS`+%Ph zZa;W#J9-N$JZK-PL3KegTIMih;s&%Za?QEgfb>a%=$&MxEFyE1v>HPV*}wPa-~#Rbr7r4)*T8wN`T_ zlM%ny0%_XEN|s8o)~W90kaPCsrcawFrHSA>KneZHecok~BOnMmy?Ec7SQ0M_5t!AQ z4srFoqJuFY8c|$ZO%?CcE9D3%puLEw#LcK{+wP~37dka&DFQiL!r0FUH680Sp`hE4 z?u^QE*X~tBE=0=j49GN@Z{ZnFLoui!kn(fFMx=~%i>kz+ZmfgPC#_N*(BK+6Hqv6_=|N> z#@xGmk&hy~y==MvT<8Zh^gMNb>acoR?O*+2&F!?9zY#@;)=~OvTPzd)4E3vLy{7dJ z`)^DHpQaoo1PDrlfPw{>&F>#{(0iBBVQ4z%P21y_JF1EpT)@Fy#;nLh0^fITDO~M3U?`Mw*E$Mw&MA@)l3bj zQV$Qui-am5W=0M}2ZQ)*y@m^Y|a4 z{4g)nkq%lNuE`}=fC<^I(y#WRYug+4Kas;^C!oMv*kF(3rA-N_F=xb7(`FYUpZih~ z@nPm4B-3z__+yC~luT%Y0NDLA^D=Ap)qUF!ucfjfX9qG(tHWF=Cb76sV4xH;9tL@^ z@WIlP4%IsudN3aK_uo~@N5343)cwivAo+~SaevV3C7Cuqy4CxObh;7mpV(!u+)0Ov z^>dD$s#N;^rXp28Y#ZdXjMx$xAlmzYwU)zK&2Sxiw|*k*hl6Aru)OjxasJ=^kQ{Vcm+`<#ReEenTkt#_iuWmdm6 zp&iCE$Xwy8w*xKDo{epd1pE1QuTj%P@46MPXDM^3`37rL!1#Nd+Tvu_+%f7$NZr9C zH#Vf7G0}JatZt>l$&#MCmRny|X%#VI08-dR<}x_HDw;hJy;Rlr+b7-BQ^;%G_)7S^ zF;4&EOoU8mech}A+)LwRVLL-KX~6@Y%dxXCQ5L!;r$KvQnOJmsay!?%v%P%Xz0Ou( zQy{Ht>6-Vdlz_XX;FmqoRHsbai64(x6ALdse+Q#GYadTY?#xEN2i65*?d)c~M&4$K z73a*6+=Ev^uZ}*%q#qEINri-sh=~9LFXA;lCyJq`DLU|DBOp*%cQyfxlv%mHG!Y zPk4FFzW45I>QO|DmpR0|v#rL;BA87KD7HyImYcWtvlvWdH`c{N_MUGmXlqH1&O6b2 zdlF;Mlp8#ZFqa{MfQ|ekfP3ktw*zLx-!`@dRVES^BAv8+D7`0F_UL)ACO(j-;`HUS zb?Jtk>}0?)J>>Y}|BJ15-0l@dBZMcj%ytG<;YYd7TC(u3{b7N9PK)l3f{L$o8cuK} z`=MaFI~&fPcjvP>gT8t4tQFOgOrBp8@Imrv%Kxs*9y&3Cbl1Zn+L)hM(r)l1$4g_` z(8ndTvGeX|hS4M7b(d%`l$jr6eDLudAw@_;2Ti^H=cFsC&zxi(-^*0;L|xQnN{~mO zG!vm|)xvSc-p+?D4pbL2h3=cUIZo_3m^lK7J$Ua**SyH?7R5-OZIPWZA%$U+FXgB+ zAw}=g>$O)rDnFCDHG!`!6<=UdTnwcn#bzn0ou6vzMz}REv-C_%)-ac|)mh&18GGv}uFH z;O;?=3S*CF!ZiDW!tjH)9Qey|1x6KnDv^mI5hg>j9{hrf9 zXuOtHa*LPSX<2PzM0I2*^sjJMTkMFJQWg_~EOV-umU1StmH`&(RQ-@p=3}abC zyLZOyJU^%$D;cFxvh>8cL8aYeDN@_}wd!oMYHeSq_&$AfgW|c$&&5A>{-4AL=w1z< zq3)^CEgH%MIu%_IHNJ{znXNJuDsxb14p`mr^;CwQn&JkO-?{?KVtpbPM$0y*%YY6y zTV*#}RVafi0K51CljpW&83_=psY;f%xyW2L&{e+e?881{sVOVT6YIzn=GR}K&Cy-#(VDuZzOa8|XD9$N@ERD4_BX9e;RhmX8%#S%(M!F;eH>H1i_JR#pZ_W7f<)ruKw2wg(6#?W38}_6iEdwedv7hMyf1rAx!F$wmjm8{U z@zTqgrbG_YHCkCN0S9M9)QhyIHr7W3$&wwvUiu3g?uu@UI0*p&EQohgmmrM=^)hEx6jFw=dWhx(z+i|(H2$CT$>c5D}hZ-{V z-Zy6sB!MSu;!X#7-XUd@?>f(a%^J>1t=<}-QMJ5OWEu9KEh zjyo5VP^!UFO$DBgY6xKe@wLjKBV&SHQ_GdPQ1iycGe^WFKvO>g>oh}~`L-P6Ef3J5 zBz%I7gxgA$U2*%{q$I)Ln91KG!ZDo1 zA?2^f%NT-swmd}sg(9y zh21M;Mnz;BK!Kg)yW2hg;pdN<->#laEkNT`9pZxiK7Tx%)9Cow;z;7NYrXN(Zo5c-PU^EyPrFZWA0{OumIVYQbK{BYA*#C6%?W`_vMZ42-oh`9tK+SqYctFI}viEGG{hc5! zLT+7U$t-=ZQ*VE8)Fc!iuP7h_2-wSLpH~ zrFa56A7X`Q@-V;YEe&M{Ha8`Qru5abte4Bk0_KPm4otE9W%puG98JHv7v2cwaFoqM zdz3zZ2w}XUN!rqC>NPqJQy&oo2?lMF;$ll?`fMxVgBTwa>A|7%Z1 zs%uFC)0i~}`#;EPHF_eZZ0?%jI zz(RJ3^DV>T!M7?G9?+Tn&RLkrOUfm#2VQwqqP@h&&~udLv}N9FeA-um^_p=}Ej)j> zM3cax%oY>-aGbr(VFGxP%szk3r4Ihi?D7idB*k4rl(#5B`^e&O)y$rs5*0nxr-xnF zcD#BS_n>;|<=W=6mh@bm6Mzx7rZaayN7I5-A*~$o-(k>t^UD6j$`p`$rNJ*zEZJ+V zRSf9ok0#vbr+}0Q33*7-{M$Q(cu464Z3W@)ce2O|xQoM`1m=Zf8WR3UmQ1t?y> zlYb4S7hY}>e6XC%z*$+gxi7W&?c)>6RNXWgS1Tx<^~mrZ59FIZ^{5Spi&`tKA5$w| zie~%{r|lb>JqcJquAYd{2fU)YNJIwg^~^4zoBu<#^>BEB>=}Ibiw?7Cw)uW_w_hhv z@JSvR(V`3D58b<-kKx>T{8!h-$NgiLY>G7((4~>G0`~lSphqIsgsPwVKBmAXvX8&w zI|1UlrZ7T85QDp9b2?R(gDs5lbY+4-mMW_w~p>k;c&<1F?K2LylEBh5DOM_=5bo`ABw0VNvlLrQs_mleq zq|!bKmSI9~H!f_GKJ#RU$cBH{&eAIxifStNy z2%}mO$^-AU%Sbt_a_2Q!dnKm08(dyaqVZ2 z{N;5CorB<>IWY}*cZT~LC~pDxls&-qePDE-eTlBG92!dfYCi=`N4I3_n`CDAwa7{X zB}bq`5gG5^zNdpRaC!-P;$?f^7vO+ZmG7k8`EZ|$Z%E;T8PFq-WD|&j2|NVgT&IxD zg;XR(9N6a;YxVVgMW5v#V8GY7KmjO;2r*5 z_NRY(*XmkmR$rI9{y-lQWD$h;LbY#Yhc`>os9f{Y4Pkeg7d~YMKI)rXo#8VZ&PKn4fXoGDOOc-&fZ8aL`V8XM z#IIbOoQVhLQns4jZ8~-0JjiWu8{Vam60qa%+qG}iEj^{C)Eisg6k~w?Xku=^pMs(4 zDY|C5*lfpcb-+wAVi+mV6wQQ(sR>4OKQ}gN7!1uutYr1^kRlZIUm*VWENiZS4nw;8 z$Jko@cMO+rSD+^ zHBe81M*YEbEMv$|E$(g7Wp{sF?%402o?Z5=#V3p9KRY1fD4y)Ghyyx@xg#SG*NuS> zcQdVkXWS-tcU;Wwfn>0P?5wH>X5_|1znc8u+HnUNPlw`Lo8=#9U~gW-%)1DT@-%wx z@AVb)f;xprx%zVGvTg2oUF7ow=|Dybxv>uItrbBSqYeBC&FcNJF7e#rX(@_*q76aP z`7uGby_A`p#HgqNf&d06l#vBD1L(N=foEBigfo`}LadN);4JB5i6UOOwJrOi1A*iF z8oOI@2^K+HNq;`sDcn@#(b}2JLJM(w5qiz_>=V=G{F^K$aQi1d(4K#HpnoqX$U%KZ zOv|9avT9j*oTo|fKlMnlJr@_(?XByWk7ANOhIB$7Xv}G^q~wE00h=hpgdJ5CZVL&{ zTD}kIoNX7&1&eP=E{Q(Ld*W#UtbR*X-oy-ovk>;&;2N>sIf(J|sh(^V5q@90OzR+6 zRTx@$+47XWsOD+tuAa4RhRF*W+sUV+lG2A*pT03I!lXC}fRwOr!y&E=mH_(Vfkz(i zfk$oeWGQ@F(H!hFa48jV`ipPa5_5)FWX|*9_p`hMmfS_BP?kOs#P>nEs@<$+#Mau@ zbiou8kcSKvB@J??&NphQT7h`4E7oDB#*W}}6EyM>Z3a$KYt5``OHGSDN@986ehbYgQvUdW0X{~~vE3NOQ zfp0bqrOEqS%`g|;fz|GLfgRu|CZY=~%GFl`V1_Fbq2k&ZRZ8A3B$Txi_Xn*+YrtPX z^S~4~cHc(AmI1s1I#_-c_H{m9$A@R@lalBzFng4YE-9z(YR+Td5$*LSbuAcsj@xi` zKyO?JdzxoOUfuJGQL0uhOTwr=yORec|BYMLwJRX8EJiL(U)4s-aO=gkSrX z|21}gjJ5FnaBqhXvkRokJu1+}RZW+pU47iO`}1sy3Wo?;Q(+LJ{yym5Ez~n|U?55D z!pede-2M*m?m)vP79LHcG;8rm!NX#|VXG}@G@- z=Y5HeEMd4;X+GE{u><8wFyt>BzD1_BRoGn(y$GPGuglK#stagVNs{mWb#`i>8)Q0FXKz4yF?*l&B&0w@ zOk0rLfT*fu&f|Ah(#+WM&fA-BCHYRTijTJKNQF#B$=MVUF+{$kxtg<67h%jF>?I|n zBZ(cAOj#;t82|M`?)^0NrHysH@Izw?F<3Ln6eoydsVjR4+kUWf;FX>A441^+1pUTn zYdV*#WFsViyOQ@1h`|vpN)$77!}14%ZGAb~%I|Oe+Bt4wJDn#>ZO=tQfaaZP$|7%R zH+7-51&_X@32?~`)RjuUNz0`H6*DS2HSKSH{Inn#GphhYYS`_f?x*yY8KefX&V9bs zNxEsQeme*{1?Ub|)z!hfAiRO|@T$VqP(~-Gr4Wcg;@Ex4{lTiYtSAZX=f=jNnjfN? zbD&Z%=$o9mTKLf2<6$Bv0IzQxPIcJX70(k6(;!0;7?s$)t&W-h;N|N)nme{$)4$RU znH9^AN3V*@Fd6*osE=i#ku|^|aY}TGFtGfeF9gd0fRb(Cf$lDVQiB{{x7|th0 z?hn{r)&5~ny%oYr@(1!hI%OsoecdwfT6?*}msWJ#vs^>bDbd1El>05~#8<4k+>d&HrLeo&jQUgR}tKSBd2E2d3p>^Vn)zgpAKf>)fH1x{a$ zI|3^(xjp614A6)_EA8;QxV7gTCTD!*UZru+Qk5!C#oL|ci#99ofb;^!$3$cUjxLG2 z)wYH=tvK8Zqn-G9zDzT2BV?Bc#LvLSdRsXW$Ts`Ea&mT{nRq`it z#TQ%J{hOna|8-hEP>WzxRp?`<|7H0?viEsi>Xw77_*8SjnT;gA^JhceCXvU%W4Sd-b#{S7%##(LB}vZq1rYRN4NZY2*&aO+|(uSOc3 z59l2~sn+TLHhu#?emj5OnX3HmH=1c0u#b&0c@fQ@?JMTi4x`SGqQ` z!waMx`nh#+Tfi9m%iPYuC|!5&6b1)UdNb^6+(A7HcIX?FKDRq74eECXG=2Np0V}pl zeR=c5^)sSKiKIn9jO`129p=fet-8UP1f&3A%N~JJp*IFUax=%tqZ1p1!rrQj{!TLi z7YmX)E&K8Aa2Co<5)PJ6UaT1IkN~r&w)7OcA96L?9C=wR^I&mMw_tB|^ZB}JviH_7 z$&8(V*4QUi{};=q+%A*BY6i6X$M>Ovq35;=Dd&%+cZwLy!Ap zIzc?l%@4^D;GAoU_rZ0gh2>u(i%4t-#!B#Vd>h9Q$Qd+a7h$mPKD4^T;Me`BUGb(2 zvA?M9AiH8`_CU%ek!Tl0b_XiAGfjpjSc51TIfugqGf9n%Xh`XS(mM^`FMza-fcAZV zVfOE?Bd_P#ps8EQT}sK(=5;nJHkZi=H$9}+8=Jf`Qk@Fp(r`4dqt?yARh6&tI3DV! zpk8im#Wh+2IXp3n$S%UT(xt<=&FcESTQF@6QZCyBhj@b0Z280zHP=gab5x%@3W9g; z@yOmq^bhpu-rg>!FjYJ3FIkRE$_9-YKrB3^r%vM8doSJA6OJI=d0XQBGDmOe_{H%X zRVF`-AJQJL`!^>r15P`{|3N>A+zp1aj4`?E5z z(D^loAA77M@O0(}cG{QR{yy%&@WiTtn)L62w(s+3Er}E1GtMrZ#~h|x;jN>$oUm<; zDD>fC(>siC0T^s|~1ros@8X?-6sJ2A6UtBt#b)va+2s>9VEKT(1o~ z{;fM?D3%UyyhD=2Mze!CvkY+RIM`chL9#6b|04^ZRm80aHwmPW1^W_pG)nouejURV zn=dt&&?71ZxsPzxIU^9>WDL*&S6x z2}O*LN8NJ7Q{?uK2^}kq9T}ShZ6rFu?ne;GG9PSRkXdgRCYB&*6(gNFCu=AdTw=Y8 z#>qUCS*L@`9BJ`yRlD+mZWAm+k_rnT4Fke6-Zyab{{yZ@5}SL3cv><$_iI+e?E7B? za~D|abD1Yo8PpHnh_h(asd_@*1SKwVg(<+PcSB|H%&m6GL~}h%LVSC1vq}d0{jHc@ zUF1ECLFc?(NA7Z|Ed#QYLD-8ny&(kFqDTMBvDn@IBlawV1i5!0Su8OcMMNR?9Sl~VwG;+rfwqU}H9O3kCOra0dHZ=CJPvC-rGMV4KJ#@lXB zLJ~H)@D5Z37$p(1Y%sg<4tyto=h|_tI>|cW_ufJ1n8+KDdWs)eHle*CNCuskBV!5V z{DQod4Qi&(s;+mwW#DM`#yk8&dH62hJeK`C>)DeomAgKu?D+$h$upZu+~2|Sx~FIcfcZZ~9}poUF*DvBI0OwK zvNCtJU40}IKi=G?E7s_<>ELoW7e3U-+5oA>^pFk{CLg!n++B->ZDyR!!#?2#KbGIHrw*? zAf#+T(@A?%V5kp4*_D2aL?l|af850I4%&F$x}3Oly471_5YRz0At%LINyG>#99SE^ z&Y~+)BcKUMW4BnAY%A7zzE@A$wGL*=0H*iuDpDTDt+$^A%@bm#f(y`wDpT0v>55^_ zd2RYy`M4jT97GPwp%9?e()xQ|2Qxp%)J%~h7R!dJFw_eHM~vsze+WIx_|?*6vAD+{ znv@;iX7Sj;G)fmoRMb_2Aro`%dKN4&0Ie;@UoAUczgwWw-;957!~*dAV1VqiudoKi!cJ z#**1#)m-u@RFk>x^>Iv(1gF`f8VXqh{z*$#yW0^lim4r=?r`tbFx&)PA97VgR0HT% zM8I>fNt87mzg_LRTmmqA0O@WBU66T$10#4l%qk6sD$~8Jn*un!=<9TaArQ|f#Am4( zp&5qIXAo}w0XsKlBmoZ~r>l8tf3Q2(g-<2k)B{-}FQgYs#b|8fNsm5`o)eofvC?Or5*qYc_DSiqL-zm|^O{Hz zjGA3)z0{#;E|uy|MI@+uD)r6B$Lz;bcrWDB?N%hTo}=hBxwfo(ImJ4LkI{om)q~?_ z0RMVS_ZW;)2Md-Ang}V?Y1P^Z5$EX&2cQzk98e*Jy)O{UY%j6S(k?FRa}x8k*Pm z9_kSd&@|-EVdF0>k6U{tpKZ>xBZIP79?cab8L)DJivpW6HAsPy0`@kS|xS5o{P~}>3oSpew4G~O@G!RX)-CG3>y1;c74G*7@rt~IO6=F z$}c0^{{^U-s)=-FolW-PJF=7n_a4bQu=dE#*#?Qqi6)cSLr3?akTU~|h1A>|epjjR zD#r?Un25j3`^H-er5&31ecX1dSuI_N!*icZ1r#dGKR*sV%U z`ZU(v8jrQ@nk1x7dt&Cg-to)Np8>%<>qe7vm3WCd z0k*jbdnydrRC?8eOwi&nsfA)>WJ886adhKq;6eF6FNV*BH%|>Lt{Ca4RE|VnKYo_y z478;{OeCa|&U&Q6%8{8mb>|%~57r{`6J<qcNclRH2X0JomXRszQ1~HS5*2q7v#b zp~@?pin)r4`$C~6e6(4yfkf2(cl@1;UTGgi1qOhsK(BVAE~NKMuY%KGQ&nn)`^|cGtZA4@>29Q0x1I(V? zLL5&y z>z%~voYuppKW^^y{m1p!d`9h3m2FR$CwgGOWD`7TVWy%8@BVS9H$je%hll(zfvoJO zKub9~WGV7roz_bvMW&$7eHQ(rh$htQSp+w@W1|*446E&!K8keY@Sg$#KB$s!9(;b1 zn?u8L0kS@7oCl>?dw3=FqdOF6ui{75P3@{bG(3l#39?T?@9WfKDo{y6*Qyb`ts?dR zalC@q2n`1K4~`*pQI@JegPg4qbzLp0##I0*&*V^o6JmK48E4@ws)O_kw4#d=uPHld zg3R^Q39+%Kk=>=rk9fa8Xw4@8_cfuU0xaRcl;?qz)wODx0JNP&NYERkCJ)j7vl(yt z1XNdi7Q2%Vq+%>+gZ_Beb3A(8m3AYj!njILEWB@ZK(m}QmDF@7mB{X%sSyx#) zzu+@3)~OTphJ!6p`&AR9UF1I=rJFB^1}oe=bQZauE6LD!=#ZW)B*qbPoy8Vm=BM@2+(3-IRm-W>)(=Iracitg=;O z%a}n(l-TL{EJ)Me&Cw()r3fge=+;Z_Uw{{3JX*a?RRRJzsf{dp+N+F~z?xV6e6dNN zKIo;$z?L@)1afF>DJAz)JLrJM6y7Vy+nlpWX6{Gzj>}W3ec8*_NPN-ua8o`m?5}E3 zNfj&oS!lwTa8kLke4St41u-Zg;Kx8-NCTr2PPWAJ?Qn^#A4hGduwXk!vA+&N5HJ&C z+>cngUYB$^6`F$_w21Lk@#OujdN&CimM%wC7uIwG>s)20WOfdOcdG)rGCV*Ew~5Gv&)!hpP(yj3VU- zk2w1B6Cwzmci&|P;#ORPe+lnm{H)pQO6+UNln((%pMVJTq)M4DH!x|hUNtl~2ZHgg z%Oenrg9;f9R*h9V2Fa|ddIngwu|GUHT?oiArgj@_kZrztJAxY@xH*C!B!2@k^Z?>`F5|N!lud%I$K*1h(T=te8yms2TKxoxDfg}wh7?%C zLBX;0*&06z2CS}*F%hm7_Y=LevW89}SSLrvsr3%X`%~$f;`IH*Sf5l_1QK!3+A#7R zogjD)*E5bBEU(;{0Oc+;%jR!KUiIfyFye4l_dgjVBvmldA-{inkf&ha*b(fFRy!o` z88#6lrhElVRm;UNiUqt&zyeCDOPzqEHY$$dl#Wid2Cnf~>E5w0ivPA4j6KCmG!mp; zee??8sS{Fx=HFs27BhibPl^JNzqGNq5jONrjQKIOcnxK?a9Q?0m+x4;`XSJ<&uJWR3?T3!%R-_%JR0 zf1WxUD+S3>NJ-@031F8a!qdbLfbZ`rx|y$fY+_vsri6ijuYF030{TZM#PZG_nnTtP zwq3t61<8Q%Xfv8RSQqWr{;u5O+ujQk?xXDvd-qNv=hwMnQz;$8v&Kz=x;%@JBVYQ~ z9S&09o`zCVZm_VkDm^*m0>`+o#1x)dt>yQRf{eQd|GS+Y=-=$dcj(R>-{9@-xceu)L#NFzIp3L zak0F&Q+&JsXR$I2A*Zk^HCl6{o%euD2PYLa`+v3e<>65F@7pSot)h)k*6eFj*&<6A z3_@fWQnnF_jN-9{?8z=$V_&k>*g~jBWtj>?gyea$&m_-deT?wCZhd>c?;r1b9KS!_ z4u5nw#yvCl{rP;Z>pHLVJg>W;fb?VW3*PbThLgK&V5602(JJ~0oP0Z~nC3${^fNq1 z&%z_>_r~<=6!e8mC}EBSo%yp+x&yRfs*e01eH%GC=ld1hphUG-?=#vQTo{5#ByBTp z9|o%XpG)yQ5Q<=*a9uG4(Ko5Z!lA$u!^0#<dsbs^}Z=|3RV)@rnd34vp=&e zc*nj9IIg;+WNyClk@{WHy$BqM4V%PgU+tM=Qix)GP7kUR)ofc92v0!^TsyU>upIX5 z(4mfAJcbMB;600KA#q7MCk`@)fp+=b-AN8VfO>Q;4Ghho-)Km6gqfC@i1qe%$)}(H z@0-@D0d9=}vA`ilS26bX_*eJr^p8aulii2jI3zIj`CohV?~`Z^r1ulwBJNhjL#3#{ z(*rp01DAN4L~+LlKC|Nsn2!m zLN@DT)*56Ic93Vp7n0S1Ynq`vi`xc8p&&#>fuO$YAd$q*5>Sdh{ormmfOckUGXj(O z;xWrwNF;+QEEqEACwEPZjrThaFu?*hdvy%RRArSMl)*R9wQ@qtFAP*_#kR~8(hUgD zt{ADCY{!8O!GL5QKVDg2;lKoADwRG9^+S;vD0tDm;60C+3kFv0fO?pzM%$_~LuenN zZR~lj4MZh-b&o5miKlEP)onBX8&dN*!6gT>Fhg?mrMF0Z{PV8-_e_28w8Hji0SKuf zUotk`FkPi=BhC@7DPX$^-xf{}ZYEHu^?J*)A4hH~Y8Cv=CzNjC&UJmUl}%@+U*~(< z?#`i_z{&WkR|CmfTu<%fH`hT~xQpS*M@i6^R9f=h)eGD+A2Pytn6t{z%JBcJ4Y zq%x#xt`&2EJduel!(Q2W4@_i4b{soUR?PwU0-Kdc(5gAOxf$4WS1B~tOh1!mLoEMyF#x2j3sL9wYer;a{{~~(8GDwH z(dB0BhAVjr@>FCXJbDOlbFWu>kNwdE^Uh%YHWN?*)`Y&B zdBg9J)Cv0cKnK7N$4}J0h(x~gEIjgj1^8ZwZedpxgAPalAR~y$7Y6|tXGkaNy?9K% zqdEk4B)W=dL?X5hCj!?&-0;!uj9)bv)VVw=<{7Tfc_&5Vlpx|{(+kA6j35gdy!_<` zx2^@WTzmS4xOp{TmosmoKp*z*N}s!&&VO26_R9|2nUqt73ssQX!342BK;Ud zmeyN|J`Z=S!$kae3Z1}L286kUb!!&>@{f`TQK>;CAaz11gtWgHf0mHf3ubQ4W^mmS zFt`SLE-3!;jrvgzqrQv7RC+lLnIeHYglzA8WppSj*Dcd#>wVB4*k;QxiTLy_VvEOA zT});T;(Vc9iotHSsPTd-$Ruu1oBmYV4+J;Rxw2lhD8CG{h`dwHNTf7qN6yao%5I_> zmELl%7|;J%h)p{xKC9RK^8Rtt*`GchaLPev@j|wUdFcKf)o};!B+yu+#y)VlNO)6!{wS0x{kz*(^6p7ws>~x(J_OI+-ETrzW8= z9m3fn3W`x_bKFVm$0HAPk_@^RVHrg^cWVt1a4iY2qWE@-m))rc#s`=?V1)5WeSL0f%vhD zmkN&|1rpr-AEMfNhGY%cQfc?Xe%kCA(k`%+l1+!)ZUrU4O;jgFdW>QY2s3qoL{9Wt zqwRP9*QHQTm47@aDI=Fx8}ZFoxgg|M8Y2#1WRh?vVi6AJLd+xQO?Go9eS3aIIv+X= zkbJ{``dnX(tg0Z5v2j_u}z9A%1p0-s;jZU z!w4D~voKf1OAv1ES#xN zKHAvQcTB(jAitik1H3A^eO8{~*qFr|k1{~;O;|kJ%12-6rxp%P$3STLc zIRIj;IZCAnyEY;Rw`5a%7~p(Ph7VRPj32nM9pT00VJ$|`K)mB_h*;fMY3Wzbmi_mg z;zuq&ioD3d)%lcS%B5)Vftlwl7JGFEV1SqhktrcgcW+S^g48%sO$T7w3xmFHY?5EN zkKYBg<$5H|$tmp;{0-(&z&k?N1~3W`Cuy?5YFMJ;K*RAk*M zbModZk%*pjPnLx^ARO2`9xg~GoKMz9QVuvpp3#A#nqKwNZ+`1@ZV(ZR34?H%_jIq zM{wN9%~?H6HBi^F#7Z&Bfc^ZmfM$3tYfWbRD|vRL!(eic{n@Zoi_Lw5!RwHR%^pnY z6rwQ|m1{S^1NW9+&d&-|`Zr^|1QTHfUT|=Qa~IaKRzs5yP>Lz(5H1E&z3UhFvJ+w& zV44J~AC;x+b|I_hy<`J0bT`_OwBg2+zcyOOYvU@W<{T&4c$D*b9So%I`x0e>BfF!x z=U4v}_!mpLV9WDtVT;KG20&e2rV=1ePDD&J%^g;|>ba+^mlR zx4>{H9OyE+j43 z`FC9D7hnUtmI? z5}3C+a8uvNW6b3;`ymz-neAc#L*_tl7`p-i0P*c=m!BR}&tjW=&a4@N`1XCb86*^) zX5SOxF#G1T+m(%ws3@i|86>fEM>~N*7f2bwFp730H{;Y1S;T&>3|t9-Q~+@lJi{Eg zA$oSsEa5II;(cH@jju5Wecpkj;rBa%Rj#oW!1bNO49Zu!jf<_^j5;Slj$Iq7*Bxx+ z55LFig=dA2XP>dhx;BBv~YYD7bLzqF6QXRW$0w)4-?EV^9eSs>42D z(D~?gnRpm1?>gn<2JmA~J#6HTa}N;;!2kJJk0CwiEao&~2j9LEArN1Q&rHN;f?Xvh%z+!mKEY4d&2Iu#>!`)i6kmR@r~P@^N8i`jFejl$3HT)M zeM^q^yOhtBp7k2++E>!{j@Fwun*dnXE>HK5ISg5R{dsn6^;hvl^E00Dfz>f@^6dt$ z9UQxL*F9P{8CK*^iZGj3Y`!KNf-$^8y&v$r1c`^OY(u~rK^u$j{k8)SOiQ0Z2IUyu zmtk)&%qaaa4(7FU5HnlkD&N7$U97tuhoZQnJnXDZ2!SNLGcK>v`hf!9<{+XO7qZdLMRVQyJlD0{!zcpWHo677~*{$bP@X{Vq}%MKj0K%r)@oH1DNWQ>Z?ht~XFB;OL%x~sbv@JB zJe?IEDuKD+F={-W=*dh~TMS;e$t*Rj@v4)+x7a;|OHeZf#yE}o)n z#gJYat;v&OWIJf+aa3{v+RqqTf0f~Vcbn$vBL5o8mU7*ZtptyA?7}d*s_&zfeGcl? zP}W4Z8cDh(qryEMH8->t{bGnPR_wRDS&%$o4@zbzEC#g)?d#E7)!#?m>PCFxYd0#O zq|Iht3|=`;cPY%t!Y^-1&pi4z4hZbGW7@Ri!n8ck$Zrol3jujA>GJQ|=CqIXK2ef$ z&bq^LwbA|hxeeQ@BLw?yI%08|eS_Yl!!wg&fD;?;m>-OGtx8;|&`i2%J zj@tQi6Ra(_R3~@ot3IZlEjwbLA!g>OfBzfNcf;Dg!QbHOMfK{TB{2^T8ab)QeCU#W zN$rgS{szoT!R{-O`lAbrcVKq%g1S|@@WurW)>(4zT2br6ft7(FiJ0_~lWh5oZdX#h zH_tskMqP*us3sF^tU~FFs)gGToeL2GZ7VhV2l%0l0TmF&W49Z(uA8~rT#YBUu42Wc zJpIX0DC-@hR-A_^ju!SV;!6gH`fB|~=A2ImlD|)##~by#W79lF_cil*KjZH@BZ%XJ zB5af(1H2|oI7hoHmf5eGk9U2PFjmOZ3L~37x_M+BHZ%|1`~>XF0S=*gCBr*uLEHxB zE}H>~-&r5%7SQHr^+r@6COfu8I(#CH_b)f~+#i`7k&i3vX$Y``sYmfQN0dv6=2U_t z&eSmn^qbAm=lYabzE;^%a`f9gi^11#7XzfgE9el))H_ner07@<>?isM1unb! zla4=#P-ih!x8Ge)Z2}Sf7?fXr=3NLewwA$uwHK_jh zyYE9@Lff6&tqr5W)*CMlDSt65>%^wl-ZtQ{EW8{=$BE^8U3*)g(@98P%@K!zEig=E zHEjq+iRAW-6pkrA_(w40hrNERGLXH>Ed1#~GDl-E!na=Inf$i`bHDH8G)s2aH<~Bk zYgXPWd_T^?L%*RixjFCIOL^qA7U6s3_W0WA$oG7QSEe}wzZ(uS`BB_qC5j#~I2Rf8 zIDz+Xb1wpayZpOvp2DN41b6>N!4rq^2H%v7$-cF!$0*C*6^?Cd$6r8bDZ<}w4t#v^ zJkyT3cs(u2*VxZ%O|)G;=|Z?9P3h44i;Mi3M+3~y|M;t2rdrw?Lmjd_7y7(8vR0is zy>VH+gxFKCwFApvmz#gv){jXvXjxB)iBx;vCj}D7Q>iv9fx4eLgpGKh9%Mz#oZSmqM zvN}L~4qQ&E1DkTRJf>YD^0~}W;V;_29<;@3(>8BG)i2c9`sjdcO}=uJE0(rc_{?uzFe^<&6G^;N@pM{o#N6 zN#5JD&Zq3Ma{Ge4hTlWm8;irxb0EDs#U6p(T83|Vr!gVn9 z##2!2y8}5faBpY7dU^j%c}&>ji-4oK+Lme$pwqV9)Ru`n2lT`*E>0$ z=c3vTt>0C;x;@^u@9W;KQS>2`K2DM8WS1`AvQgPlhYs!Q&Vdo%4h-+Nj?d)8O(!ZU z$CJMKZm5~*S|{`*1Y#>K3%u2A^IIQ&T%*!JP-}1fIFJ1m55C=|dn;jm#g9PUidIY$ zb)s*5C#dG%wHfH?P3I|Mk1aDX*Dy&d3aYP2ZW^OD<;UjzGwGiYSP-l`@TRHrtYeNj zwXVE6+4d0qx2S2ZMpmbTmd~rCD%Kr~*xbVpXXRZTK0jgciWS1I=b9v28aiaVyfM`+*f(+1a7T|0bP5r4 z1+58cANh&C5RtYx&IV`JuRmnlKb_tE^FTA?P;lcSp;+_)q%v9WGY;Ux-A}BI|G0&Y_iVH`pa1rBI4JBGcf$w#e`WDM0oATvYsvHL;CKt4n<3=dvrnk=^PBi(r})!6 zK}P)T0ZopN6nJO>h$|mYdka5LmB*(I-SZ&Lb5YINrYkP5ugq!jQ)C~s@+F?xc88oC z`g<0V4+i`RLbUgS_-NO|PwxY6_n>a4y#xinWNM5$4;0c9?#!ADfb0 zod9PgKy;v?5)8%Ma;q2kMv_b*omT7VcbT(3$OzC?70uaFWJU;tEgL~rqkN3JfP(;S zDMi^p$kWw{3(xJm#{e*USKbPHzED)IA6w+t)j_fD^wHniHBWz_qIP6no$=O=(goOf4e~>iPI< z)J~L|51Evb$^IwlY8oSZIVPVw?HvIh)DN*Jux!ZIdjib@XcdCnntPl zav7&)5)HbW$I`8d!c!nFGP0UUY+LC`o2Sei^WG}Lgslb%!Cj<=9;@t9L@Fr6m z00d%zUto~ygc4NYWiQK2Ha`Xpxb1~h8CXb#TB4gn9PJ?h8Z1&eU3~qEy0FEwC|v0% zfOnbo7rXg7(DVr!?a=GMcMp^z0ZlZKyh3R0idiRayq9MOmdewImaL176pCi&y1n?}q$ML<{jnL_ z#M_1%Ym?jv)ijIV%IW=OcB#`ZV}(a+^`5_x)A(%1aO|GF_qFnzGQ770)VC*X4Qxt4 z8x0Z$Kriq*IyJ@%*TGL5KLlib%Ud}(6D%J0E`L~QX% z8$aUK2|_dUSX54U|M=MUaV6EmLc>a34;TM!&3DEK*#;-x%;=T;U{S=?TaDeJD)g>J zN=iKp1uOF<#phL70gnSL2kG5gX`4=?zTPslQT{ihRkAkYCGduoSwY`E`qeDzl~O*O zcOAMlSkJQgNLOO~$b4E^bwf}KV5tMw3m+H2SB*wW_$Yq?lcgz^XKMFkfKQsC*7Z#% z(}aA;z|VVxcIeMb`RbdGv$Gjf1`knQ*pJ)gpytxIg*}xMFQ(;94=1dZIV}cTpG3GM z`%}DOP$zR~QxFqm-LQ_!xSvmws(p{C_` zrqsO2kOTtjS5I#1(L--#4mqyte%8{&K6O;U9Gdc0Zj7)U^1v&E-yJWfJ+7Clc;0-V zHc)De(^12~ka=_SN$ou9t0mK7t0w#TGC6bHO?+i;b%QMoMym6<*|!b!E6fgL$Pr=@ z4ghG8*&MZ418M%gDa;2ZAe4)MNv~MnyE|<8IWTuj3WJ?Oe|%3jw&K3iVcO6DpvL7t zR4ZqE3T43b0k|V_MOon%)MW$wh&Pj_z31S1;?ci4^Y7dM?p>MIXguRX$8kuf>c=ER zaz)0D3>efJ848ce?1}*;KtIYRJ{Ag-lqxBevL7WlW-cYlnp5?zFu;BR^7>Pa8P0!E za*C{==?B#nzJ!o;zIIa@1(u>4a_k5lXPC2>Jrd?wY9iSzcx2hQ8f7NYr(v?l`p%;0 zV)cvIibV(3iLGUj&7m!OMLotS-AwQ0<0k)EEkBoq%0b?VbU(bff1fqVx1i22SRIUg1ToEQ%a z&RAi~!==3+cPnn*O^NWO3CyH2s*3dH8xO4&1^cAB`cX-Auxg9{<`?qU4RwGe%`cGnrHoIg9-15(Y`_I~ z-N~3x0+!?H`k>7nSOtOz@0O>hxW(9qYAc3Z3!{%>zx&_U3B(GeiZ$L#D9M_ayxcFC z3Il+xP3o*0E40sTXyPoJ*FVmm`MO25Q`>zXB5|ecz-EDg)FFI7rJGVZnrmh3R0y`5 ze~Jqy&sOTFpl#B9$qDO4q;cAMp8WyRc+KwJcp3;-iMlZRZI8LJn!nERy1G!hwKb7f z@`&58WtlBm!JxTyAg0k`&?%?zW7v_Iy9t3Abta_JQGtzuF~v9~lkgXln~S0KUr0tP z;{qRC1Znerb|e$gs}f|>wMhRHo$H1J)eaVJQz~iSwH9p-zkYQ|B6%TK6AT;X&y-7b ztUEb}ayHb5QHqSQvQIXksUvc)Cu3ibN!EXIZKMbuVCZs}<5?5k)0Wb}J9VLSK>@ld z=b^`{$S!PMto}9pg^A({&HAsY7OFRH2-Cf?9;vCT4xoi1rp76Jw{^NPz#@i7)f^W0 z2@1a~UA~0CqVL+S!m^#}$$?j@hq6}7TdOXqVebXRze~0(u6n3KtcZIKmcg(SA$v#P4$+crkJbNvy$(Y*YjOLZmF5SSkopNmQe*g~Xfc+KE zCk1+T2w(2pUU8b+=mX;Kb%@X|@1>Br}v3Kr9_<%bp1n8qE*x zh>P*vz_wWT_j z>6@*mO+KP3hUba8rL;L+v4M2JuU!dUDcz33=zpE|+`~+bIx-gj?lf< z-ua|IqEpJq(h}|3wYpGFujddy-rQNi9agWW} zHRwO-msPB3V_9lPO&Yo$T^m68GXOKJQOrr$k4#r}@FK8L-yj2rK)yE8~0_$jP|f7L-Z* zh)N%zt?3)yIp9}BxRJoxmRSPR*)NnmHEjG*j?69_n9}ll2t*+#oNflh3urNUNc~N9 zajG*j?eZ@r(qM)GG#<^uy@vJ$4sZzWVj{FSYWO6d%TKmFd@<2Z|7VHGSVX}B!1g%l zK&Z6+88F_IgQ}k|^p;s?T6RK*>cC_H5_^{u3=e{tFg1dpNQkdO~+o9}u!knOtHA$$4Kcu*>%%yQs7v<9p-c9I9=bol?% ct=igNQcg_rD6WQP6$1nMqOneu*45De0Y}m-3jhEB literal 0 HcmV?d00001 diff --git a/docs/GridKit/Model/PhasorDynamics/Governor/HYGOV/README.md b/docs/GridKit/Model/PhasorDynamics/Governor/HYGOV/README.md new file mode 100644 index 000000000..0e4dd411a --- /dev/null +++ b/docs/GridKit/Model/PhasorDynamics/Governor/HYGOV/README.md @@ -0,0 +1,6 @@ +# HYGOV + +```{include} ../../../../../../GridKit/Model/PhasorDynamics/Governor/HYGOV/README.md +:start-line: 1 +:relative-images: +``` diff --git a/docs/GridKit/Model/PhasorDynamics/Governor/README.md b/docs/GridKit/Model/PhasorDynamics/Governor/README.md index fe3d5b5b7..4e38abd24 100644 --- a/docs/GridKit/Model/PhasorDynamics/Governor/README.md +++ b/docs/GridKit/Model/PhasorDynamics/Governor/README.md @@ -8,6 +8,7 @@ TGOV1 IEEEG1 GGOV1 +HYGOV ``` ```{include} ../../../../../GridKit/Model/PhasorDynamics/Governor/README.md diff --git a/tests/UnitTests/Math/SmoothnessIndicatorTests.hpp b/tests/UnitTests/Math/SmoothnessIndicatorTests.hpp index 1ce2d8b1f..9c85c6388 100644 --- a/tests/UnitTests/Math/SmoothnessIndicatorTests.hpp +++ b/tests/UnitTests/Math/SmoothnessIndicatorTests.hpp @@ -55,7 +55,31 @@ namespace GridKit return success.report(__func__); } - TestOutcome deadband() + TestOutcome deadband1() + { + TestStatus success = true; + + const ScalarT lower = scalar(-0.05); + const ScalarT upper = scalar(0.10); + const ScalarT tol = scalar(kSmoothTolerance); + + // Type 1 (no-offset) deadband returns ~0 inside the band ... + const ScalarT inside_input = scalar(0.02); + success *= (std::abs(Math::deadband1(inside_input, lower, upper)) < tol * tol); + + // ... and passes the input through unchanged (no offset) outside the band. + const ScalarT far_above_input = scalar(4.0); + const ScalarT far_below_input = scalar(-4.0); + + success *= std::isfinite(Math::deadband1(far_above_input, lower, upper)); + success *= within(Math::deadband1(far_above_input, lower, upper), far_above_input, tol); + success *= std::isfinite(Math::deadband1(far_below_input, lower, upper)); + success *= within(Math::deadband1(far_below_input, lower, upper), far_below_input, tol); + + return success.report(__func__); + } + + TestOutcome deadband2() { TestStatus success = true; @@ -69,12 +93,12 @@ namespace GridKit const ScalarT expected_below = below_input - lower; const ScalarT expected_above = above_input - upper; - success *= within(Math::deadband(below_input, lower, upper), expected_below, tol); - success *= (std::abs(Math::deadband(inside_input, lower, upper)) < tol * tol); - success *= within(Math::deadband(above_input, lower, upper), expected_above, tol); + success *= within(Math::deadband2(below_input, lower, upper), expected_below, tol); + success *= (std::abs(Math::deadband2(inside_input, lower, upper)) < tol * tol); + success *= within(Math::deadband2(above_input, lower, upper), expected_above, tol); - const ScalarT lower_breakpoint = Math::deadband(lower, lower, upper); - const ScalarT upper_breakpoint = Math::deadband(upper, lower, upper); + const ScalarT lower_breakpoint = Math::deadband2(lower, lower, upper); + const ScalarT upper_breakpoint = Math::deadband2(upper, lower, upper); success *= (lower_breakpoint < scalar(0.0)); success *= (upper_breakpoint > scalar(0.0)); @@ -82,7 +106,7 @@ namespace GridKit success *= (std::abs(upper_breakpoint) < tol); const ScalarT x = scalar(-0.4); - success *= (std::abs(Math::deadband(x, lower, upper) + success *= (std::abs(Math::deadband2(x, lower, upper) - (x - Math::clamp(x, lower, upper))) < scalar(kRoundoffTolerance)); @@ -91,17 +115,17 @@ namespace GridKit const ScalarT expected_far_above = far_above_input - upper; const ScalarT expected_far_below = far_below_input - lower; - success *= std::isfinite(Math::deadband(far_above_input, lower, upper)); - success *= within(Math::deadband(far_above_input, lower, upper), expected_far_above, tol); - success *= std::isfinite(Math::deadband(far_below_input, lower, upper)); - success *= within(Math::deadband(far_below_input, lower, upper), expected_far_below, tol); + success *= std::isfinite(Math::deadband2(far_above_input, lower, upper)); + success *= within(Math::deadband2(far_above_input, lower, upper), expected_far_above, tol); + success *= std::isfinite(Math::deadband2(far_below_input, lower, upper)); + success *= within(Math::deadband2(far_below_input, lower, upper), expected_far_below, tol); const ScalarT point = scalar(0.25); const ScalarT above_point = scalar(0.75); const ScalarT below_point = scalar(-0.25); - success *= (std::abs(Math::deadband(above_point, point, point) - (above_point - point)) + success *= (std::abs(Math::deadband2(above_point, point, point) - (above_point - point)) < scalar(kRoundoffTolerance)); - success *= (std::abs(Math::deadband(below_point, point, point) - (below_point - point)) + success *= (std::abs(Math::deadband2(below_point, point, point) - (below_point - point)) < scalar(kRoundoffTolerance)); return success.report(__func__); diff --git a/tests/UnitTests/Math/runSmoothnessIndicatorTests.cpp b/tests/UnitTests/Math/runSmoothnessIndicatorTests.cpp index 233c36ddd..aecdc7ae8 100644 --- a/tests/UnitTests/Math/runSmoothnessIndicatorTests.cpp +++ b/tests/UnitTests/Math/runSmoothnessIndicatorTests.cpp @@ -7,7 +7,8 @@ int main() GridKit::Testing::SmoothnessIndicatorTests test; result += test.clamp(); - result += test.deadband(); + result += test.deadband1(); + result += test.deadband2(); result += test.limitIndicators(); result += test.slew(); result += test.linseg(); diff --git a/tests/UnitTests/PhasorDynamics/CMakeLists.txt b/tests/UnitTests/PhasorDynamics/CMakeLists.txt index 3d2973b62..42d64bbf2 100644 --- a/tests/UnitTests/PhasorDynamics/CMakeLists.txt +++ b/tests/UnitTests/PhasorDynamics/CMakeLists.txt @@ -88,6 +88,14 @@ target_link_libraries( GridKit::phasor_dynamics_bus_dependency_tracking GridKit::testing) +add_executable(test_phasor_governor_hygov runGovernorHygovTests.cpp) +target_link_libraries( + test_phasor_governor_hygov + GridKit::definitions + GridKit::phasor_dynamics_systemmodel + GridKit::phasor_dynamics_systemmodel_dependency_tracking + GridKit::testing) + add_executable(test_phasor_exciter_ieeet1 runExciterIeeet1Tests.cpp) target_link_libraries( test_phasor_exciter_ieeet1 @@ -146,6 +154,7 @@ add_test(NAME PhasorDynamicsBusToSignalAdapterTest COMMAND test_phasor_bustosign add_test(NAME PhasorDynamicsBranchTest COMMAND test_phasor_branch) add_test(NAME PhasorDynamicsGenrouTest COMMAND test_phasor_genrou) add_test(NAME PhasorDynamicsGovernorTgov1Test COMMAND test_phasor_governor_tgov1) +add_test(NAME PhasorDynamicsGovernorHygovTest COMMAND test_phasor_governor_hygov) add_test(NAME PhasorDynamicsExciterIeeet1Test COMMAND test_phasor_exciter_ieeet1) add_test(NAME PhasorDynamicsGensalTest COMMAND test_phasor_gensal) add_test(NAME PhasorDynamicsExciterSexsPtiTest COMMAND test_phasor_exciter_sexspti) @@ -168,6 +177,7 @@ install( test_phasor_loadzip test_phasor_genrou test_phasor_governor_tgov1 + test_phasor_governor_hygov test_phasor_exciter_ieeet1 test_phasor_gensal test_phasor_exciter_sexspti diff --git a/tests/UnitTests/PhasorDynamics/GovernorHygovTests.hpp b/tests/UnitTests/PhasorDynamics/GovernorHygovTests.hpp new file mode 100644 index 000000000..261298b47 --- /dev/null +++ b/tests/UnitTests/PhasorDynamics/GovernorHygovTests.hpp @@ -0,0 +1,415 @@ +#pragma once + +#include +#include +#include +#include + +#include +#include +#include +#include +#include +#include +#include + +namespace GridKit +{ + namespace Testing + { + template + class GovernorHygovTests + { + public: + using ScalarT = scalar_type; + using IdxT = index_type; + using RealT = typename PhasorDynamics::Component::RealT; + + static constexpr ScalarT kTol = static_cast(1.0e-8); + + TestOutcome constructionAndValidation() + { + using Params = PhasorDynamics::Governor::HygovParameters; + + TestStatus success = true; + + PhasorDynamics::Governor::Hygov hygov(makeHygovData()); + success *= (hygov.size() + == static_cast( + PhasorDynamics::Governor::HygovInternalVariables::MAXIMUM)); + success *= (hygov.getMonitor() != nullptr); + success *= (hygov.verify() == 0); + + auto source_default_hygov = makeHygovSourceDefaultData(); + PhasorDynamics::Governor::Hygov source_default_hygov_model( + source_default_hygov); + success *= (source_default_hygov_model.verify() == 0); + + auto bad_hygov_curve = makeHygovData(); + bad_hygov_curve.parameters[Params::Gv2] = static_cast(0.2); + PhasorDynamics::Governor::Hygov bad_hygov_curve_model(bad_hygov_curve); + success *= (bad_hygov_curve_model.verify() > 0); + + auto bad_hygov_trate = makeHygovData(); + bad_hygov_trate.parameters[Params::Trate] = static_cast(0.0); + PhasorDynamics::Governor::Hygov bad_hygov_trate_model(bad_hygov_trate); + success *= (bad_hygov_trate_model.verify() > 0); + + auto missing_hygov_trate = makeHygovData(); + missing_hygov_trate.parameters.erase(Params::Trate); + PhasorDynamics::Governor::Hygov missing_hygov_trate_model( + missing_hygov_trate); + success *= (missing_hygov_trate_model.verify() > 0); + + auto bad_hygov_hdam = makeHygovData(); + bad_hygov_hdam.parameters[Params::Hdam] = static_cast(0.0); + PhasorDynamics::Governor::Hygov bad_hygov_hdam_model( + bad_hygov_hdam); + success *= (bad_hygov_hdam_model.verify() > 0); + + auto bad_hygov_rtemp = makeHygovData(); + bad_hygov_rtemp.parameters[Params::Rtemp] = static_cast(0.0); + PhasorDynamics::Governor::Hygov bad_hygov_rtemp_model( + bad_hygov_rtemp); + success *= (bad_hygov_rtemp_model.verify() > 0); + + auto bad_hygov_limits = makeHygovData(); + bad_hygov_limits.parameters[Params::Velm] = static_cast(-0.1); + bad_hygov_limits.parameters[Params::db1] = static_cast(-0.1); + PhasorDynamics::Governor::Hygov bad_hygov_limits_model( + bad_hygov_limits); + success *= (bad_hygov_limits_model.verify() > 0); + + auto bad_hygov_gate_limits = makeHygovData(); + bad_hygov_gate_limits.parameters[Params::Gmin] = static_cast(1.1); + PhasorDynamics::Governor::Hygov bad_hygov_gate_limits_model( + bad_hygov_gate_limits); + success *= (bad_hygov_gate_limits_model.verify() > 0); + + return success.report(__func__); + } + + TestOutcome signals() + { + using Var = PhasorDynamics::Governor::HygovInternalVariables; + + TestStatus success = true; + + PhasorDynamics::SignalNode pmech_node; + PhasorDynamics::SignalNode omega_node; + ScalarT pmech_value{0.40}; + ScalarT omega_value{0.0}; + IdxT pmech_index = 4; + IdxT omega_index = 5; + pmech_node.set(&pmech_value, &pmech_index); + omega_node.set(&omega_value, &omega_index); + + PhasorDynamics::Governor::Hygov hygov(makeHygovData()); + auto& hygov_signals = hygov.getSignals(); + hygov_signals + .template assignSignalNode( + &pmech_node); + hygov_signals + .template attachSignalNode( + &omega_node); + + success *= (hygov.allocate() == 0); + hygov.y()[index(Var::PMECH)] = pmech_value; + success *= (hygov.verify() == 0); + success *= (hygov.initialize() == 0); + success *= (hygov.tagDifferentiable() == 0); + success *= (hygov.evaluateResidual() == 0); + + success *= isEqual(hygov.y()[index(Var::Q)], static_cast(0.50), kTol); + success *= isEqual(hygov.y()[index(Var::G)], static_cast(0.50), kTol); + success *= isEqual(hygov.y()[index(Var::C)], static_cast(0.50), kTol); + success *= (hygov.tag()[index(Var::G)] == true); + + for (size_t i = 0; i < hygov.getResidual().size(); ++i) + { + success *= isEqual(hygov.getResidual()[i], static_cast(0.0), kTol); + success *= isEqual(hygov.yp()[i], static_cast(0.0), kTol); + } + + return success.report(__func__); + } + + TestOutcome sourceDefault() + { + using Var = PhasorDynamics::Governor::HygovInternalVariables; + + TestStatus success = true; + + PhasorDynamics::SignalNode pmech_node; + ScalarT pmech_value{0.40}; + IdxT pmech_index = 4; + pmech_node.set(&pmech_value, &pmech_index); + + PhasorDynamics::Governor::Hygov hygov(makeHygovSourceDefaultData()); + hygov.getSignals() + .template assignSignalNode( + &pmech_node); + + success *= (hygov.allocate() == 0); + hygov.y()[index(Var::PMECH)] = pmech_value; + success *= (hygov.verify() == 0); + success *= (hygov.initialize() == 0); + success *= (hygov.tagDifferentiable() == 0); + success *= (hygov.evaluateResidual() == 0); + + success *= isEqual(hygov.y()[index(Var::Q)], static_cast(0.50), kTol); + success *= isEqual(hygov.y()[index(Var::PGV)], static_cast(0.50), kTol); + success *= isEqual(hygov.y()[index(Var::G)], static_cast(0.50), kTol); + success *= (hygov.tag()[index(Var::XN)] == true); + + for (size_t i = 0; i < hygov.getResidual().size(); ++i) + { + success *= isEqual(hygov.getResidual()[i], static_cast(0.0), kTol); + success *= isEqual(hygov.yp()[i], static_cast(0.0), kTol); + } + + return success.report(__func__); + } + + TestOutcome zeroTimeConstants() + { + using Params = PhasorDynamics::Governor::HygovParameters; + using Var = PhasorDynamics::Governor::HygovInternalVariables; + + TestStatus success = true; + + PhasorDynamics::SignalNode pmech_node; + ScalarT pmech_value{0.40}; + IdxT pmech_index = 4; + pmech_node.set(&pmech_value, &pmech_index); + + auto data = makeHygovData(); + data.parameters[Params::Tr] = static_cast(0.0); + data.parameters[Params::Tf] = static_cast(0.0); + data.parameters[Params::Tg] = static_cast(0.0); + data.parameters[Params::Tw] = static_cast(0.0); + data.parameters[Params::Tnp] = static_cast(0.0); + + PhasorDynamics::Governor::Hygov hygov(data); + hygov.getSignals() + .template assignSignalNode( + &pmech_node); + + success *= (hygov.allocate() == 0); + hygov.y()[index(Var::PMECH)] = pmech_value; + success *= (hygov.verify() == 0); + success *= (hygov.initialize() == 0); + success *= (hygov.tagDifferentiable() == 0); + success *= (hygov.evaluateResidual() == 0); + + success *= (hygov.tag()[index(Var::XN)] == true); + success *= (hygov.tag()[index(Var::XF)] == true); + success *= (hygov.tag()[index(Var::C)] == true); + success *= (hygov.tag()[index(Var::G)] == true); + success *= (hygov.tag()[index(Var::Q)] == true); + + for (size_t i = 0; i < hygov.getResidual().size(); ++i) + { + success *= isEqual(hygov.getResidual()[i], static_cast(0.0), kTol); + success *= isEqual(hygov.yp()[i], static_cast(0.0), kTol); + } + + return success.report(__func__); + } + + TestOutcome baseConversion() + { + using Var = PhasorDynamics::Governor::HygovInternalVariables; + using Params = PhasorDynamics::Governor::HygovParameters; + + TestStatus success = true; + + PhasorDynamics::SignalNode pmech_node; + ScalarT pmech_value{0.40}; + IdxT pmech_index = 6; + pmech_node.set(&pmech_value, &pmech_index); + + auto data = makeHygovData(); + data.parameters[Params::Trate] = static_cast(50.0); + + PhasorDynamics::Governor::Hygov hygov(data); + hygov.getSignals().template assignSignalNode(&pmech_node); + + success *= (hygov.allocate() == 0); + hygov.y()[index(Var::PMECH)] = pmech_value; + success *= (hygov.verify() == 0); + success *= (hygov.initialize() == 0); + success *= (hygov.evaluateResidual() == 0); + + success *= isEqual(hygov.y()[index(Var::Q)], static_cast(0.90), kTol); + success *= isEqual(hygov.y()[index(Var::G)], static_cast(0.90), kTol); + success *= isEqual(hygov.y()[index(Var::PMECH)], static_cast(0.40), kTol); + success *= isEqual(pmech_node.read(), static_cast(0.40), kTol); + + for (size_t i = 0; i < hygov.getResidual().size(); ++i) + { + success *= isEqual(hygov.getResidual()[i], static_cast(0.0), kTol); + success *= isEqual(hygov.yp()[i], static_cast(0.0), kTol); + } + + return success.report(__func__); + } + + TestOutcome jsonParseAndSystemAssembly() + { + TestStatus success = true; + + std::istringstream input(R"json( +{ + "header": { + "format_version": 0, + "format_revision": 1, + "case_name": "hydro governor", + "case_description": "HYGOV parser test", + "case_comments": "", + "freq_base": 60.0, + "va_base": 100000000.0 + }, + "buses": [ + { + "number": 1, + "class": "bus", + "name": "Bus 1", + "init": { "Vr": 1.0, "Vi": 0.0 }, + "v_base": 1.0 + } + ], + "signals": [ + { "signal_id": 10, "name": "Pmech" } + ], + "devices": [ + { + "class": "Genrou", + "ports": { "bus": 1, "pmech": 10 }, + "id": "GEN1", + "params": { + "p0": 0.3, "q0": 0.0, "H": 3.0, "D": 0.0, "Ra": 0.0, + "Tdop": 7.0, "Tdopp": 0.04, "Tqop": 0.75, "Tqopp": 0.05, + "Xd": 2.1, "Xdp": 0.2, "Xdpp": 0.18, "Xq": 0.5, "Xqp": 0.5, + "Xqpp": 0.18, "Xl": 0.15, "S10": 0.0, "S12": 0.0, + "mva": 100.0 + } + }, + { + "class": "Hygov", + "ports": { "pmech": 10 }, + "id": "HYG1", + "params": { + "Trate": 50.0, "Rperm": 0.05, "Rtemp": 0.4, "Tr": 5.0, + "Tf": 0.2, "Tg": 0.0, "Velm": 0.5, "Gmax": 1.0, "Gmin": 0.0, + "Tw": 1.0, "At": 1.0, "Dturb": 0.0, "Qnl": 0.1, + "Tn": 0.0, "Tnp": 1.0, "db1": 0.0, "db2": 0.0, "Hdam": 1.0, + "Gv0": 0.0, "Gv1": 0.2, "Gv2": 0.4, "Gv3": 0.6, "Gv4": 0.8, "Gv5": 1.0, + "Pgv0": 0.0, "Pgv1": 0.2, "Pgv2": 0.4, "Pgv3": 0.6, "Pgv4": 0.8, "Pgv5": 1.0 + } + } + ] +} +)json"); + + auto data = PhasorDynamics::parseSystemModelData(input); + success *= (data.hygov.size() == 1); + const auto trate_param = + data.hygov[0].parameters.at(PhasorDynamics::Governor::HygovParameters::Trate); + success *= (std::get(trate_param) == static_cast(50.0)); + PhasorDynamics::SystemModel system(data); + success *= (system.allocate() == 0); + success *= (system.initialize() == 0); + const auto hygov_size = + static_cast(PhasorDynamics::Governor::HygovInternalVariables::MAXIMUM); + const auto hygov_offset = static_cast(system.size() - hygov_size); + success *= isEqual( + system.y()[hygov_offset + index(PhasorDynamics::Governor::HygovInternalVariables::Q)], + static_cast(0.70), + kTol); + success *= isEqual( + system.y()[hygov_offset + index(PhasorDynamics::Governor::HygovInternalVariables::G)], + static_cast(0.70), + kTol); + success *= (system.evaluateResidual() == 0); + success *= (system.size() == 33); + + return success.report(__func__); + } + + private: + static size_t index(PhasorDynamics::Governor::HygovInternalVariables variable) + { + return static_cast(variable); + } + + auto makeHygovData() -> PhasorDynamics::Governor::HygovData + { + using Params = PhasorDynamics::Governor::HygovParameters; + using Mon = PhasorDynamics::Governor::HygovMonitorableVariables; + + PhasorDynamics::Governor::HygovData data; + data.device_class = "Hygov"; + data.disambiguation_string = "hygov_test"; + data.monitored_variables.insert(Mon::pmech); + data.monitored_variables.insert(Mon::gate); + + data.parameters[Params::Trate] = static_cast(100.0); + data.parameters[Params::Rperm] = static_cast(0.05); + data.parameters[Params::Rtemp] = static_cast(0.4); + data.parameters[Params::Tr] = static_cast(5.0); + data.parameters[Params::Tf] = static_cast(0.2); + data.parameters[Params::Tg] = static_cast(0.0); + data.parameters[Params::Velm] = static_cast(0.5); + data.parameters[Params::Gmax] = static_cast(1.0); + data.parameters[Params::Gmin] = static_cast(0.0); + data.parameters[Params::Tw] = static_cast(1.0); + data.parameters[Params::At] = static_cast(1.0); + data.parameters[Params::Dturb] = static_cast(0.0); + data.parameters[Params::Qnl] = static_cast(0.1); + data.parameters[Params::Tn] = static_cast(0.0); + data.parameters[Params::Tnp] = static_cast(1.0); + data.parameters[Params::db1] = static_cast(0.0); + data.parameters[Params::db2] = static_cast(0.0); + data.parameters[Params::Hdam] = static_cast(1.0); + data.parameters[Params::Gv0] = static_cast(0.0); + data.parameters[Params::Gv1] = static_cast(0.2); + data.parameters[Params::Gv2] = static_cast(0.4); + data.parameters[Params::Gv3] = static_cast(0.6); + data.parameters[Params::Gv4] = static_cast(0.8); + data.parameters[Params::Gv5] = static_cast(1.0); + data.parameters[Params::Pgv0] = static_cast(0.0); + data.parameters[Params::Pgv1] = static_cast(0.2); + data.parameters[Params::Pgv2] = static_cast(0.4); + data.parameters[Params::Pgv3] = static_cast(0.6); + data.parameters[Params::Pgv4] = static_cast(0.8); + data.parameters[Params::Pgv5] = static_cast(1.0); + + return data; + } + + auto makeHygovSourceDefaultData() -> PhasorDynamics::Governor::HygovData + { + using Params = PhasorDynamics::Governor::HygovParameters; + + auto data = makeHygovData(); + data.parameters[Params::Tn] = static_cast(0.0); + data.parameters[Params::Tnp] = static_cast(0.0); + data.parameters[Params::Gv0] = static_cast(0.0); + data.parameters[Params::Gv1] = static_cast(0.0); + data.parameters[Params::Gv2] = static_cast(0.0); + data.parameters[Params::Gv3] = static_cast(0.0); + data.parameters[Params::Gv4] = static_cast(0.0); + data.parameters[Params::Gv5] = static_cast(0.0); + data.parameters[Params::Pgv0] = static_cast(0.0); + data.parameters[Params::Pgv1] = static_cast(0.0); + data.parameters[Params::Pgv2] = static_cast(0.0); + data.parameters[Params::Pgv3] = static_cast(0.0); + data.parameters[Params::Pgv4] = static_cast(0.0); + data.parameters[Params::Pgv5] = static_cast(0.0); + + return data; + } + }; + } // namespace Testing +} // namespace GridKit diff --git a/tests/UnitTests/PhasorDynamics/runGovernorHygovTests.cpp b/tests/UnitTests/PhasorDynamics/runGovernorHygovTests.cpp new file mode 100644 index 000000000..8eb6e82f2 --- /dev/null +++ b/tests/UnitTests/PhasorDynamics/runGovernorHygovTests.cpp @@ -0,0 +1,17 @@ +#include "GovernorHygovTests.hpp" + +int main() +{ + GridKit::Testing::TestingResults result; + + GridKit::Testing::GovernorHygovTests test; + + result += test.constructionAndValidation(); + result += test.signals(); + result += test.sourceDefault(); + result += test.zeroTimeConstants(); + result += test.baseConversion(); + result += test.jsonParseAndSystemAssembly(); + + return result.summary(); +} From 78dd07c0a2615c8271d5c29e98b4c0188995de81 Mon Sep 17 00:00:00 2001 From: lukelowry Date: Fri, 26 Jun 2026 22:08:26 +0000 Subject: [PATCH 2/2] Apply pre-commit fixes --- .../Model/PhasorDynamics/SystemModelData.hpp | 30 +++++++++---------- 1 file changed, 15 insertions(+), 15 deletions(-) diff --git a/GridKit/Model/PhasorDynamics/SystemModelData.hpp b/GridKit/Model/PhasorDynamics/SystemModelData.hpp index 32c1a5c5c..7ce549029 100644 --- a/GridKit/Model/PhasorDynamics/SystemModelData.hpp +++ b/GridKit/Model/PhasorDynamics/SystemModelData.hpp @@ -108,21 +108,21 @@ namespace GridKit std::vector stabilizer; ///< Stabilizers within the model std::vector constant_source; ///< Constant signal sources within the model std::vector signal; ///< Signal nodes - std::vector bus; ///< Buses within the model - std::vector adapter; ///< bus-to-signal adapters within the model - std::vector branch; ///< Branches within the model - std::vector bus_fault; ///< Bus faults within the model - std::vector genrou; ///< GENROU instances within the model - std::vector gensal; ///< GENSAL instances within the model - std::vector genclassical; ///< Classical generator instances within the model - std::vector loadz; ///< LoadZ instances within the model - std::vector loadzip; ///< LoadZIP instances within the model - std::vector gov; ///< Governors within the model - std::vector hygov; ///< HYGOV governors within the model - std::vector exciter; ///< Exciters within the model - std::vector sexspti; ///< SEXS-PTI exciters within the model - std::vector stabilizer; ///< Stabilizers within the model - std::vector signal; ///< Signal nodes + std::vector bus; ///< Buses within the model + std::vector adapter; ///< bus-to-signal adapters within the model + std::vector branch; ///< Branches within the model + std::vector bus_fault; ///< Bus faults within the model + std::vector genrou; ///< GENROU instances within the model + std::vector gensal; ///< GENSAL instances within the model + std::vector genclassical; ///< Classical generator instances within the model + std::vector loadz; ///< LoadZ instances within the model + std::vector loadzip; ///< LoadZIP instances within the model + std::vector gov; ///< Governors within the model + std::vector hygov; ///< HYGOV governors within the model + std::vector exciter; ///< Exciters within the model + std::vector sexspti; ///< SEXS-PTI exciters within the model + std::vector stabilizer; ///< Stabilizers within the model + std::vector signal; ///< Signal nodes /// Monitor sink specs std::vector monitor_sink;