From 04002ed99c861d2879888a9445d56f0092147f62 Mon Sep 17 00:00:00 2001 From: "martin.moraga" Date: Wed, 15 May 2024 10:38:58 +0200 Subject: [PATCH] WIP Signed-off-by: martin.moraga --- .../dpsim-models/Signal/VSIControlType1.h | 276 ++++---- .../Base/Base_VSIVoltageSourceInverterDQ.cpp | 226 +++---- .../src/DP/DP_Ph1_VSIVoltageControlDQ.cpp | 304 +++++---- .../src/DP/DP_Ph1_VSIVoltageControlVBR.cpp | 118 ++-- .../src/EMT/EMT_Ph3_VSIVoltageControlDQ.cpp | 19 +- .../src/EMT/EMT_Ph3_VSIVoltageControlVBR.cpp | 631 ++++++++++-------- dpsim-models/src/SP/SP_Ph1_Capacitor.cpp | 2 +- .../src/SP/SP_Ph1_VSIVoltageControlVBR.cpp | 423 ++++++------ dpsim-models/src/Signal/VSIControlType1.cpp | 460 +++++++------ dpsim/examples/cxx/CMakeLists.txt | 1 + .../GridFormingInverter_ControllerType1.ipynb | 37 +- 11 files changed, 1336 insertions(+), 1161 deletions(-) diff --git a/dpsim-models/include/dpsim-models/Signal/VSIControlType1.h b/dpsim-models/include/dpsim-models/Signal/VSIControlType1.h index 46d23bdafd..c38b9f0afe 100644 --- a/dpsim-models/include/dpsim-models/Signal/VSIControlType1.h +++ b/dpsim-models/include/dpsim-models/Signal/VSIControlType1.h @@ -1,144 +1,148 @@ #pragma once -#include #include #include +#include namespace CPS { namespace Signal { - class VSIControlType1Parameters : - public Base::VSIControlParameters, - public SharedFactory { - - public: - /// - Real Kpv; - /// - Real Kiv; - /// - Real Kpc; - /// - Real Kic; - /// - Real VdRef; - /// - Real VqRef; - /// Equivalent time constant of the PT1 block - Real tau; - }; - - /// DOC: Controller for a grid-forming power converter, - /// which is implemented by using two cascaded - /// PI controllers working on the dq reference frame. - /// The output of the inverter is used as the reference voltage - /// for the voltage source of the inverter - /// Optional: if mModelAsCurrentSource == true, the PI controller - /// modelling the current loop is replaced by one equivalent PT1 block - /// In this case the output of the inverter is the equivalent current flowing - /// into the RLC-Filter, e.g. the reference for the current source - /// of the inverter - /// *** This controller does not consider the feedforward terms - /// Ref.: Joan Rocabert, Control of Power Converters in AC Microgrids - /// & BA of Matthias Mees - class VSIControlType1 : - public SimSignalComp, - public Base::VSIControlDQ, - public SharedFactory { - - public: - /// @brief VBR auxiliar variables - Real mA_VBR; - Real mB_VBR; - Real mB2_VBR; - Real mC_VBR; - Real mD_VBR; - Real mE_VBR; - - private: - /// Controller Parameters - std::shared_ptr mParameters; - - /// Controller variables - /// state variable of the outer loop (d-component) - const Attribute::Ptr mPhi_d; - /// state variable of the outer loop (q-component) - const Attribute::Ptr mPhi_q; - /// state variable of the inner loop (d-component) - const Attribute::Ptr mGamma_d; - /// state variable of the inner loop (q-component) - const Attribute::Ptr mGamma_q; - - // State space matrices - /// matrix A of state space model - Matrix mA = Matrix::Zero(4, 4); - /// matrix B of state space model - Matrix mB = Matrix::Zero(4, 6); - /// matrix C of state space model - Matrix mC = Matrix::Zero(2, 4); - /// matrix D of state space model - Matrix mD = Matrix::Zero(2, 6); - - /// Trapedzoidal based state space matrix A - Matrix mATrapezoidal; - /// Trapedzoidal based state space matrix B - Matrix mBTrapezoidal; - /// Trapedzoidal based state space matrix Matrix::Zero(2,1) - Matrix mCTrapezoidal; - - // input, state and output vectors - // matrixes of the seperate control loops previous and current - /// Previous Input - const Attribute::Ptr mInputPrev; - /// Current Input - const Attribute::Ptr mInputCurr; - /// Previous State - const Attribute::Ptr mStatePrev; - /// Current State - const Attribute::Ptr mStateCurr; - /// Output - const Attribute::Ptr mOutput; - - /// ### VBR Auxiliar Variables ### - Complex mVcap_dq; - Complex mIfilter_dq; - Real mIref_d; - Real mIref_q; - - public: - /// - explicit VSIControlType1(const String & name, CPS::Logger::Level logLevel) : - SimSignalComp(name, name, logLevel), - mPhi_d(mAttributes->create("Phi_d", 0)), - mPhi_q(mAttributes->create("Phi_q", 0)), - mGamma_d(mAttributes->create("Gamma_d", 0)), - mGamma_q(mAttributes->create("Gamma_q", 0)), - mInputPrev(mAttributes->create("input_prev", Matrix::Zero(6,1))), - mInputCurr(mAttributes->create("input_curr", Matrix::Zero(6,1))), - mStatePrev(mAttributes->create("state_prev", Matrix::Zero(4,1))), - mStateCurr(mAttributes->create("state_curr", Matrix::Zero(4,1))), - mOutput(mAttributes->create("output", Matrix::Zero(2,1))) { } - - /// Constructor with log level - VSIControlType1(const String & name); - - /// Sets Parameters of the vsi controller - void setParameters(std::shared_ptr parameters) final; - - /// - void calculateVBRconstants() final; - - /// Initialises the initial state of the vsi controller - void initialize(const Complex& Vsref_dq, const Complex& Vcap_dq, - const Complex& Ifilter_dq, Real time_step, Bool modelAsCurrentSource) final; - - /// Performs a step to update all state variables and the output - Complex step(const Complex& Vcap_dq, const Complex& Ifilter_dq) final; - - /// - Complex stepVBR(const Complex& Vcap_dq, const Complex& Ifilter_dq) final; - - }; - -} -} \ No newline at end of file +class VSIControlType1Parameters + : public Base::VSIControlParameters, + public SharedFactory { + +public: + /// + Real Kpv; + /// + Real Kiv; + /// + Real Kpc; + /// + Real Kic; + /// + Real VdRef; + /// + Real VqRef; + /// Equivalent time constant of the PT1 block + Real tau; +}; + +/// DOC: Controller for a grid-forming power converter, +/// which is implemented by using two cascaded +/// PI controllers working on the dq reference frame. +/// The output of the inverter is used as the reference voltage +/// for the voltage source of the inverter +/// Optional: if mModelAsCurrentSource == true, the PI controller +/// modelling the current loop is replaced by one equivalent PT1 block +/// In this case the output of the inverter is the equivalent current flowing +/// into the RLC-Filter, e.g. the reference for the current source +/// of the inverter +/// *** This controller does not consider the feedforward terms +/// Ref.: Joan Rocabert, Control of Power Converters in AC Microgrids +/// & BA of Matthias Mees +class VSIControlType1 : public SimSignalComp, + public Base::VSIControlDQ, + public SharedFactory { + +public: + /// @brief VBR auxiliar variables + Real mA_VBR; + Real mB_VBR; + Real mB2_VBR; + Real mC_VBR; + Real mD_VBR; + Real mE_VBR; + +private: + /// Controller Parameters + std::shared_ptr mParameters; + + /// Controller variables + /// state variable of the outer loop (d-component) + const Attribute::Ptr mPhi_d; + /// state variable of the outer loop (q-component) + const Attribute::Ptr mPhi_q; + /// state variable of the inner loop (d-component) + const Attribute::Ptr mGamma_d; + /// state variable of the inner loop (q-component) + const Attribute::Ptr mGamma_q; + + // State space matrices + /// matrix A of state space model + Matrix mA = Matrix::Zero(4, 4); + /// matrix B of state space model + Matrix mB = Matrix::Zero(4, 6); + /// matrix C of state space model + Matrix mC = Matrix::Zero(2, 4); + /// matrix D of state space model + Matrix mD = Matrix::Zero(2, 6); + + /// Trapedzoidal based state space matrix A + Matrix mATrapezoidal; + /// Trapedzoidal based state space matrix B + Matrix mBTrapezoidal; + /// Trapedzoidal based state space matrix Matrix::Zero(2,1) + Matrix mCTrapezoidal; + + // input, state and output vectors + // matrixes of the seperate control loops previous and current + /// Previous Input + const Attribute::Ptr mInputPrev; + /// Current Input + const Attribute::Ptr mInputCurr; + /// Previous State + const Attribute::Ptr mStatePrev; + /// Current State + const Attribute::Ptr mStateCurr; + /// Output + const Attribute::Ptr mOutput; + + /// ### VBR Auxiliar Variables ### + Complex mVcap_dq; + Complex mIfilter_dq; + Real mIref_d; + Real mIref_q; + +public: + /// + explicit VSIControlType1(const String &name, CPS::Logger::Level logLevel) + : SimSignalComp(name, name, logLevel), + mPhi_d(mAttributes->create("Phi_d", 0)), + mPhi_q(mAttributes->create("Phi_q", 0)), + mGamma_d(mAttributes->create("Gamma_d", 0)), + mGamma_q(mAttributes->create("Gamma_q", 0)), + mInputPrev( + mAttributes->create("input_prev", Matrix::Zero(6, 1))), + mInputCurr( + mAttributes->create("input_curr", Matrix::Zero(6, 1))), + mStatePrev( + mAttributes->create("state_prev", Matrix::Zero(4, 1))), + mStateCurr( + mAttributes->create("state_curr", Matrix::Zero(4, 1))), + mOutput(mAttributes->create("output", Matrix::Zero(2, 1))) {} + + /// Constructor with log level + VSIControlType1(const String &name); + + /// Sets Parameters of the vsi controller + void + setParameters(std::shared_ptr parameters) final; + + /// + void calculateVBRconstants() final; + + /// Initialises the initial state of the vsi controller + void initialize(const Complex &Vsref_dq, const Complex &Vcap_dq, + const Complex &Ifilter_dq, Real time_step, + Bool modelAsCurrentSource) final; + + /// Performs a step to update all state variables and the output + Complex step(const Complex &Vcap_dq, const Complex &Ifilter_dq) final; + + /// + Complex stepVBR(const Complex &Vcap_dq, const Complex &Ifilter_dq) final; +}; + +} // namespace Signal +} // namespace CPS \ No newline at end of file diff --git a/dpsim-models/src/Base/Base_VSIVoltageSourceInverterDQ.cpp b/dpsim-models/src/Base/Base_VSIVoltageSourceInverterDQ.cpp index 0a3dc7ff68..2b93b9f484 100644 --- a/dpsim-models/src/Base/Base_VSIVoltageSourceInverterDQ.cpp +++ b/dpsim-models/src/Base/Base_VSIVoltageSourceInverterDQ.cpp @@ -11,137 +11,137 @@ using namespace CPS; template -void Base::VSIVoltageSourceInverterDQ::setParameters( - Real sysOmega, Real VdRef, Real VqRef) -{ - - mOmegaNom = sysOmega; - mVdRef = VdRef; - mVqRef = VqRef; - - SPDLOG_LOGGER_INFO(mLogger, - "\nGeneral Parameters:" - "\n\tNominal Omega = {} [1/s]" - "\n\tVdRef = {} [V] " - "\n\tVqRef = {} [V]", - mOmegaNom, mVdRef, mVqRef); +void Base::VSIVoltageSourceInverterDQ::setParameters(Real sysOmega, + Real VdRef, + Real VqRef) { + + mOmegaNom = sysOmega; + mVdRef = VdRef; + mVqRef = VqRef; + + SPDLOG_LOGGER_INFO(mLogger, + "\nGeneral Parameters:" + "\n\tNominal Omega = {} [1/s]" + "\n\tVdRef = {} [V] " + "\n\tVqRef = {} [V]", + mOmegaNom, mVdRef, mVqRef); } template -void Base::VSIVoltageSourceInverterDQ::setFilterParameters( - Real Lf, Real Cf, Real Rf, Real Rc) -{ - - mLf = Lf; - mCf = Cf; - mRf = Rf; - mRc = Rc; - - createSubComponents(); - - SPDLOG_LOGGER_INFO(mLogger, - "\nFilter Parameters:" - "\n\tInductance Lf = {} [H]" - "\n\tCapacitance Cf = {} [F]" - "\n\tResistance Rf = {} [H]" - "\n\tResistance Rc = {} [F]", - mLf, mCf, mRf, mRc); - mLogger->flush(); +void Base::VSIVoltageSourceInverterDQ::setFilterParameters(Real Lf, + Real Cf, + Real Rf, + Real Rc) { + + mLf = Lf; + mCf = Cf; + mRf = Rf; + mRc = Rc; + + createSubComponents(); + + SPDLOG_LOGGER_INFO(mLogger, + "\nFilter Parameters:" + "\n\tInductance Lf = {} [H]" + "\n\tCapacitance Cf = {} [F]" + "\n\tResistance Rf = {} [H]" + "\n\tResistance Rc = {} [F]", + mLf, mCf, mRf, mRc); + mLogger->flush(); } template -int Base::VSIVoltageSourceInverterDQ::determineNumberOfVirtualNodes() -{ - // first virtual node is the second node of the rl element - int numberOfVirtualNodes = 1; +int Base::VSIVoltageSourceInverterDQ::determineNumberOfVirtualNodes() { + // first virtual node is the second node of the rl element + int numberOfVirtualNodes = 1; - if (mWithInterfaceResistor) - numberOfVirtualNodes += 1; + if (mWithInterfaceResistor) + numberOfVirtualNodes += 1; - return numberOfVirtualNodes; + return numberOfVirtualNodes; } template -void Base::VSIVoltageSourceInverterDQ::addVSIController(std::shared_ptr VSIController) -{ - mVSIController = VSIController; - mWithControl = true; +void Base::VSIVoltageSourceInverterDQ::addVSIController( + std::shared_ptr VSIController) { + mVSIController = VSIController; + mWithControl = true; } template -void Base::VSIVoltageSourceInverterDQ::addDroopController(std::shared_ptr DroopController) -{ - mDroopController = DroopController; - mWithDroop = true; +void Base::VSIVoltageSourceInverterDQ::addDroopController( + std::shared_ptr DroopController) { + mDroopController = DroopController; + mWithDroop = true; } template void Base::VSIVoltageSourceInverterDQ::initializeFilterVariables( - const Complex &interfaceVoltage, const Complex &interfaceCurrent, - typename SimNode::List virtualNodesList) -{ - - // derive initialization quantities of filter - /// initial filter capacitor voltage - Complex vcInit; - if (mWithInterfaceResistor) - vcInit = interfaceVoltage + interfaceCurrent * mRc; - else - vcInit = interfaceVoltage; - - /// initial filter capacitor current - Complex icfInit = vcInit * Complex(0., mOmegaNom * mCf); - - /// initial voltage/current equivalent source - Complex filterCurrentInit = interfaceCurrent + icfInit; - Complex sourceInitialValue; - if (mModelAsCurrentSource) - sourceInitialValue = filterCurrentInit; - else - sourceInitialValue = vcInit + filterCurrentInit * Complex(mRf, mOmegaNom * mLf); - - // initialize angles - **mThetaSys = 0; - **mThetaInv = std::arg(vcInit); - - // Initialie filter variables in dq domain - **mVcap_dq = Math::rotatingFrame2to1(vcInit, **mThetaInv, **mThetaSys); - **mIfilter_dq = Math::rotatingFrame2to1(filterCurrentInit, **mThetaInv, **mThetaSys); - **mSourceValue_dq = Math::rotatingFrame2to1(sourceInitialValue, **mThetaInv, **mThetaSys); - - String inverter_type = mModelAsCurrentSource ? "current source" : "voltage source"; - String unit = mModelAsCurrentSource ? "[A]" : "[V]"; - SPDLOG_LOGGER_INFO(mLogger, - "\nInverter will be modelled as {}" - "\nInitialize Filter Variables:" - "\n\tInitial capacitor voltage: {}[V]" - "\n\tInitial capacitor voltage d-axis: {}[V]" - "\n\tInitial capacitor voltage q-axis: {}[V]" - "\n\tInitial filter current: {}[A]" - "\n\tInitial filter d-axis: {}[A]" - "\n\tInitial filter q-axis: {}[A]" - "\n\tInitial equivalent source: {}{}" - "\n\tInverter equivalent source d-axis value: {}{}" - "\n\tInverter equivalent source q-axis value: {}{}", - inverter_type, - Logger::phasorToString(vcInit), - (**mVcap_dq).real(), (**mVcap_dq).imag(), - Logger::phasorToString(filterCurrentInit), - (**mIfilter_dq).real(), (**mIfilter_dq).imag(), - sourceInitialValue, unit, - (**mSourceValue_dq).real(), unit, - (**mSourceValue_dq).imag(), unit); - mLogger->flush(); - - // TODO: MOVE - // initialize voltage of virtual nodes - virtualNodesList[0]->setInitialVoltage(vcInit + filterCurrentInit * Complex(mRf, mOmegaNom * mLf)); - if (mWithInterfaceResistor) - { - // filter capacitor is connected to mVirtualNodes[1], the second - // node of the interface resistor is mTerminals[0] - virtualNodesList[1]->setInitialVoltage(vcInit); - } + const Complex &interfaceVoltage, const Complex &interfaceCurrent, + typename SimNode::List virtualNodesList) { + + // derive initialization quantities of filter + /// initial filter capacitor voltage + Complex vcInit; + if (mWithInterfaceResistor) + vcInit = interfaceVoltage + interfaceCurrent * mRc; + else + vcInit = interfaceVoltage; + + /// initial filter capacitor current + Complex icfInit = vcInit * Complex(0., mOmegaNom * mCf); + + /// initial voltage/current equivalent source + Complex filterCurrentInit = interfaceCurrent + icfInit; + Complex sourceInitialValue; + if (mModelAsCurrentSource) + sourceInitialValue = filterCurrentInit; + else + sourceInitialValue = + vcInit + filterCurrentInit * Complex(mRf, mOmegaNom * mLf); + + // initialize angles + **mThetaSys = 0; + **mThetaInv = std::arg(vcInit); + + // Initialie filter variables in dq domain + **mVcap_dq = Math::rotatingFrame2to1(vcInit, **mThetaInv, **mThetaSys); + **mIfilter_dq = + Math::rotatingFrame2to1(filterCurrentInit, **mThetaInv, **mThetaSys); + **mSourceValue_dq = + Math::rotatingFrame2to1(sourceInitialValue, **mThetaInv, **mThetaSys); + + String inverter_type = + mModelAsCurrentSource ? "current source" : "voltage source"; + String unit = mModelAsCurrentSource ? "[A]" : "[V]"; + SPDLOG_LOGGER_INFO( + mLogger, + "\nInverter will be modelled as {}" + "\nInitialize Filter Variables:" + "\n\tInitial capacitor voltage: {}[V]" + "\n\tInitial capacitor voltage d-axis: {}[V]" + "\n\tInitial capacitor voltage q-axis: {}[V]" + "\n\tInitial filter current: {}[A]" + "\n\tInitial filter d-axis: {}[A]" + "\n\tInitial filter q-axis: {}[A]" + "\n\tInitial equivalent source: {}{}" + "\n\tInverter equivalent source d-axis value: {}{}" + "\n\tInverter equivalent source q-axis value: {}{}", + inverter_type, Logger::phasorToString(vcInit), (**mVcap_dq).real(), + (**mVcap_dq).imag(), Logger::phasorToString(filterCurrentInit), + (**mIfilter_dq).real(), (**mIfilter_dq).imag(), sourceInitialValue, unit, + (**mSourceValue_dq).real(), unit, (**mSourceValue_dq).imag(), unit); + mLogger->flush(); + + // TODO: MOVE + // initialize voltage of virtual nodes + virtualNodesList[0]->setInitialVoltage( + vcInit + filterCurrentInit * Complex(mRf, mOmegaNom * mLf)); + if (mWithInterfaceResistor) { + // filter capacitor is connected to mVirtualNodes[1], the second + // node of the interface resistor is mTerminals[0] + virtualNodesList[1]->setInitialVoltage(vcInit); + } } // Declare specializations to move definitions to .cpp diff --git a/dpsim-models/src/DP/DP_Ph1_VSIVoltageControlDQ.cpp b/dpsim-models/src/DP/DP_Ph1_VSIVoltageControlDQ.cpp index d23802e3ba..d99c877eaf 100644 --- a/dpsim-models/src/DP/DP_Ph1_VSIVoltageControlDQ.cpp +++ b/dpsim-models/src/DP/DP_Ph1_VSIVoltageControlDQ.cpp @@ -10,164 +10,200 @@ using namespace CPS; -DP::Ph1::VSIVoltageControlDQ::VSIVoltageControlDQ(String uid, String name, Logger::Level logLevel, - Bool modelAsCurrentSource, Bool withInterfaceResistor) : - CompositePowerComp(uid, name, true, true, logLevel), - VSIVoltageSourceInverterDQ(this->mSLog, mAttributes, modelAsCurrentSource, withInterfaceResistor) { - - setTerminalNumber(1); - setVirtualNodeNumber(this->determineNumberOfVirtualNodes()); - - **mIntfVoltage = MatrixComp::Zero(1, 1); - **mIntfCurrent = MatrixComp::Zero(1, 1); +DP::Ph1::VSIVoltageControlDQ::VSIVoltageControlDQ(String uid, String name, + Logger::Level logLevel, + Bool modelAsCurrentSource, + Bool withInterfaceResistor) + : CompositePowerComp(uid, name, true, true, logLevel), + VSIVoltageSourceInverterDQ(this->mSLog, mAttributes, + modelAsCurrentSource, + withInterfaceResistor) { + + setTerminalNumber(1); + setVirtualNodeNumber(this->determineNumberOfVirtualNodes()); + + **mIntfVoltage = MatrixComp::Zero(1, 1); + **mIntfCurrent = MatrixComp::Zero(1, 1); } -void DP::Ph1::VSIVoltageControlDQ::createSubComponents() { - // voltage source - if (!mModelAsCurrentSource) { - mSubCtrledVoltageSource = DP::Ph1::VoltageSource::make(**mName + "_src", mLogLevel); - addMNASubComponent(mSubCtrledVoltageSource, MNA_SUBCOMP_TASK_ORDER::NO_TASK, MNA_SUBCOMP_TASK_ORDER::TASK_BEFORE_PARENT, true); - } - - // RL Element as part of the LC filter - mSubFilterRL = DP::Ph1::ResIndSeries::make(**mName + "_FilterRL", mLogLevel); - addMNASubComponent(mSubFilterRL, MNA_SUBCOMP_TASK_ORDER::TASK_BEFORE_PARENT, MNA_SUBCOMP_TASK_ORDER::TASK_BEFORE_PARENT, true); - mSubFilterRL->setParameters(mRf, mLf); - - // Capacitor as part of the LC filter - mSubCapacitorF = DP::Ph1::Capacitor::make(**mName + "_CapF", mLogLevel); - addMNASubComponent(mSubCapacitorF, MNA_SUBCOMP_TASK_ORDER::TASK_BEFORE_PARENT, MNA_SUBCOMP_TASK_ORDER::TASK_BEFORE_PARENT, true); - mSubCapacitorF->setParameters(mCf); - - // optinal: interface resistor - if (mWithInterfaceResistor) { - mSubResistorC = DP::Ph1::Resistor::make(**mName + "_ResC", mLogLevel); - mSubResistorC->setParameters(mRc); - addMNASubComponent(mSubResistorC, MNA_SUBCOMP_TASK_ORDER::TASK_BEFORE_PARENT, MNA_SUBCOMP_TASK_ORDER::TASK_BEFORE_PARENT, false); - } +void DP::Ph1::VSIVoltageControlDQ::createSubComponents() { + // voltage source + if (!mModelAsCurrentSource) { + mSubCtrledVoltageSource = + DP::Ph1::VoltageSource::make(**mName + "_src", mLogLevel); + addMNASubComponent(mSubCtrledVoltageSource, MNA_SUBCOMP_TASK_ORDER::NO_TASK, + MNA_SUBCOMP_TASK_ORDER::TASK_BEFORE_PARENT, true); + } + + // RL Element as part of the LC filter + mSubFilterRL = DP::Ph1::ResIndSeries::make(**mName + "_FilterRL", mLogLevel); + addMNASubComponent(mSubFilterRL, MNA_SUBCOMP_TASK_ORDER::TASK_BEFORE_PARENT, + MNA_SUBCOMP_TASK_ORDER::TASK_BEFORE_PARENT, true); + mSubFilterRL->setParameters(mRf, mLf); + + // Capacitor as part of the LC filter + mSubCapacitorF = DP::Ph1::Capacitor::make(**mName + "_CapF", mLogLevel); + addMNASubComponent(mSubCapacitorF, MNA_SUBCOMP_TASK_ORDER::TASK_BEFORE_PARENT, + MNA_SUBCOMP_TASK_ORDER::TASK_BEFORE_PARENT, true); + mSubCapacitorF->setParameters(mCf); + + // optinal: interface resistor + if (mWithInterfaceResistor) { + mSubResistorC = DP::Ph1::Resistor::make(**mName + "_ResC", mLogLevel); + mSubResistorC->setParameters(mRc); + addMNASubComponent(mSubResistorC, + MNA_SUBCOMP_TASK_ORDER::TASK_BEFORE_PARENT, + MNA_SUBCOMP_TASK_ORDER::TASK_BEFORE_PARENT, false); + } } void DP::Ph1::VSIVoltageControlDQ::connectSubComponents() { - // TODO: COULD WE MOVE THIS FUNCTION TO THE BASE CLASS? - if (!mModelAsCurrentSource) - mSubCtrledVoltageSource->connect({ SimNode::GND, mVirtualNodes[0] }); - if (mWithInterfaceResistor) { - // with interface resistor - mSubFilterRL->connect({ mVirtualNodes[1], mVirtualNodes[0] }); - mSubCapacitorF->connect({ SimNode::GND, mVirtualNodes[1] }); - mSubResistorC->connect({ mTerminals[0]->node(), mVirtualNodes[1]}); - } else { - // without interface resistor - mSubFilterRL->connect({ mTerminals[0]->node(), mVirtualNodes[0] }); - mSubCapacitorF->connect({ SimNode::GND, mTerminals[0]->node() }); - } + // TODO: COULD WE MOVE THIS FUNCTION TO THE BASE CLASS? + if (!mModelAsCurrentSource) + mSubCtrledVoltageSource->connect({SimNode::GND, mVirtualNodes[0]}); + if (mWithInterfaceResistor) { + // with interface resistor + mSubFilterRL->connect({mVirtualNodes[1], mVirtualNodes[0]}); + mSubCapacitorF->connect({SimNode::GND, mVirtualNodes[1]}); + mSubResistorC->connect({mTerminals[0]->node(), mVirtualNodes[1]}); + } else { + // without interface resistor + mSubFilterRL->connect({mTerminals[0]->node(), mVirtualNodes[0]}); + mSubCapacitorF->connect({SimNode::GND, mTerminals[0]->node()}); + } } -void DP::Ph1::VSIVoltageControlDQ::initializeFromNodesAndTerminals(Real frequency) { - // terminal powers in consumer system -> convert to generator system - **mPower = -terminal(0)->singlePower(); - - // set initial interface quantities --> Current flowing into the inverter is positive - (**mIntfVoltage)(0, 0) = initialSingleVoltage(0); - (**mIntfCurrent)(0, 0) = std::conj(**mPower / (**mIntfVoltage)(0,0)); - - // initialize filter variables and set initial voltage of virtual nodes - initializeFilterVariables((**mIntfVoltage)(0, 0), (**mIntfCurrent)(0, 0), mVirtualNodes); - - // calculate initial source value - (**mSourceValue)(0,0) = Math::rotatingFrame2to1(**mSourceValue_dq, **mThetaSys, **mThetaInv); - - // Create & Initialize electrical subcomponents - this->connectSubComponents(); - for (auto subcomp: mSubComponents) { - subcomp->initialize(mFrequencies); - subcomp->initializeFromNodesAndTerminals(frequency); - } - - // TODO: droop - **mOmega = mOmegaNom; - - SPDLOG_LOGGER_INFO(mSLog, - "\n--- Initialization from powerflow ---" - "\nTerminal 0 connected to {} = sim node {}" - "\nInverter terminal voltage: {}[V]" - "\nInverter output current: {}[A]", - mTerminals[0]->node()->name(), mTerminals[0]->node()->matrixNodeIndex(), - Logger::phasorToString((**mIntfVoltage)(0, 0)), - Logger::phasorToString((**mIntfCurrent)(0, 0))); - mSLog->flush(); +void DP::Ph1::VSIVoltageControlDQ::initializeFromNodesAndTerminals( + Real frequency) { + // terminal powers in consumer system -> convert to generator system + **mPower = -terminal(0)->singlePower(); + + // set initial interface quantities --> Current flowing into the inverter is positive + (**mIntfVoltage)(0, 0) = initialSingleVoltage(0); + (**mIntfCurrent)(0, 0) = std::conj(**mPower / (**mIntfVoltage)(0, 0)); + + // initialize filter variables and set initial voltage of virtual nodes + initializeFilterVariables((**mIntfVoltage)(0, 0), (**mIntfCurrent)(0, 0), + mVirtualNodes); + + // calculate initial source value + (**mSourceValue)(0, 0) = + Math::rotatingFrame2to1(**mSourceValue_dq, **mThetaSys, **mThetaInv); + + // Create & Initialize electrical subcomponents + this->connectSubComponents(); + for (auto subcomp : mSubComponents) { + subcomp->initialize(mFrequencies); + subcomp->initializeFromNodesAndTerminals(frequency); + } + + // TODO: droop + **mOmega = mOmegaNom; + + SPDLOG_LOGGER_INFO(mSLog, + "\n--- Initialization from powerflow ---" + "\nTerminal 0 connected to {} = sim node {}" + "\nInverter terminal voltage: {}[V]" + "\nInverter output current: {}[A]", + mTerminals[0]->node()->name(), + mTerminals[0]->node()->matrixNodeIndex(), + Logger::phasorToString((**mIntfVoltage)(0, 0)), + Logger::phasorToString((**mIntfCurrent)(0, 0))); + mSLog->flush(); } -void DP::Ph1::VSIVoltageControlDQ::mnaParentInitialize(Real omega, Real timeStep, Attribute::Ptr leftVector) { - this->updateMatrixNodeIndices(); - mTimeStep = timeStep; - if (mWithControl) - mVSIController->initialize(**mSourceValue_dq, **mVcap_dq, **mIfilter_dq, mTimeStep, mModelAsCurrentSource); +void DP::Ph1::VSIVoltageControlDQ::mnaParentInitialize( + Real omega, Real timeStep, Attribute::Ptr leftVector) { + this->updateMatrixNodeIndices(); + mTimeStep = timeStep; + if (mWithControl) + mVSIController->initialize(**mSourceValue_dq, **mVcap_dq, **mIfilter_dq, + mTimeStep, mModelAsCurrentSource); } -void DP::Ph1::VSIVoltageControlDQ::mnaParentAddPreStepDependencies(AttributeBase::List &prevStepDependencies, AttributeBase::List &attributeDependencies, AttributeBase::List &modifiedAttributes) { - modifiedAttributes.push_back(mRightVector); - prevStepDependencies.push_back(mIntfVoltage); +void DP::Ph1::VSIVoltageControlDQ::mnaParentAddPreStepDependencies( + AttributeBase::List &prevStepDependencies, + AttributeBase::List &attributeDependencies, + AttributeBase::List &modifiedAttributes) { + modifiedAttributes.push_back(mRightVector); + prevStepDependencies.push_back(mIntfVoltage); } -void DP::Ph1::VSIVoltageControlDQ::mnaParentPreStep(Real time, Int timeStepCount) { - // get measurements - **mVcap_dq = Math::rotatingFrame2to1((**mSubCapacitorF->mIntfVoltage)(0,0), **mThetaInv, **mThetaSys); - **mIfilter_dq = Math::rotatingFrame2to1((**mSubFilterRL->mIntfCurrent)(0, 0), **mThetaInv, **mThetaSys); - - // TODO: droop - //if (mWithDroop) - // mDroop->signalStep(time, timeStepCount); - - // VCO Step - **mThetaInv = **mThetaInv + mTimeStep * **mOmega; - - // Update nominal system angle - **mThetaSys = **mThetaSys + mTimeStep * mOmegaNom; - - // - if (mWithControl) - **mSourceValue_dq = mVSIController->step(**mVcap_dq, **mIfilter_dq); - - // Transformation interface backward - (**mSourceValue)(0,0) = Math::rotatingFrame2to1(**mSourceValue_dq, **mThetaSys, **mThetaInv); - - // set reference voltage of voltage source - if (!mModelAsCurrentSource) { - // pre-step of voltage source - **mSubCtrledVoltageSource->mVoltageRef = (**mSourceValue)(0,0); - std::dynamic_pointer_cast(mSubCtrledVoltageSource)->mnaPreStep(time, timeStepCount); - } - - // stamp right side vector - mnaApplyRightSideVectorStamp(**mRightVector); +void DP::Ph1::VSIVoltageControlDQ::mnaParentPreStep(Real time, + Int timeStepCount) { + // get measurements + **mVcap_dq = Math::rotatingFrame2to1((**mSubCapacitorF->mIntfVoltage)(0, 0), + **mThetaInv, **mThetaSys); + **mIfilter_dq = Math::rotatingFrame2to1((**mSubFilterRL->mIntfCurrent)(0, 0), + **mThetaInv, **mThetaSys); + + // TODO: droop + //if (mWithDroop) + // mDroop->signalStep(time, timeStepCount); + + // VCO Step + **mThetaInv = **mThetaInv + mTimeStep * **mOmega; + + // Update nominal system angle + **mThetaSys = **mThetaSys + mTimeStep * mOmegaNom; + + // + if (mWithControl) + **mSourceValue_dq = mVSIController->step(**mVcap_dq, **mIfilter_dq); + + // Transformation interface backward + (**mSourceValue)(0, 0) = + Math::rotatingFrame2to1(**mSourceValue_dq, **mThetaSys, **mThetaInv); + + // set reference voltage of voltage source + if (!mModelAsCurrentSource) { + // pre-step of voltage source + **mSubCtrledVoltageSource->mVoltageRef = (**mSourceValue)(0, 0); + std::dynamic_pointer_cast(mSubCtrledVoltageSource) + ->mnaPreStep(time, timeStepCount); + } + + // stamp right side vector + mnaApplyRightSideVectorStamp(**mRightVector); } -void DP::Ph1::VSIVoltageControlDQ::mnaParentApplyRightSideVectorStamp(Matrix& rightVector) { - if (mModelAsCurrentSource) - Math::addToVectorElement(**mRightVector, mVirtualNodes[0]->matrixNodeIndex(), (**mSourceValue)(0,0)); +void DP::Ph1::VSIVoltageControlDQ::mnaParentApplyRightSideVectorStamp( + Matrix &rightVector) { + if (mModelAsCurrentSource) + Math::addToVectorElement(**mRightVector, + mVirtualNodes[0]->matrixNodeIndex(), + (**mSourceValue)(0, 0)); } -void DP::Ph1::VSIVoltageControlDQ::mnaParentAddPostStepDependencies(AttributeBase::List &prevStepDependencies, AttributeBase::List &attributeDependencies, AttributeBase::List &modifiedAttributes, Attribute::Ptr &leftVector) { - attributeDependencies.push_back(leftVector); - modifiedAttributes.push_back(mIntfVoltage); - modifiedAttributes.push_back(mIntfCurrent); +void DP::Ph1::VSIVoltageControlDQ::mnaParentAddPostStepDependencies( + AttributeBase::List &prevStepDependencies, + AttributeBase::List &attributeDependencies, + AttributeBase::List &modifiedAttributes, + Attribute::Ptr &leftVector) { + attributeDependencies.push_back(leftVector); + modifiedAttributes.push_back(mIntfVoltage); + modifiedAttributes.push_back(mIntfCurrent); } -void DP::Ph1::VSIVoltageControlDQ::mnaParentPostStep(Real time, Int timeStepCount, Attribute::Ptr &leftVector) { - mnaCompUpdateCurrent(**leftVector); - mnaCompUpdateVoltage(**leftVector); - updatePower(); +void DP::Ph1::VSIVoltageControlDQ::mnaParentPostStep( + Real time, Int timeStepCount, Attribute::Ptr &leftVector) { + mnaCompUpdateCurrent(**leftVector); + mnaCompUpdateVoltage(**leftVector); + updatePower(); } -void DP::Ph1::VSIVoltageControlDQ::mnaCompUpdateCurrent(const Matrix& leftvector) { - **mIntfCurrent = mSubCapacitorF->mIntfCurrent->get() + mSubFilterRL->mIntfCurrent->get(); +void DP::Ph1::VSIVoltageControlDQ::mnaCompUpdateCurrent( + const Matrix &leftvector) { + **mIntfCurrent = + mSubCapacitorF->mIntfCurrent->get() + mSubFilterRL->mIntfCurrent->get(); } -void DP::Ph1::VSIVoltageControlDQ::mnaCompUpdateVoltage(const Matrix& leftVector) { - (**mIntfVoltage)(0,0) = Math::complexFromVectorElement(leftVector, matrixNodeIndex(0)); +void DP::Ph1::VSIVoltageControlDQ::mnaCompUpdateVoltage( + const Matrix &leftVector) { + (**mIntfVoltage)(0, 0) = + Math::complexFromVectorElement(leftVector, matrixNodeIndex(0)); } void DP::Ph1::VSIVoltageControlDQ::updatePower() { - **mPower = (**mIntfVoltage)(0,0) * std::conj((**mIntfCurrent)(0,0)); + **mPower = (**mIntfVoltage)(0, 0) * std::conj((**mIntfCurrent)(0, 0)); } diff --git a/dpsim-models/src/DP/DP_Ph1_VSIVoltageControlVBR.cpp b/dpsim-models/src/DP/DP_Ph1_VSIVoltageControlVBR.cpp index e98d927401..6fb8b12d61 100644 --- a/dpsim-models/src/DP/DP_Ph1_VSIVoltageControlVBR.cpp +++ b/dpsim-models/src/DP/DP_Ph1_VSIVoltageControlVBR.cpp @@ -16,8 +16,7 @@ DP::Ph1::VSIVoltageControlVBR::VSIVoltageControlVBR(String uid, String name, Bool modelAsCurrentSource) : CompositePowerComp(uid, name, true, true, logLevel), VSIVoltageSourceInverterDQ(this->mSLog, mAttributes, - modelAsCurrentSource, false) -{ + modelAsCurrentSource, false) { setTerminalNumber(1); setVirtualNodeNumber(this->determineNumberOfVirtualNodes()); @@ -27,8 +26,7 @@ DP::Ph1::VSIVoltageControlVBR::VSIVoltageControlVBR(String uid, String name, mDqToComplexA = Matrix::Zero(2, 2); } -void DP::Ph1::VSIVoltageControlVBR::createSubComponents() -{ +void DP::Ph1::VSIVoltageControlVBR::createSubComponents() { // Capacitor as part of the LC filter mSubCapacitorF = DP::Ph1::Capacitor::make(**mName + "_CapF", mLogLevel); mSubCapacitorF->setParameters(mCf); @@ -37,8 +35,7 @@ void DP::Ph1::VSIVoltageControlVBR::createSubComponents() } void DP::Ph1::VSIVoltageControlVBR::initializeFromNodesAndTerminals( - Real frequency) -{ + Real frequency) { // terminal powers in consumer system -> convert to generator system **mPower = -terminal(0)->singlePower(); @@ -57,8 +54,7 @@ void DP::Ph1::VSIVoltageControlVBR::initializeFromNodesAndTerminals( // Connect & Initialize electrical subcomponents mSubCapacitorF->connect({SimNode::GND, mTerminals[0]->node()}); - for (auto subcomp : mSubComponents) - { + for (auto subcomp : mSubComponents) { subcomp->initialize(mFrequencies); subcomp->initializeFromNodesAndTerminals(frequency); } @@ -68,7 +64,8 @@ void DP::Ph1::VSIVoltageControlVBR::initializeFromNodesAndTerminals( /// initialize filter current in dp domain mFilterCurrent = Matrix::Zero(1, 1); - mFilterCurrent(0, 0) = Math::rotatingFrame2to1(**mIfilter_dq, **mThetaSys, **mThetaInv); + mFilterCurrent(0, 0) = + Math::rotatingFrame2to1(**mIfilter_dq, **mThetaSys, **mThetaInv); SPDLOG_LOGGER_INFO(mSLog, "\n--- Initialization from powerflow ---" @@ -82,23 +79,26 @@ void DP::Ph1::VSIVoltageControlVBR::initializeFromNodesAndTerminals( mSLog->flush(); } -void DP::Ph1::VSIVoltageControlVBR::initVars(Real timeStep) -{ - for (Int freq = 0; freq < mNumFreqs; freq++) - { +void DP::Ph1::VSIVoltageControlVBR::initVars(Real timeStep) { + for (Int freq = 0; freq < mNumFreqs; freq++) { Real a = timeStep / (2. * mLf); Real b = timeStep * 2. * PI * mFrequencies(freq, 0) / 2.; - Real equivCondReal = (a + mRf * std::pow(a, 2)) / (std::pow(1. + mRf * a, 2) + std::pow(b, 2)); + Real equivCondReal = (a + mRf * std::pow(a, 2)) / + (std::pow(1. + mRf * a, 2) + std::pow(b, 2)); Real equivCondImag = -a * b / (std::pow(1. + mRf * a, 2) + std::pow(b, 2)); mEquivCond(freq, 0) = {equivCondReal, equivCondImag}; - Real preCurrFracReal = (1. - std::pow(b, 2) + -std::pow(mRf * a, 2)) / (std::pow(1. + mRf * a, 2) + std::pow(b, 2)); - Real preCurrFracImag = (-2. * b) / (std::pow(1. + mRf * a, 2) + std::pow(b, 2)); + Real preCurrFracReal = (1. - std::pow(b, 2) + -std::pow(mRf * a, 2)) / + (std::pow(1. + mRf * a, 2) + std::pow(b, 2)); + Real preCurrFracImag = + (-2. * b) / (std::pow(1. + mRf * a, 2) + std::pow(b, 2)); mPrevCurrFac(freq, 0) = {preCurrFracReal, preCurrFracImag}; // TODO: implement it for mNumFreqs>1 - mEquivCurrent(freq, 0) = mEquivCond(freq, 0) * ((**mSourceValue)(0, 0) - (**mIntfVoltage)(0, 0)) + mPrevCurrFac(freq, 0) * mFilterCurrent(0, 0); + mEquivCurrent(freq, 0) = mEquivCond(freq, 0) * ((**mSourceValue)(0, 0) - + (**mIntfVoltage)(0, 0)) + + mPrevCurrFac(freq, 0) * mFilterCurrent(0, 0); } // initialize auxiliar @@ -117,12 +117,12 @@ void DP::Ph1::VSIVoltageControlVBR::initVars(Real timeStep) ->mE_VBR; mZ = Matrix::Zero(2, 2); - mZ << mEquivCond(0, 0).real(), -mEquivCond(0, 0).imag(), mEquivCond(0, 0).imag(), mEquivCond(0, 0).real(); + mZ << mEquivCond(0, 0).real(), -mEquivCond(0, 0).imag(), + mEquivCond(0, 0).imag(), mEquivCond(0, 0).real(); } void DP::Ph1::VSIVoltageControlVBR::mnaParentInitialize( - Real omega, Real timeStep, Attribute::Ptr leftVector) -{ + Real omega, Real timeStep, Attribute::Ptr leftVector) { mTimeStep = timeStep; // initialize RL Element matrix @@ -130,14 +130,12 @@ void DP::Ph1::VSIVoltageControlVBR::mnaParentInitialize( mEquivCond = MatrixComp::Zero(mNumFreqs, 1); mPrevCurrFac = MatrixComp::Zero(mNumFreqs, 1); - if (mWithControl) - { + if (mWithControl) { mVSIController->initialize(**mSourceValue_dq, **mVcap_dq, **mIfilter_dq, mTimeStep, mModelAsCurrentSource); mVSIController->calculateVBRconstants(); } - if (mWithDroopControl) - { + if (mWithDroopControl) { mDroopController->initialize(mTimeStep, **mOmega); } initVars(timeStep); @@ -177,15 +175,13 @@ void DP::Ph1::VSIVoltageControlVBR::mnaParentInitialize( SPDLOG_LOGGER_INFO(mSLog, "({}, {})", indexPair.first, indexPair.second); } -void DP::Ph1::VSIVoltageControlVBR::calculateResistanceMatrix() -{ +void DP::Ph1::VSIVoltageControlVBR::calculateResistanceMatrix() { mAMatrix = mDqToComplexA * mA * mComplexAToDq * mZ; mEMatrix = mDqToComplexA * mE * mComplexAToDq; } void DP::Ph1::VSIVoltageControlVBR::mnaParentApplySystemMatrixStamp( - SparseMatrixRow &systemMatrix) -{ + SparseMatrixRow &systemMatrix) { // TODO update_DqToComplexATransformMatrix(); mComplexAToDq = mDqToComplexA.transpose(); @@ -198,30 +194,28 @@ void DP::Ph1::VSIVoltageControlVBR::mnaParentApplySystemMatrixStamp( Math::addToMatrixElement(systemMatrix, matrixNodeIndex(0), matrixNodeIndex(0), mEquivCond(0, 0)); Math::addToMatrixElement(systemMatrix, matrixNodeIndex(0), - mVirtualNodes[0]->matrixNodeIndex(), -mEquivCond(0, 0)); + mVirtualNodes[0]->matrixNodeIndex(), + -mEquivCond(0, 0)); // Stamp voltage source equation Math::addToMatrixElement(systemMatrix, mVirtualNodes[0]->matrixNodeIndex(), mVirtualNodes[0]->matrixNodeIndex(), Matrix::Identity(2, 2) + mAMatrix); Math::addToMatrixElement(systemMatrix, mVirtualNodes[0]->matrixNodeIndex(), - matrixNodeIndex(0), - -mEMatrix - mAMatrix); + matrixNodeIndex(0), -mEMatrix - mAMatrix); } void DP::Ph1::VSIVoltageControlVBR::mnaParentAddPreStepDependencies( AttributeBase::List &prevStepDependencies, AttributeBase::List &attributeDependencies, - AttributeBase::List &modifiedAttributes) -{ + AttributeBase::List &modifiedAttributes) { modifiedAttributes.push_back(mRightVector); prevStepDependencies.push_back(mIntfVoltage); prevStepDependencies.push_back(mIntfCurrent); } void DP::Ph1::VSIVoltageControlVBR::mnaParentPreStep(Real time, - Int timeStepCount) -{ + Int timeStepCount) { // get measurements at time t=k **mVcap_dq = Math::rotatingFrame2to1((**mSubCapacitorF->mIntfVoltage)(0, 0), **mThetaInv, **mThetaSys); @@ -229,8 +223,7 @@ void DP::Ph1::VSIVoltageControlVBR::mnaParentPreStep(Real time, Math::rotatingFrame2to1(mFilterCurrent(0, 0), **mThetaInv, **mThetaSys); // Droop contoller step - if (mWithDroopControl) - { + if (mWithDroopControl) { mDroopController->step((**mPower).real()); } @@ -248,17 +241,21 @@ void DP::Ph1::VSIVoltageControlVBR::mnaParentPreStep(Real time, } void DP::Ph1::VSIVoltageControlVBR::mnaParentApplyRightSideVectorStamp( - Matrix &rightVector) -{ + Matrix &rightVector) { // TODO: implement it for mNumFreqs>1 - mEquivCurrent(0, 0) = mEquivCond(0, 0) * ((**mSourceValue)(0, 0) - (**mIntfVoltage)(0, 0)) + mPrevCurrFac(0, 0) * mFilterCurrent(0, 0); + mEquivCurrent(0, 0) = + mEquivCond(0, 0) * ((**mSourceValue)(0, 0) - (**mIntfVoltage)(0, 0)) + + mPrevCurrFac(0, 0) * mFilterCurrent(0, 0); // in complex: Complex dqToComplexA = Complex(mDqToComplexA(0, 0), -mDqToComplexA(0, 1)); - Complex eqVoltageSource = dqToComplexA * mVhist - std::dynamic_pointer_cast(mVSIController) - ->mA_VBR * - mEquivCurrent(0, 0); - Math::setVectorElement(rightVector, mVirtualNodes[0]->matrixNodeIndex(), eqVoltageSource); + Complex eqVoltageSource = + dqToComplexA * mVhist - + std::dynamic_pointer_cast(mVSIController) + ->mA_VBR * + mEquivCurrent(0, 0); + Math::setVectorElement(rightVector, mVirtualNodes[0]->matrixNodeIndex(), + eqVoltageSource); // stamp equivalent current RL element // TODO: Implement function addToVectorElement for multiple frequencies @@ -267,15 +264,15 @@ void DP::Ph1::VSIVoltageControlVBR::mnaParentApplyRightSideVectorStamp( // Math::addToVectorElement(rightVector, matrixNodeIndex(0), mEquivCurrent(freq, 0), 1, 0, freq); //} // - Math::addToVectorElement(rightVector, matrixNodeIndex(0), mEquivCurrent(0, 0)); + Math::addToVectorElement(rightVector, matrixNodeIndex(0), + mEquivCurrent(0, 0)); } void DP::Ph1::VSIVoltageControlVBR::mnaParentAddPostStepDependencies( AttributeBase::List &prevStepDependencies, AttributeBase::List &attributeDependencies, AttributeBase::List &modifiedAttributes, - Attribute::Ptr &leftVector) -{ + Attribute::Ptr &leftVector) { attributeDependencies.push_back(leftVector); attributeDependencies.push_back(mRightVector); modifiedAttributes.push_back(mIntfVoltage); @@ -283,16 +280,14 @@ void DP::Ph1::VSIVoltageControlVBR::mnaParentAddPostStepDependencies( } void DP::Ph1::VSIVoltageControlVBR::mnaParentPostStep( - Real time, Int timeStepCount, Attribute::Ptr &leftVector) -{ + Real time, Int timeStepCount, Attribute::Ptr &leftVector) { mnaCompUpdateVoltage(**leftVector); mnaCompUpdateCurrent(**leftVector); updatePower(); } void DP::Ph1::VSIVoltageControlVBR::mnaCompUpdateVoltage( - const Matrix &leftVector) -{ + const Matrix &leftVector) { (**mIntfVoltage)(0, 0) = Math::complexFromVectorElement(leftVector, matrixNodeIndex(0)); (**mSourceValue)(0, 0) = Math::complexFromVectorElement( @@ -300,21 +295,20 @@ void DP::Ph1::VSIVoltageControlVBR::mnaCompUpdateVoltage( } void DP::Ph1::VSIVoltageControlVBR::mnaCompUpdateCurrent( - const Matrix &leftvector) -{ + const Matrix &leftvector) { // TODO: calculations for number_frequencies>1 - mFilterCurrent(0, 0) = mEquivCond(0, 0) * ((**mSourceValue)(0, 0) - (**mIntfVoltage)(0, 0)) + mEquivCurrent(0, 0); - **mIntfCurrent = - mSubCapacitorF->mIntfCurrent->get() + mFilterCurrent; + mFilterCurrent(0, 0) = + mEquivCond(0, 0) * ((**mSourceValue)(0, 0) - (**mIntfVoltage)(0, 0)) + + mEquivCurrent(0, 0); + **mIntfCurrent = mSubCapacitorF->mIntfCurrent->get() + mFilterCurrent; } -void DP::Ph1::VSIVoltageControlVBR::updatePower() -{ +void DP::Ph1::VSIVoltageControlVBR::updatePower() { **mPower = (**mIntfVoltage)(0, 0) * std::conj((**mIntfCurrent)(0, 0)); } -void DP::Ph1::VSIVoltageControlVBR::update_DqToComplexATransformMatrix() -{ - mDqToComplexA << cos(**mThetaInv - **mThetaSys), -sin(**mThetaInv - **mThetaSys), - sin(**mThetaInv - **mThetaSys), cos(**mThetaInv - **mThetaSys); +void DP::Ph1::VSIVoltageControlVBR::update_DqToComplexATransformMatrix() { + mDqToComplexA << cos(**mThetaInv - **mThetaSys), + -sin(**mThetaInv - **mThetaSys), sin(**mThetaInv - **mThetaSys), + cos(**mThetaInv - **mThetaSys); } \ No newline at end of file diff --git a/dpsim-models/src/EMT/EMT_Ph3_VSIVoltageControlDQ.cpp b/dpsim-models/src/EMT/EMT_Ph3_VSIVoltageControlDQ.cpp index f3d4753722..6e8b425b1f 100644 --- a/dpsim-models/src/EMT/EMT_Ph3_VSIVoltageControlDQ.cpp +++ b/dpsim-models/src/EMT/EMT_Ph3_VSIVoltageControlDQ.cpp @@ -21,6 +21,9 @@ EMT::Ph3::VSIVoltageControlDQ::VSIVoltageControlDQ(String uid, String name, mPhaseType = PhaseType::ABC; setTerminalNumber(1); setVirtualNodeNumber(this->determineNumberOfVirtualNodes()); + SPDLOG_LOGGER_INFO(mSLog, "\nthis->determineNumberOfVirtualNodes()={}", + this->determineNumberOfVirtualNodes()); + mSLog->flush(); } void EMT::Ph3::VSIVoltageControlDQ::createSubComponents() { @@ -96,6 +99,13 @@ void EMT::Ph3::VSIVoltageControlDQ::initializeFromNodesAndTerminals( **mSourceValue = inverseParkTransformPowerInvariant(**mThetaInv, **mSourceValue_dq).real(); + // Create & Initialize electrical subcomponents + this->connectSubComponents(); + for (auto subcomp : mSubComponents) { + subcomp->initialize(mFrequencies); + subcomp->initializeFromNodesAndTerminals(frequency); + } + // droop **mOmega = mOmegaNom; @@ -112,7 +122,9 @@ void EMT::Ph3::VSIVoltageControlDQ::initializeFromNodesAndTerminals( void EMT::Ph3::VSIVoltageControlDQ::mnaParentInitialize( Real omega, Real timeStep, Attribute::Ptr leftVector) { - this->updateMatrixNodeIndices(); + + SPDLOG_LOGGER_INFO(mSLog, "\n--- TEST ---"); + mSLog->flush(); mTimeStep = timeStep; if (mWithControl) mVSIController->initialize(**mSourceValue_dq, **mVcap_dq, **mIfilter_dq, @@ -123,8 +135,13 @@ void EMT::Ph3::VSIVoltageControlDQ::mnaParentAddPreStepDependencies( AttributeBase::List &prevStepDependencies, AttributeBase::List &attributeDependencies, AttributeBase::List &modifiedAttributes) { + modifiedAttributes.push_back(mRightVector); prevStepDependencies.push_back(mIntfVoltage); + prevStepDependencies.push_back(mIntfCurrent); + + SPDLOG_LOGGER_INFO(mSLog, "\n--- TEST1 ---"); + mSLog->flush(); } void EMT::Ph3::VSIVoltageControlDQ::mnaParentPreStep(Real time, diff --git a/dpsim-models/src/EMT/EMT_Ph3_VSIVoltageControlVBR.cpp b/dpsim-models/src/EMT/EMT_Ph3_VSIVoltageControlVBR.cpp index 2f8ae21b35..bbc1dfcc6f 100644 --- a/dpsim-models/src/EMT/EMT_Ph3_VSIVoltageControlVBR.cpp +++ b/dpsim-models/src/EMT/EMT_Ph3_VSIVoltageControlVBR.cpp @@ -11,314 +11,417 @@ using namespace CPS; -EMT::Ph3::VSIVoltageControlVBR::VSIVoltageControlVBR(String uid, String name, Logger::Level logLevel, - Bool modelAsCurrentSource) : CompositePowerComp(uid, name, true, true, logLevel), - VSIVoltageSourceInverterDQ(this->mSLog, mAttributes, modelAsCurrentSource, false) -{ - - mPhaseType = PhaseType::ABC; - setTerminalNumber(1); - setVirtualNodeNumber(this->determineNumberOfVirtualNodes()); - - **mIntfVoltage = Matrix::Zero(3, 1); - **mIntfCurrent = Matrix::Zero(3, 1); +EMT::Ph3::VSIVoltageControlVBR::VSIVoltageControlVBR(String uid, String name, + Logger::Level logLevel, + Bool modelAsCurrentSource) + : CompositePowerComp(uid, name, true, true, logLevel), + VSIVoltageSourceInverterDQ(this->mSLog, mAttributes, modelAsCurrentSource, + false) { + + mPhaseType = PhaseType::ABC; + setTerminalNumber(1); + setVirtualNodeNumber(this->determineNumberOfVirtualNodes()); + + **mIntfVoltage = Matrix::Zero(3, 1); + **mIntfCurrent = Matrix::Zero(3, 1); } -void EMT::Ph3::VSIVoltageControlVBR::createSubComponents() -{ - // Capacitor as part of the LC filter - mSubCapacitorF = EMT::Ph3::Capacitor::make(**mName + "_CapF", mLogLevel); - addMNASubComponent(mSubCapacitorF, MNA_SUBCOMP_TASK_ORDER::TASK_BEFORE_PARENT, MNA_SUBCOMP_TASK_ORDER::TASK_BEFORE_PARENT, true); - mSubCapacitorF->setParameters(Math::singlePhaseParameterToThreePhase(mCf)); +void EMT::Ph3::VSIVoltageControlVBR::createSubComponents() { + // Capacitor as part of the LC filter + mSubCapacitorF = EMT::Ph3::Capacitor::make(**mName + "_CapF", mLogLevel); + addMNASubComponent(mSubCapacitorF, MNA_SUBCOMP_TASK_ORDER::TASK_BEFORE_PARENT, + MNA_SUBCOMP_TASK_ORDER::TASK_BEFORE_PARENT, true); + mSubCapacitorF->setParameters(Math::singlePhaseParameterToThreePhase(mCf)); } -void EMT::Ph3::VSIVoltageControlVBR::initializeFromNodesAndTerminals(Real frequency) -{ - // terminal powers in consumer system -> convert to generator system - **mPower = -terminal(0)->singlePower(); - - // set initial interface quantities --> Current flowing into the inverter is positive - Complex intfVoltageComplex = initialSingleVoltage(0); - Complex intfCurrentComplex = std::conj(**mPower / intfVoltageComplex); - **mIntfVoltage = Math::singlePhaseVariableToThreePhase(RMS3PH_TO_PEAK1PH * intfVoltageComplex).real(); - **mIntfCurrent = Math::singlePhaseVariableToThreePhase(RMS3PH_TO_PEAK1PH * intfCurrentComplex).real(); - - // initialize filter variables and set initial voltage of virtual nodes - initializeFilterVariables(intfVoltageComplex, intfCurrentComplex, mVirtualNodes); - - // calculate initial source value - // FIXME: later is mSourceValue used to store the history term - **mSourceValue = inverseParkTransformPowerInvariant(**mThetaInv, **mSourceValue_dq).real(); - - // Connect & Initialize electrical subcomponents - mSubCapacitorF->connect({SimNode::GND, mTerminals[0]->node()}); - for (auto subcomp : mSubComponents) - { - subcomp->initialize(mFrequencies); - subcomp->initializeFromNodesAndTerminals(frequency); - } - - // droop - **mOmega = mOmegaNom; - - /// initialize filter current in dp domain - mFilterCurrent = inverseParkTransformPowerInvariant(**mThetaInv, **mIfilter_dq); - - SPDLOG_LOGGER_INFO(mSLog, - "\n--- Initialization from powerflow ---" - "\nTerminal 0 connected to {} = sim node {}" - "\nInverter terminal voltage: {}[V]" - "\nInverter output current: {}[A]", - mTerminals[0]->node()->name(), mTerminals[0]->node()->matrixNodeIndex(), - (**mIntfVoltage)(0, 0), - (**mIntfCurrent)(0, 0)); - mSLog->flush(); +void EMT::Ph3::VSIVoltageControlVBR::initializeFromNodesAndTerminals( + Real frequency) { + // terminal powers in consumer system -> convert to generator system + **mPower = -terminal(0)->singlePower(); + + // set initial interface quantities --> Current flowing into the inverter is positive + Complex intfVoltageComplex = initialSingleVoltage(0); + Complex intfCurrentComplex = std::conj(**mPower / intfVoltageComplex); + **mIntfVoltage = Math::singlePhaseVariableToThreePhase(RMS3PH_TO_PEAK1PH * + intfVoltageComplex) + .real(); + **mIntfCurrent = Math::singlePhaseVariableToThreePhase(RMS3PH_TO_PEAK1PH * + intfCurrentComplex) + .real(); + + // initialize filter variables and set initial voltage of virtual nodes + initializeFilterVariables(intfVoltageComplex, intfCurrentComplex, + mVirtualNodes); + + // calculate initial source value + // FIXME: later is mSourceValue used to store the history term + **mSourceValue = + inverseParkTransformPowerInvariant(**mThetaInv, **mSourceValue_dq).real(); + + // Connect & Initialize electrical subcomponents + mSubCapacitorF->connect({SimNode::GND, mTerminals[0]->node()}); + for (auto subcomp : mSubComponents) { + subcomp->initialize(mFrequencies); + subcomp->initializeFromNodesAndTerminals(frequency); + } + + // droop + **mOmega = mOmegaNom; + + /// initialize filter current in dp domain + mFilterCurrent = + inverseParkTransformPowerInvariant(**mThetaInv, **mIfilter_dq); + + SPDLOG_LOGGER_INFO(mSLog, + "\n--- Initialization from powerflow ---" + "\nTerminal 0 connected to {} = sim node {}" + "\nInverter terminal voltage: {}[V]" + "\nInverter output current: {}[A]", + mTerminals[0]->node()->name(), + mTerminals[0]->node()->matrixNodeIndex(), + (**mIntfVoltage)(0, 0), (**mIntfCurrent)(0, 0)); + mSLog->flush(); } -void EMT::Ph3::VSIVoltageControlVBR::initVars(Real timeStep) -{ - // Assumption: symmetric R and L matrix - Real a = timeStep * mRf / (2. * mLf); - Real b = timeStep / (2. * mLf); - - mEquivCond = b / (1. + a); - mPrevCurrFac = (1. - a) / (1. + a); - - // - mEquivCurrent = mEquivCond * ((**mSourceValue).real() - **mIntfVoltage) + mPrevCurrFac * mFilterCurrent; - - // initialize auxiliar - mA = Matrix::Zero(2, 2); - mA << std::dynamic_pointer_cast(mVSIController) - ->mA_VBR, - 0, 0, - std::dynamic_pointer_cast(mVSIController) - ->mA_VBR; - - mE = Matrix::Zero(2, 2); - mE << std::dynamic_pointer_cast(mVSIController) - ->mE_VBR, - 0, 0, - std::dynamic_pointer_cast(mVSIController) - ->mE_VBR; +void EMT::Ph3::VSIVoltageControlVBR::initVars(Real timeStep) { + // Assumption: symmetric R and L matrix + Real a = timeStep * mRf / (2. * mLf); + Real b = timeStep / (2. * mLf); + + mEquivCond = b / (1. + a); + mPrevCurrFac = (1. - a) / (1. + a); + + // + mEquivCurrent = mEquivCond * ((**mSourceValue).real() - **mIntfVoltage) + + mPrevCurrFac * mFilterCurrent; + + // initialize auxiliar + mA = Matrix::Zero(2, 2); + mA << std::dynamic_pointer_cast(mVSIController) + ->mA_VBR, + 0, 0, + std::dynamic_pointer_cast(mVSIController) + ->mA_VBR; + + mE = Matrix::Zero(2, 2); + mE << std::dynamic_pointer_cast(mVSIController) + ->mE_VBR, + 0, 0, + std::dynamic_pointer_cast(mVSIController) + ->mE_VBR; } -void EMT::Ph3::VSIVoltageControlVBR::mnaParentInitialize(Real omega, Real timeStep, Attribute::Ptr leftVector) -{ - mTimeStep = timeStep; - if (mWithControl) - { - mVSIController->initialize(**mSourceValue_dq, **mVcap_dq, **mIfilter_dq, mTimeStep, mModelAsCurrentSource); - mVSIController->calculateVBRconstants(); - } - initVars(timeStep); - - // get matrix dimension to properly set variable entries - // upper left block - mVariableSystemMatrixEntries.push_back(std::make_pair(mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), mVirtualNodes[0]->matrixNodeIndex(PhaseType::A))); - mVariableSystemMatrixEntries.push_back(std::make_pair(mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), mVirtualNodes[0]->matrixNodeIndex(PhaseType::B))); - mVariableSystemMatrixEntries.push_back(std::make_pair(mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), mVirtualNodes[0]->matrixNodeIndex(PhaseType::C))); - mVariableSystemMatrixEntries.push_back(std::make_pair(mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), mVirtualNodes[0]->matrixNodeIndex(PhaseType::A))); - mVariableSystemMatrixEntries.push_back(std::make_pair(mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), mVirtualNodes[0]->matrixNodeIndex(PhaseType::B))); - mVariableSystemMatrixEntries.push_back(std::make_pair(mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), mVirtualNodes[0]->matrixNodeIndex(PhaseType::C))); - mVariableSystemMatrixEntries.push_back(std::make_pair(mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), mVirtualNodes[0]->matrixNodeIndex(PhaseType::A))); - mVariableSystemMatrixEntries.push_back(std::make_pair(mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), mVirtualNodes[0]->matrixNodeIndex(PhaseType::B))); - mVariableSystemMatrixEntries.push_back(std::make_pair(mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), mVirtualNodes[0]->matrixNodeIndex(PhaseType::C))); - - // off diagonal blocks - mVariableSystemMatrixEntries.push_back(std::make_pair(mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), matrixNodeIndex(0, 0))); - mVariableSystemMatrixEntries.push_back(std::make_pair(mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), matrixNodeIndex(0, 1))); - mVariableSystemMatrixEntries.push_back(std::make_pair(mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), matrixNodeIndex(0, 2))); - mVariableSystemMatrixEntries.push_back(std::make_pair(mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), matrixNodeIndex(0, 0))); - mVariableSystemMatrixEntries.push_back(std::make_pair(mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), matrixNodeIndex(0, 1))); - mVariableSystemMatrixEntries.push_back(std::make_pair(mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), matrixNodeIndex(0, 2))); - mVariableSystemMatrixEntries.push_back(std::make_pair(mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), matrixNodeIndex(0, 0))); - mVariableSystemMatrixEntries.push_back(std::make_pair(mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), matrixNodeIndex(0, 1))); - mVariableSystemMatrixEntries.push_back(std::make_pair(mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), matrixNodeIndex(0, 2))); - - SPDLOG_LOGGER_INFO(mSLog, "List of index pairs of varying matrix entries: "); - for (auto indexPair : mVariableSystemMatrixEntries) - SPDLOG_LOGGER_INFO(mSLog, "({}, {})", indexPair.first, indexPair.second); - mSLog->flush(); +void EMT::Ph3::VSIVoltageControlVBR::mnaParentInitialize( + Real omega, Real timeStep, Attribute::Ptr leftVector) { + mTimeStep = timeStep; + if (mWithControl) { + mVSIController->initialize(**mSourceValue_dq, **mVcap_dq, **mIfilter_dq, + mTimeStep, mModelAsCurrentSource); + mVSIController->calculateVBRconstants(); + } + initVars(timeStep); + + // get matrix dimension to properly set variable entries + // upper left block + mVariableSystemMatrixEntries.push_back(std::make_pair( + mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), + mVirtualNodes[0]->matrixNodeIndex(PhaseType::A))); + mVariableSystemMatrixEntries.push_back(std::make_pair( + mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), + mVirtualNodes[0]->matrixNodeIndex(PhaseType::B))); + mVariableSystemMatrixEntries.push_back(std::make_pair( + mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), + mVirtualNodes[0]->matrixNodeIndex(PhaseType::C))); + mVariableSystemMatrixEntries.push_back(std::make_pair( + mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), + mVirtualNodes[0]->matrixNodeIndex(PhaseType::A))); + mVariableSystemMatrixEntries.push_back(std::make_pair( + mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), + mVirtualNodes[0]->matrixNodeIndex(PhaseType::B))); + mVariableSystemMatrixEntries.push_back(std::make_pair( + mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), + mVirtualNodes[0]->matrixNodeIndex(PhaseType::C))); + mVariableSystemMatrixEntries.push_back(std::make_pair( + mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), + mVirtualNodes[0]->matrixNodeIndex(PhaseType::A))); + mVariableSystemMatrixEntries.push_back(std::make_pair( + mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), + mVirtualNodes[0]->matrixNodeIndex(PhaseType::B))); + mVariableSystemMatrixEntries.push_back(std::make_pair( + mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), + mVirtualNodes[0]->matrixNodeIndex(PhaseType::C))); + + // off diagonal blocks + mVariableSystemMatrixEntries.push_back(std::make_pair( + mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), matrixNodeIndex(0, 0))); + mVariableSystemMatrixEntries.push_back(std::make_pair( + mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), matrixNodeIndex(0, 1))); + mVariableSystemMatrixEntries.push_back(std::make_pair( + mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), matrixNodeIndex(0, 2))); + mVariableSystemMatrixEntries.push_back(std::make_pair( + mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), matrixNodeIndex(0, 0))); + mVariableSystemMatrixEntries.push_back(std::make_pair( + mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), matrixNodeIndex(0, 1))); + mVariableSystemMatrixEntries.push_back(std::make_pair( + mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), matrixNodeIndex(0, 2))); + mVariableSystemMatrixEntries.push_back(std::make_pair( + mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), matrixNodeIndex(0, 0))); + mVariableSystemMatrixEntries.push_back(std::make_pair( + mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), matrixNodeIndex(0, 1))); + mVariableSystemMatrixEntries.push_back(std::make_pair( + mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), matrixNodeIndex(0, 2))); + + SPDLOG_LOGGER_INFO(mSLog, "List of index pairs of varying matrix entries: "); + for (auto indexPair : mVariableSystemMatrixEntries) + SPDLOG_LOGGER_INFO(mSLog, "({}, {})", indexPair.first, indexPair.second); + mSLog->flush(); } -void EMT::Ph3::VSIVoltageControlVBR::calculateResistanceMatrix() -{ - mAMatrix = mDqToABC * mA * mABCToDq * mEquivCond; - mAMatrix_ = Matrix::Identity(3, 3) + mAMatrix; - mEMatrix = mDqToABC * mE * mABCToDq; +void EMT::Ph3::VSIVoltageControlVBR::calculateResistanceMatrix() { + mAMatrix = mDqToABC * mA * mABCToDq * mEquivCond; + mAMatrix_ = Matrix::Identity(3, 3) + mAMatrix; + mEMatrix = mDqToABC * mE * mABCToDq; } void EMT::Ph3::VSIVoltageControlVBR::mnaParentApplySystemMatrixStamp( - SparseMatrixRow &systemMatrix) -{ - // - mABCToDq = getParkTransformMatrixPowerInvariant(**mThetaInv); - mDqToABC = getInverseParkTransformMatrixPowerInvariant(**mThetaInv); - - // calculate resistance matrix at t=k+1 - calculateResistanceMatrix(); - - // Stamp filter current: I = (Vc - Vsource) * Admittance - Math::addToMatrixElement(systemMatrix, matrixNodeIndex(0, 0), matrixNodeIndex(0, 0), mEquivCond); - Math::addToMatrixElement(systemMatrix, matrixNodeIndex(0, 1), matrixNodeIndex(0, 1), mEquivCond); - Math::addToMatrixElement(systemMatrix, matrixNodeIndex(0, 2), matrixNodeIndex(0, 2), mEquivCond); - Math::addToMatrixElement(systemMatrix, matrixNodeIndex(0, 0), mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), -mEquivCond); - Math::addToMatrixElement(systemMatrix, matrixNodeIndex(0, 1), mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), -mEquivCond); - Math::addToMatrixElement(systemMatrix, matrixNodeIndex(0, 2), mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), -mEquivCond); - - // Stamp voltage source equation - Math::addToMatrixElement(systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), mAMatrix_(0, 0)); - Math::addToMatrixElement(systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), mAMatrix_(0, 1)); - Math::addToMatrixElement(systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), mAMatrix_(0, 2)); - Math::addToMatrixElement(systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), mAMatrix_(1, 0)); - Math::addToMatrixElement(systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), mAMatrix_(1, 1)); - Math::addToMatrixElement(systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), mAMatrix_(1, 2)); - Math::addToMatrixElement(systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), mAMatrix_(2, 0)); - Math::addToMatrixElement(systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), mAMatrix_(2, 1)); - Math::addToMatrixElement(systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), mAMatrix_(2, 2)); - - Math::addToMatrixElement(systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), matrixNodeIndex(0, 0), -mEMatrix(0, 0) - mAMatrix(0, 0)); - Math::addToMatrixElement(systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), matrixNodeIndex(0, 1), -mEMatrix(0, 1) - mAMatrix(0, 1)); - Math::addToMatrixElement(systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), matrixNodeIndex(0, 2), -mEMatrix(0, 2) - mAMatrix(0, 2)); - Math::addToMatrixElement(systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), matrixNodeIndex(0, 0), -mEMatrix(1, 0) - mAMatrix(1, 0)); - Math::addToMatrixElement(systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), matrixNodeIndex(0, 1), -mEMatrix(1, 1) - mAMatrix(1, 1)); - Math::addToMatrixElement(systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), matrixNodeIndex(0, 2), -mEMatrix(1, 2) - mAMatrix(1, 2)); - Math::addToMatrixElement(systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), matrixNodeIndex(0, 0), -mEMatrix(2, 0) - mAMatrix(2, 0)); - Math::addToMatrixElement(systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), matrixNodeIndex(0, 1), -mEMatrix(2, 1) - mAMatrix(2, 1)); - Math::addToMatrixElement(systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), matrixNodeIndex(0, 2), -mEMatrix(2, 2) - mAMatrix(2, 2)); + SparseMatrixRow &systemMatrix) { + // + mABCToDq = getParkTransformMatrixPowerInvariant(**mThetaInv); + mDqToABC = getInverseParkTransformMatrixPowerInvariant(**mThetaInv); + + // calculate resistance matrix at t=k+1 + calculateResistanceMatrix(); + + // Stamp filter current: I = (Vc - Vsource) * Admittance + Math::addToMatrixElement(systemMatrix, matrixNodeIndex(0, 0), + matrixNodeIndex(0, 0), mEquivCond); + Math::addToMatrixElement(systemMatrix, matrixNodeIndex(0, 1), + matrixNodeIndex(0, 1), mEquivCond); + Math::addToMatrixElement(systemMatrix, matrixNodeIndex(0, 2), + matrixNodeIndex(0, 2), mEquivCond); + Math::addToMatrixElement(systemMatrix, matrixNodeIndex(0, 0), + mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), + -mEquivCond); + Math::addToMatrixElement(systemMatrix, matrixNodeIndex(0, 1), + mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), + -mEquivCond); + Math::addToMatrixElement(systemMatrix, matrixNodeIndex(0, 2), + mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), + -mEquivCond); + + // Stamp voltage source equation + Math::addToMatrixElement( + systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), + mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), mAMatrix_(0, 0)); + Math::addToMatrixElement( + systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), + mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), mAMatrix_(0, 1)); + Math::addToMatrixElement( + systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), + mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), mAMatrix_(0, 2)); + Math::addToMatrixElement( + systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), + mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), mAMatrix_(1, 0)); + Math::addToMatrixElement( + systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), + mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), mAMatrix_(1, 1)); + Math::addToMatrixElement( + systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), + mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), mAMatrix_(1, 2)); + Math::addToMatrixElement( + systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), + mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), mAMatrix_(2, 0)); + Math::addToMatrixElement( + systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), + mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), mAMatrix_(2, 1)); + Math::addToMatrixElement( + systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), + mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), mAMatrix_(2, 2)); + + Math::addToMatrixElement( + systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), + matrixNodeIndex(0, 0), -mEMatrix(0, 0) - mAMatrix(0, 0)); + Math::addToMatrixElement( + systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), + matrixNodeIndex(0, 1), -mEMatrix(0, 1) - mAMatrix(0, 1)); + Math::addToMatrixElement( + systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::A), + matrixNodeIndex(0, 2), -mEMatrix(0, 2) - mAMatrix(0, 2)); + Math::addToMatrixElement( + systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), + matrixNodeIndex(0, 0), -mEMatrix(1, 0) - mAMatrix(1, 0)); + Math::addToMatrixElement( + systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), + matrixNodeIndex(0, 1), -mEMatrix(1, 1) - mAMatrix(1, 1)); + Math::addToMatrixElement( + systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::B), + matrixNodeIndex(0, 2), -mEMatrix(1, 2) - mAMatrix(1, 2)); + Math::addToMatrixElement( + systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), + matrixNodeIndex(0, 0), -mEMatrix(2, 0) - mAMatrix(2, 0)); + Math::addToMatrixElement( + systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), + matrixNodeIndex(0, 1), -mEMatrix(2, 1) - mAMatrix(2, 1)); + Math::addToMatrixElement( + systemMatrix, mVirtualNodes[0]->matrixNodeIndex(PhaseType::C), + matrixNodeIndex(0, 2), -mEMatrix(2, 2) - mAMatrix(2, 2)); } -void EMT::Ph3::VSIVoltageControlVBR::mnaParentAddPreStepDependencies(AttributeBase::List &prevStepDependencies, AttributeBase::List &attributeDependencies, AttributeBase::List &modifiedAttributes) -{ - modifiedAttributes.push_back(mRightVector); - prevStepDependencies.push_back(mIntfVoltage); - prevStepDependencies.push_back(mIntfCurrent); +void EMT::Ph3::VSIVoltageControlVBR::mnaParentAddPreStepDependencies( + AttributeBase::List &prevStepDependencies, + AttributeBase::List &attributeDependencies, + AttributeBase::List &modifiedAttributes) { + modifiedAttributes.push_back(mRightVector); + prevStepDependencies.push_back(mIntfVoltage); + prevStepDependencies.push_back(mIntfCurrent); } // TODO -void EMT::Ph3::VSIVoltageControlVBR::mnaParentPreStep(Real time, Int timeStepCount) -{ - // Transformation interface forward - **mVcap_dq = parkTransformPowerInvariant(**mThetaInv, **mSubCapacitorF->mIntfVoltage); - **mIfilter_dq = parkTransformPowerInvariant(**mThetaInv, mFilterCurrent); +void EMT::Ph3::VSIVoltageControlVBR::mnaParentPreStep(Real time, + Int timeStepCount) { + // Transformation interface forward + **mVcap_dq = + parkTransformPowerInvariant(**mThetaInv, **mSubCapacitorF->mIntfVoltage); + **mIfilter_dq = parkTransformPowerInvariant(**mThetaInv, mFilterCurrent); - // TODO: droop - // if (mWithDroop) - // mDroop->signalStep(time, timeStepCount); + // TODO: droop + // if (mWithDroop) + // mDroop->signalStep(time, timeStepCount); - // Update nominal system angle - **mThetaSys = **mThetaSys + mTimeStep * mOmegaNom; + // Update nominal system angle + **mThetaSys = **mThetaSys + mTimeStep * mOmegaNom; - // VCO Step - **mThetaInv = **mThetaInv + mTimeStep * **mOmega; + // VCO Step + **mThetaInv = **mThetaInv + mTimeStep * **mOmega; - // - if (mWithControl) - mVhist_dq = mVSIController->step(**mVcap_dq, **mIfilter_dq); + // + if (mWithControl) + mVhist_dq = mVSIController->step(**mVcap_dq, **mIfilter_dq); - // Transformation interface backward - mVhist = inverseParkTransformPowerInvariant(**mThetaInv, mVhist_dq); + // Transformation interface backward + mVhist = inverseParkTransformPowerInvariant(**mThetaInv, mVhist_dq); - // stamp right side vector - mnaApplyRightSideVectorStamp(**mRightVector); + // stamp right side vector + mnaApplyRightSideVectorStamp(**mRightVector); } -void EMT::Ph3::VSIVoltageControlVBR::mnaParentApplyRightSideVectorStamp(Matrix &rightVector) -{ - mEquivCurrent = mEquivCond * ((**mSourceValue).real() - **mIntfVoltage) + mPrevCurrFac * mFilterCurrent; - Matrix res = mVhist - mAMatrix * mEquivCurrent; - Math::addToVectorElement(**mRightVector, mVirtualNodes[0]->matrixNodeIndex(CPS::PhaseType::A), res(0, 0)); - Math::addToVectorElement(**mRightVector, mVirtualNodes[0]->matrixNodeIndex(CPS::PhaseType::B), res(1, 0)); - Math::addToVectorElement(**mRightVector, mVirtualNodes[0]->matrixNodeIndex(CPS::PhaseType::C), res(2, 0)); - - Math::addToVectorElement(**mRightVector, matrixNodeIndex(0, 0), mEquivCurrent(0, 0)); - Math::addToVectorElement(**mRightVector, matrixNodeIndex(0, 1), mEquivCurrent(1, 0)); - Math::addToVectorElement(**mRightVector, matrixNodeIndex(0, 2), mEquivCurrent(2, 0)); +void EMT::Ph3::VSIVoltageControlVBR::mnaParentApplyRightSideVectorStamp( + Matrix &rightVector) { + mEquivCurrent = mEquivCond * ((**mSourceValue).real() - **mIntfVoltage) + + mPrevCurrFac * mFilterCurrent; + Matrix res = mVhist - mAMatrix * mEquivCurrent; + Math::addToVectorElement(**mRightVector, + mVirtualNodes[0]->matrixNodeIndex(CPS::PhaseType::A), + res(0, 0)); + Math::addToVectorElement(**mRightVector, + mVirtualNodes[0]->matrixNodeIndex(CPS::PhaseType::B), + res(1, 0)); + Math::addToVectorElement(**mRightVector, + mVirtualNodes[0]->matrixNodeIndex(CPS::PhaseType::C), + res(2, 0)); + + Math::addToVectorElement(**mRightVector, matrixNodeIndex(0, 0), + mEquivCurrent(0, 0)); + Math::addToVectorElement(**mRightVector, matrixNodeIndex(0, 1), + mEquivCurrent(1, 0)); + Math::addToVectorElement(**mRightVector, matrixNodeIndex(0, 2), + mEquivCurrent(2, 0)); } -void EMT::Ph3::VSIVoltageControlVBR::mnaParentAddPostStepDependencies(AttributeBase::List &prevStepDependencies, AttributeBase::List &attributeDependencies, AttributeBase::List &modifiedAttributes, Attribute::Ptr &leftVector) -{ - attributeDependencies.push_back(leftVector); - attributeDependencies.push_back(mRightVector); - modifiedAttributes.push_back(mIntfVoltage); - modifiedAttributes.push_back(mIntfCurrent); +void EMT::Ph3::VSIVoltageControlVBR::mnaParentAddPostStepDependencies( + AttributeBase::List &prevStepDependencies, + AttributeBase::List &attributeDependencies, + AttributeBase::List &modifiedAttributes, + Attribute::Ptr &leftVector) { + attributeDependencies.push_back(leftVector); + attributeDependencies.push_back(mRightVector); + modifiedAttributes.push_back(mIntfVoltage); + modifiedAttributes.push_back(mIntfCurrent); } -void EMT::Ph3::VSIVoltageControlVBR::mnaParentPostStep(Real time, Int timeStepCount, Attribute::Ptr &leftVector) -{ - mnaCompUpdateVoltage(**leftVector); - mnaCompUpdateCurrent(**leftVector); - updatePower(); +void EMT::Ph3::VSIVoltageControlVBR::mnaParentPostStep( + Real time, Int timeStepCount, Attribute::Ptr &leftVector) { + mnaCompUpdateVoltage(**leftVector); + mnaCompUpdateCurrent(**leftVector); + updatePower(); } -void EMT::Ph3::VSIVoltageControlVBR::mnaCompUpdateVoltage(const Matrix &leftVector) -{ - (**mIntfVoltage)(0, 0) = Math::realFromVectorElement(leftVector, matrixNodeIndex(0, 0)); - (**mIntfVoltage)(1, 0) = Math::realFromVectorElement(leftVector, matrixNodeIndex(0, 1)); - (**mIntfVoltage)(2, 0) = Math::realFromVectorElement(leftVector, matrixNodeIndex(0, 2)); - (**mSourceValue)(0, 0) = Math::realFromVectorElement( - leftVector, mVirtualNodes[0]->matrixNodeIndex(PhaseType::A)); - (**mSourceValue)(1, 0) = Math::realFromVectorElement( - leftVector, mVirtualNodes[0]->matrixNodeIndex(PhaseType::B)); - (**mSourceValue)(2, 0) = Math::realFromVectorElement( - leftVector, mVirtualNodes[0]->matrixNodeIndex(PhaseType::C)); +void EMT::Ph3::VSIVoltageControlVBR::mnaCompUpdateVoltage( + const Matrix &leftVector) { + (**mIntfVoltage)(0, 0) = + Math::realFromVectorElement(leftVector, matrixNodeIndex(0, 0)); + (**mIntfVoltage)(1, 0) = + Math::realFromVectorElement(leftVector, matrixNodeIndex(0, 1)); + (**mIntfVoltage)(2, 0) = + Math::realFromVectorElement(leftVector, matrixNodeIndex(0, 2)); + (**mSourceValue)(0, 0) = Math::realFromVectorElement( + leftVector, mVirtualNodes[0]->matrixNodeIndex(PhaseType::A)); + (**mSourceValue)(1, 0) = Math::realFromVectorElement( + leftVector, mVirtualNodes[0]->matrixNodeIndex(PhaseType::B)); + (**mSourceValue)(2, 0) = Math::realFromVectorElement( + leftVector, mVirtualNodes[0]->matrixNodeIndex(PhaseType::C)); } -void EMT::Ph3::VSIVoltageControlVBR::mnaCompUpdateCurrent(const Matrix &leftvector) -{ - mFilterCurrent = mEquivCond * ((**mSourceValue).real() - **mIntfVoltage) + mEquivCurrent; - **mIntfCurrent = mSubCapacitorF->mIntfCurrent->get() + mFilterCurrent; +void EMT::Ph3::VSIVoltageControlVBR::mnaCompUpdateCurrent( + const Matrix &leftvector) { + mFilterCurrent = + mEquivCond * ((**mSourceValue).real() - **mIntfVoltage) + mEquivCurrent; + **mIntfCurrent = mSubCapacitorF->mIntfCurrent->get() + mFilterCurrent; } -void EMT::Ph3::VSIVoltageControlVBR::updatePower() -{ - Complex intfVoltageDQ = parkTransformPowerInvariant(**mThetaInv, **mIntfVoltage); - Complex intfCurrentDQ = parkTransformPowerInvariant(**mThetaInv, **mIntfCurrent); - **mPower = Complex(intfVoltageDQ.real() * intfCurrentDQ.real() + intfVoltageDQ.imag() * intfCurrentDQ.imag(), - intfVoltageDQ.imag() * intfCurrentDQ.real() - intfVoltageDQ.real() * intfCurrentDQ.imag()); +void EMT::Ph3::VSIVoltageControlVBR::updatePower() { + Complex intfVoltageDQ = + parkTransformPowerInvariant(**mThetaInv, **mIntfVoltage); + Complex intfCurrentDQ = + parkTransformPowerInvariant(**mThetaInv, **mIntfCurrent); + **mPower = Complex(intfVoltageDQ.real() * intfCurrentDQ.real() + + intfVoltageDQ.imag() * intfCurrentDQ.imag(), + intfVoltageDQ.imag() * intfCurrentDQ.real() - + intfVoltageDQ.real() * intfCurrentDQ.imag()); } -Complex EMT::Ph3::VSIVoltageControlVBR::parkTransformPowerInvariant(Real theta, const Matrix &fabc) -{ - // Calculates fdq = Tdq * fabc - // Assumes that d-axis starts aligned with phase a - Matrix Tdq = getParkTransformMatrixPowerInvariant(theta); - Matrix dqvector = Tdq * fabc; +Complex EMT::Ph3::VSIVoltageControlVBR::parkTransformPowerInvariant( + Real theta, const Matrix &fabc) { + // Calculates fdq = Tdq * fabc + // Assumes that d-axis starts aligned with phase a + Matrix Tdq = getParkTransformMatrixPowerInvariant(theta); + Matrix dqvector = Tdq * fabc; - return Complex(dqvector(0, 0), dqvector(1, 0)); + return Complex(dqvector(0, 0), dqvector(1, 0)); } -Matrix EMT::Ph3::VSIVoltageControlVBR::getParkTransformMatrixPowerInvariant(Real theta) -{ - // Return park matrix for theta - // Assumes that d-axis starts aligned with phase a - Matrix Tdq = Matrix::Zero(2, 3); - Real k = sqrt(2. / 3.); - Tdq << k * cos(theta), k * cos(theta - 2. * M_PI / 3.), k * cos(theta + 2. * M_PI / 3.), - -k * sin(theta), -k * sin(theta - 2. * M_PI / 3.), -k * sin(theta + 2. * M_PI / 3.); - return Tdq; +Matrix EMT::Ph3::VSIVoltageControlVBR::getParkTransformMatrixPowerInvariant( + Real theta) { + // Return park matrix for theta + // Assumes that d-axis starts aligned with phase a + Matrix Tdq = Matrix::Zero(2, 3); + Real k = sqrt(2. / 3.); + Tdq << k * cos(theta), k * cos(theta - 2. * M_PI / 3.), + k * cos(theta + 2. * M_PI / 3.), -k * sin(theta), + -k * sin(theta - 2. * M_PI / 3.), -k * sin(theta + 2. * M_PI / 3.); + return Tdq; } -Matrix EMT::Ph3::VSIVoltageControlVBR::inverseParkTransformPowerInvariant(Real theta, const Complex &fdq) -{ - // Calculates fabc = Tabc * fdq - // with d-axis starts aligned with phase a - Matrix Fdq = Matrix::Zero(2, 1); - Fdq << fdq.real(), fdq.imag(); - Matrix Tabc = getInverseParkTransformMatrixPowerInvariant(theta); +Matrix EMT::Ph3::VSIVoltageControlVBR::inverseParkTransformPowerInvariant( + Real theta, const Complex &fdq) { + // Calculates fabc = Tabc * fdq + // with d-axis starts aligned with phase a + Matrix Fdq = Matrix::Zero(2, 1); + Fdq << fdq.real(), fdq.imag(); + Matrix Tabc = getInverseParkTransformMatrixPowerInvariant(theta); - return Tabc * Fdq; + return Tabc * Fdq; } -Matrix EMT::Ph3::VSIVoltageControlVBR::getInverseParkTransformMatrixPowerInvariant(Real theta) -{ - // Return inverse park matrix for theta - /// with d-axis starts aligned with phase a - Matrix Tabc = Matrix::Zero(3, 2); - Real k = sqrt(2. / 3.); - Tabc << k * cos(theta), -k * sin(theta), - k * cos(theta - 2. * M_PI / 3.), -k * sin(theta - 2. * M_PI / 3.), - k * cos(theta + 2. * M_PI / 3.), -k * sin(theta + 2. * M_PI / 3.); - return Tabc; +Matrix +EMT::Ph3::VSIVoltageControlVBR::getInverseParkTransformMatrixPowerInvariant( + Real theta) { + // Return inverse park matrix for theta + /// with d-axis starts aligned with phase a + Matrix Tabc = Matrix::Zero(3, 2); + Real k = sqrt(2. / 3.); + Tabc << k * cos(theta), -k * sin(theta), k * cos(theta - 2. * M_PI / 3.), + -k * sin(theta - 2. * M_PI / 3.), k * cos(theta + 2. * M_PI / 3.), + -k * sin(theta + 2. * M_PI / 3.); + return Tabc; } diff --git a/dpsim-models/src/SP/SP_Ph1_Capacitor.cpp b/dpsim-models/src/SP/SP_Ph1_Capacitor.cpp index f6739394bd..0f48c24bd4 100644 --- a/dpsim-models/src/SP/SP_Ph1_Capacitor.cpp +++ b/dpsim-models/src/SP/SP_Ph1_Capacitor.cpp @@ -122,7 +122,7 @@ void SP::Ph1::Capacitor::mnaCompPostStep(Real time, Int timeStepCount, void SP::Ph1::Capacitor::mnaCompUpdateVoltage(const Matrix &leftVector) { // v1 - v0 - **mIntfVoltage = Matrix::Zero(3, 1); + **mIntfVoltage = Matrix::Zero(1, 1); if (terminalNotGrounded(1)) { (**mIntfVoltage)(0, 0) = Math::complexFromVectorElement(leftVector, matrixNodeIndex(1)); diff --git a/dpsim-models/src/SP/SP_Ph1_VSIVoltageControlVBR.cpp b/dpsim-models/src/SP/SP_Ph1_VSIVoltageControlVBR.cpp index 06b44237ed..86057251c3 100644 --- a/dpsim-models/src/SP/SP_Ph1_VSIVoltageControlVBR.cpp +++ b/dpsim-models/src/SP/SP_Ph1_VSIVoltageControlVBR.cpp @@ -12,270 +12,251 @@ using namespace CPS; SP::Ph1::VSIVoltageControlVBR::VSIVoltageControlVBR(String uid, String name, - Logger::Level logLevel, - Bool modelAsCurrentSource) - : CompositePowerComp(uid, name, true, true, logLevel), - VSIVoltageSourceInverterDQ(this->mSLog, mAttributes, - modelAsCurrentSource, - false) -{ - - setTerminalNumber(1); - setVirtualNodeNumber(this->determineNumberOfVirtualNodes()); - - **mIntfVoltage = MatrixComp::Zero(1, 1); - **mIntfCurrent = MatrixComp::Zero(1, 1); - mDqToComplexA = Matrix::Zero(2, 2); + Logger::Level logLevel, + Bool modelAsCurrentSource) + : CompositePowerComp(uid, name, true, true, logLevel), + VSIVoltageSourceInverterDQ(this->mSLog, mAttributes, + modelAsCurrentSource, false) { + + setTerminalNumber(1); + setVirtualNodeNumber(this->determineNumberOfVirtualNodes()); + + **mIntfVoltage = MatrixComp::Zero(1, 1); + **mIntfCurrent = MatrixComp::Zero(1, 1); + mDqToComplexA = Matrix::Zero(2, 2); } -void SP::Ph1::VSIVoltageControlVBR::createSubComponents() -{ - // Capacitor as part of the LC filter - mSubCapacitorF = SP::Ph1::Capacitor::make(**mName + "_CapF", mLogLevel); - mSubCapacitorF->setParameters(mCf); - addMNASubComponent(mSubCapacitorF, MNA_SUBCOMP_TASK_ORDER::TASK_BEFORE_PARENT, - MNA_SUBCOMP_TASK_ORDER::TASK_BEFORE_PARENT, false); +void SP::Ph1::VSIVoltageControlVBR::createSubComponents() { + // Capacitor as part of the LC filter + mSubCapacitorF = SP::Ph1::Capacitor::make(**mName + "_CapF", mLogLevel); + mSubCapacitorF->setParameters(mCf); + addMNASubComponent(mSubCapacitorF, MNA_SUBCOMP_TASK_ORDER::TASK_BEFORE_PARENT, + MNA_SUBCOMP_TASK_ORDER::TASK_BEFORE_PARENT, false); } void SP::Ph1::VSIVoltageControlVBR::initializeFromNodesAndTerminals( - Real frequency) -{ - // terminal powers in consumer system -> convert to generator system - **mPower = -terminal(0)->singlePower(); - - // set initial interface quantities --> Current flowing into the inverter is positive - (**mIntfVoltage)(0, 0) = initialSingleVoltage(0); - (**mIntfCurrent)(0, 0) = std::conj(**mPower / (**mIntfVoltage)(0, 0)); - - // initialize filter variables and set initial voltage of virtual nodes - initializeFilterVariables((**mIntfVoltage)(0, 0), (**mIntfCurrent)(0, 0), - mVirtualNodes); - - // calculate initial source value - (**mSourceValue)(0, 0) = - Math::rotatingFrame2to1(**mSourceValue_dq, **mThetaSys, **mThetaInv); - - // Connect & Initialize electrical subcomponents - mSubCapacitorF->connect({SimNode::GND, mTerminals[0]->node()}); - for (auto subcomp : mSubComponents) - { - subcomp->initialize(mFrequencies); - subcomp->initializeFromNodesAndTerminals(frequency); - } - - // TODO: droop - **mOmega = mOmegaNom; - - // - mAdmitance = Complex(1, 0) / Complex(mRf, **mOmega * mLf); - - /// initialize filter current in dp domain - mFilterCurrent = Matrix::Zero(1, 1); - mFilterCurrent(0, 0) = Math::rotatingFrame2to1(**mIfilter_dq, **mThetaSys, **mThetaInv); - - SPDLOG_LOGGER_INFO(mSLog, - "\n--- Initialization from powerflow ---" - "\nTerminal 0 connected to {} = sim node {}" - "\nInverter terminal voltage: {}[V]" - "\nInverter output current: {}[A]", - mTerminals[0]->node()->name(), - mTerminals[0]->node()->matrixNodeIndex(), - Logger::phasorToString((**mIntfVoltage)(0, 0)), - Logger::phasorToString((**mIntfCurrent)(0, 0))); - mSLog->flush(); + Real frequency) { + // terminal powers in consumer system -> convert to generator system + **mPower = -terminal(0)->singlePower(); + + // set initial interface quantities --> Current flowing into the inverter is positive + (**mIntfVoltage)(0, 0) = initialSingleVoltage(0); + (**mIntfCurrent)(0, 0) = std::conj(**mPower / (**mIntfVoltage)(0, 0)); + + // initialize filter variables and set initial voltage of virtual nodes + initializeFilterVariables((**mIntfVoltage)(0, 0), (**mIntfCurrent)(0, 0), + mVirtualNodes); + + // calculate initial source value + (**mSourceValue)(0, 0) = + Math::rotatingFrame2to1(**mSourceValue_dq, **mThetaSys, **mThetaInv); + + // Connect & Initialize electrical subcomponents + mSubCapacitorF->connect({SimNode::GND, mTerminals[0]->node()}); + for (auto subcomp : mSubComponents) { + subcomp->initialize(mFrequencies); + subcomp->initializeFromNodesAndTerminals(frequency); + } + + // TODO: droop + **mOmega = mOmegaNom; + + // + mAdmitance = Complex(1, 0) / Complex(mRf, **mOmega * mLf); + + /// initialize filter current in dp domain + mFilterCurrent = Matrix::Zero(1, 1); + mFilterCurrent(0, 0) = + Math::rotatingFrame2to1(**mIfilter_dq, **mThetaSys, **mThetaInv); + + SPDLOG_LOGGER_INFO(mSLog, + "\n--- Initialization from powerflow ---" + "\nTerminal 0 connected to {} = sim node {}" + "\nInverter terminal voltage: {}[V]" + "\nInverter output current: {}[A]", + mTerminals[0]->node()->name(), + mTerminals[0]->node()->matrixNodeIndex(), + Logger::phasorToString((**mIntfVoltage)(0, 0)), + Logger::phasorToString((**mIntfCurrent)(0, 0))); + mSLog->flush(); } void SP::Ph1::VSIVoltageControlVBR::mnaParentInitialize( - Real omega, Real timeStep, Attribute::Ptr leftVector) -{ - mTimeStep = timeStep; - if (mWithControl) - { - mVSIController->initialize(**mSourceValue_dq, **mVcap_dq, **mIfilter_dq, - mTimeStep, false); - mVSIController->calculateVBRconstants(); - } - - // initialize auxialiar Matrix - mA = Matrix::Zero(2, 2); - mE = Matrix::Zero(2, 2); - mY = Matrix::Zero(2, 2); - - mY << mRf, -mLf * mOmegaNom, mLf * mOmegaNom, mRf; - mY = mY.inverse(); - - mA << std::dynamic_pointer_cast(mVSIController) - ->mA_VBR, - 0, 0, - std::dynamic_pointer_cast(mVSIController) - ->mA_VBR; - - mE << std::dynamic_pointer_cast(mVSIController) - ->mE_VBR, - 0, 0, - std::dynamic_pointer_cast(mVSIController) - ->mE_VBR; - - // get matrix dimension to properly set variable entries - auto n = leftVector->asRawPointer()->rows(); - auto complexOffset = (UInt)(n / 2); - // upper left - mVariableSystemMatrixEntries.push_back( - std::make_pair(mVirtualNodes[0]->matrixNodeIndex(), - mVirtualNodes[0]->matrixNodeIndex())); - mVariableSystemMatrixEntries.push_back(std::make_pair( - mVirtualNodes[0]->matrixNodeIndex() + complexOffset, - mVirtualNodes[0]->matrixNodeIndex())); - mVariableSystemMatrixEntries.push_back(std::make_pair( - mVirtualNodes[0]->matrixNodeIndex(), - mVirtualNodes[0]->matrixNodeIndex() + complexOffset)); - mVariableSystemMatrixEntries.push_back(std::make_pair( - mVirtualNodes[0]->matrixNodeIndex() + complexOffset, - mVirtualNodes[0]->matrixNodeIndex() + complexOffset)); - - // off diagonal - mVariableSystemMatrixEntries.push_back(std::make_pair( - mVirtualNodes[0]->matrixNodeIndex(), matrixNodeIndex(0, 0))); - mVariableSystemMatrixEntries.push_back(std::make_pair( - mVirtualNodes[0]->matrixNodeIndex() + complexOffset, - matrixNodeIndex(0, 0))); - mVariableSystemMatrixEntries.push_back( - std::make_pair(mVirtualNodes[0]->matrixNodeIndex(), - matrixNodeIndex(0, 0) + complexOffset)); - mVariableSystemMatrixEntries.push_back(std::make_pair( - mVirtualNodes[0]->matrixNodeIndex() + complexOffset, - matrixNodeIndex(0, 0) + complexOffset)); - - SPDLOG_LOGGER_INFO(mSLog, "List of index pairs of varying matrix entries: "); - for (auto indexPair : mVariableSystemMatrixEntries) - SPDLOG_LOGGER_INFO(mSLog, "({}, {})", indexPair.first, indexPair.second); - mSLog->flush(); + Real omega, Real timeStep, Attribute::Ptr leftVector) { + mTimeStep = timeStep; + if (mWithControl) { + mVSIController->initialize(**mSourceValue_dq, **mVcap_dq, **mIfilter_dq, + mTimeStep, false); + mVSIController->calculateVBRconstants(); + } + + // initialize auxialiar Matrix + mA = Matrix::Zero(2, 2); + mE = Matrix::Zero(2, 2); + mY = Matrix::Zero(2, 2); + + mY << mRf, -mLf * mOmegaNom, mLf * mOmegaNom, mRf; + mY = mY.inverse(); + + mA << std::dynamic_pointer_cast(mVSIController) + ->mA_VBR, + 0, 0, + std::dynamic_pointer_cast(mVSIController) + ->mA_VBR; + + mE << std::dynamic_pointer_cast(mVSIController) + ->mE_VBR, + 0, 0, + std::dynamic_pointer_cast(mVSIController) + ->mE_VBR; + + // get matrix dimension to properly set variable entries + auto n = leftVector->asRawPointer()->rows(); + auto complexOffset = (UInt)(n / 2); + // upper left + mVariableSystemMatrixEntries.push_back( + std::make_pair(mVirtualNodes[0]->matrixNodeIndex(), + mVirtualNodes[0]->matrixNodeIndex())); + mVariableSystemMatrixEntries.push_back(std::make_pair( + mVirtualNodes[0]->matrixNodeIndex() + complexOffset, + mVirtualNodes[0]->matrixNodeIndex())); + mVariableSystemMatrixEntries.push_back(std::make_pair( + mVirtualNodes[0]->matrixNodeIndex(), + mVirtualNodes[0]->matrixNodeIndex() + complexOffset)); + mVariableSystemMatrixEntries.push_back(std::make_pair( + mVirtualNodes[0]->matrixNodeIndex() + complexOffset, + mVirtualNodes[0]->matrixNodeIndex() + complexOffset)); + + // off diagonal + mVariableSystemMatrixEntries.push_back(std::make_pair( + mVirtualNodes[0]->matrixNodeIndex(), matrixNodeIndex(0, 0))); + mVariableSystemMatrixEntries.push_back(std::make_pair( + mVirtualNodes[0]->matrixNodeIndex() + complexOffset, + matrixNodeIndex(0, 0))); + mVariableSystemMatrixEntries.push_back( + std::make_pair(mVirtualNodes[0]->matrixNodeIndex(), + matrixNodeIndex(0, 0) + complexOffset)); + mVariableSystemMatrixEntries.push_back(std::make_pair( + mVirtualNodes[0]->matrixNodeIndex() + complexOffset, + matrixNodeIndex(0, 0) + complexOffset)); + + SPDLOG_LOGGER_INFO(mSLog, "List of index pairs of varying matrix entries: "); + for (auto indexPair : mVariableSystemMatrixEntries) + SPDLOG_LOGGER_INFO(mSLog, "({}, {})", indexPair.first, indexPair.second); + mSLog->flush(); } -void SP::Ph1::VSIVoltageControlVBR::calculateResistanceMatrix() -{ - mEMatrix = mDqToComplexA * mE * mComplexAToDq; - mAMatrix = mDqToComplexA * mA * mComplexAToDq * mY; +void SP::Ph1::VSIVoltageControlVBR::calculateResistanceMatrix() { + mEMatrix = mDqToComplexA * mE * mComplexAToDq; + mAMatrix = mDqToComplexA * mA * mComplexAToDq * mY; } void SP::Ph1::VSIVoltageControlVBR::mnaParentApplySystemMatrixStamp( - SparseMatrixRow &systemMatrix) -{ - update_DqToComplexATransformMatrix(); - mComplexAToDq = mDqToComplexA.transpose(); - calculateResistanceMatrix(); - - // Stamp filter current: I = (Vc - Vsource) * Admittance - Math::addToMatrixElement(systemMatrix, matrixNodeIndex(0), matrixNodeIndex(0), - mAdmitance); - Math::addToMatrixElement(systemMatrix, matrixNodeIndex(0), - mVirtualNodes[0]->matrixNodeIndex(), -mAdmitance); - - // Stamp voltage source equation - Math::addToMatrixElement(systemMatrix, mVirtualNodes[0]->matrixNodeIndex(), - mVirtualNodes[0]->matrixNodeIndex(), - Matrix::Identity(2, 2) + mAMatrix); - Math::addToMatrixElement(systemMatrix, mVirtualNodes[0]->matrixNodeIndex(), - matrixNodeIndex(0), - -mEMatrix - mAMatrix); + SparseMatrixRow &systemMatrix) { + update_DqToComplexATransformMatrix(); + mComplexAToDq = mDqToComplexA.transpose(); + calculateResistanceMatrix(); + + // Stamp filter current: I = (Vc - Vsource) * Admittance + Math::addToMatrixElement(systemMatrix, matrixNodeIndex(0), matrixNodeIndex(0), + mAdmitance); + Math::addToMatrixElement(systemMatrix, matrixNodeIndex(0), + mVirtualNodes[0]->matrixNodeIndex(), -mAdmitance); + + // Stamp voltage source equation + Math::addToMatrixElement(systemMatrix, mVirtualNodes[0]->matrixNodeIndex(), + mVirtualNodes[0]->matrixNodeIndex(), + Matrix::Identity(2, 2) + mAMatrix); + Math::addToMatrixElement(systemMatrix, mVirtualNodes[0]->matrixNodeIndex(), + matrixNodeIndex(0), -mEMatrix - mAMatrix); } void SP::Ph1::VSIVoltageControlVBR::mnaParentAddPreStepDependencies( - AttributeBase::List &prevStepDependencies, - AttributeBase::List &attributeDependencies, - AttributeBase::List &modifiedAttributes) -{ - modifiedAttributes.push_back(mRightVector); - prevStepDependencies.push_back(mIntfVoltage); - prevStepDependencies.push_back(mIntfCurrent); + AttributeBase::List &prevStepDependencies, + AttributeBase::List &attributeDependencies, + AttributeBase::List &modifiedAttributes) { + modifiedAttributes.push_back(mRightVector); + prevStepDependencies.push_back(mIntfVoltage); + prevStepDependencies.push_back(mIntfCurrent); } void SP::Ph1::VSIVoltageControlVBR::mnaParentPreStep(Real time, - Int timeStepCount) -{ - // get measurements - **mVcap_dq = Math::rotatingFrame2to1((**mSubCapacitorF->mIntfVoltage)(0, 0), - **mThetaInv, **mThetaSys); - **mIfilter_dq = - Math::rotatingFrame2to1(mFilterCurrent(0, 0), **mThetaInv, **mThetaSys); + Int timeStepCount) { + // get measurements + **mVcap_dq = Math::rotatingFrame2to1((**mSubCapacitorF->mIntfVoltage)(0, 0), + **mThetaInv, **mThetaSys); + **mIfilter_dq = + Math::rotatingFrame2to1(mFilterCurrent(0, 0), **mThetaInv, **mThetaSys); - // TODO: droop - // if (mWithDroop) - // mDroop->signalStep(time, timeStepCount); + // TODO: droop + // if (mWithDroop) + // mDroop->signalStep(time, timeStepCount); - // VCO Step - **mThetaInv = **mThetaInv + mTimeStep * **mOmega; + // VCO Step + **mThetaInv = **mThetaInv + mTimeStep * **mOmega; - // Update nominal system angle - **mThetaSys = **mThetaSys + mTimeStep * mOmegaNom; + // Update nominal system angle + **mThetaSys = **mThetaSys + mTimeStep * mOmegaNom; - // calculat history term - if (mWithControl) - mVhist = mVSIController->stepVBR(**mVcap_dq, **mIfilter_dq); + // calculat history term + if (mWithControl) + mVhist = mVSIController->stepVBR(**mVcap_dq, **mIfilter_dq); - update_DqToComplexATransformMatrix(); - mComplexAToDq = mDqToComplexA.transpose(); + update_DqToComplexATransformMatrix(); + mComplexAToDq = mDqToComplexA.transpose(); - // calculate resistance matrix at t=k+1 - calculateResistanceMatrix(); + // calculate resistance matrix at t=k+1 + calculateResistanceMatrix(); - mnaApplyRightSideVectorStamp(**mRightVector); + mnaApplyRightSideVectorStamp(**mRightVector); } void SP::Ph1::VSIVoltageControlVBR::mnaParentApplyRightSideVectorStamp( - Matrix &rightVector) -{ - Complex dqToComplexA = Complex(mDqToComplexA(0, 0), -mDqToComplexA(0, 1)); - Complex eqVoltageSource = dqToComplexA * mVhist; - - Math::setVectorElement(rightVector, mVirtualNodes[0]->matrixNodeIndex(), - eqVoltageSource); + Matrix &rightVector) { + Complex dqToComplexA = Complex(mDqToComplexA(0, 0), -mDqToComplexA(0, 1)); + Complex eqVoltageSource = dqToComplexA * mVhist; + Math::setVectorElement(rightVector, mVirtualNodes[0]->matrixNodeIndex(), + eqVoltageSource); } void SP::Ph1::VSIVoltageControlVBR::mnaParentAddPostStepDependencies( - AttributeBase::List &prevStepDependencies, - AttributeBase::List &attributeDependencies, - AttributeBase::List &modifiedAttributes, - Attribute::Ptr &leftVector) -{ - attributeDependencies.push_back(leftVector); - attributeDependencies.push_back(mRightVector); - modifiedAttributes.push_back(mIntfVoltage); - modifiedAttributes.push_back(mIntfCurrent); + AttributeBase::List &prevStepDependencies, + AttributeBase::List &attributeDependencies, + AttributeBase::List &modifiedAttributes, + Attribute::Ptr &leftVector) { + attributeDependencies.push_back(leftVector); + attributeDependencies.push_back(mRightVector); + modifiedAttributes.push_back(mIntfVoltage); + modifiedAttributes.push_back(mIntfCurrent); } void SP::Ph1::VSIVoltageControlVBR::mnaParentPostStep( - Real time, Int timeStepCount, Attribute::Ptr &leftVector) -{ - mnaCompUpdateVoltage(**leftVector); - mnaCompUpdateCurrent(**leftVector); - updatePower(); + Real time, Int timeStepCount, Attribute::Ptr &leftVector) { + mnaCompUpdateVoltage(**leftVector); + mnaCompUpdateCurrent(**leftVector); + updatePower(); } void SP::Ph1::VSIVoltageControlVBR::mnaCompUpdateCurrent( - const Matrix &leftvector) -{ - mFilterCurrent = (**mSourceValue - **mIntfVoltage) * mAdmitance; - **mIntfCurrent = - mSubCapacitorF->mIntfCurrent->get() + mFilterCurrent; + const Matrix &leftvector) { + mFilterCurrent = (**mSourceValue - **mIntfVoltage) * mAdmitance; + **mIntfCurrent = mSubCapacitorF->mIntfCurrent->get() + mFilterCurrent; } void SP::Ph1::VSIVoltageControlVBR::mnaCompUpdateVoltage( - const Matrix &leftVector) -{ - (**mIntfVoltage)(0, 0) = - Math::complexFromVectorElement(leftVector, matrixNodeIndex(0)); - (**mSourceValue)(0, 0) = Math::complexFromVectorElement( - leftVector, mVirtualNodes[0]->matrixNodeIndex()); + const Matrix &leftVector) { + (**mIntfVoltage)(0, 0) = + Math::complexFromVectorElement(leftVector, matrixNodeIndex(0)); + (**mSourceValue)(0, 0) = Math::complexFromVectorElement( + leftVector, mVirtualNodes[0]->matrixNodeIndex()); } -void SP::Ph1::VSIVoltageControlVBR::updatePower() -{ - **mPower = (**mIntfVoltage)(0, 0) * std::conj((**mIntfCurrent)(0, 0)); +void SP::Ph1::VSIVoltageControlVBR::updatePower() { + **mPower = (**mIntfVoltage)(0, 0) * std::conj((**mIntfCurrent)(0, 0)); } -void SP::Ph1::VSIVoltageControlVBR::update_DqToComplexATransformMatrix() -{ - mDqToComplexA << cos(**mThetaInv - **mThetaSys), -sin(**mThetaInv - **mThetaSys), - sin(**mThetaInv - **mThetaSys), cos(**mThetaInv - **mThetaSys); +void SP::Ph1::VSIVoltageControlVBR::update_DqToComplexATransformMatrix() { + mDqToComplexA << cos(**mThetaInv - **mThetaSys), + -sin(**mThetaInv - **mThetaSys), sin(**mThetaInv - **mThetaSys), + cos(**mThetaInv - **mThetaSys); } diff --git a/dpsim-models/src/Signal/VSIControlType1.cpp b/dpsim-models/src/Signal/VSIControlType1.cpp index adb97e9b38..3cbda3e54f 100644 --- a/dpsim-models/src/Signal/VSIControlType1.cpp +++ b/dpsim-models/src/Signal/VSIControlType1.cpp @@ -1,227 +1,261 @@ -#include #include +#include using namespace CPS; using namespace CPS::Signal; -void VSIControlType1::setParameters(std::shared_ptr parameters) -{ - if (auto params = std::dynamic_pointer_cast(parameters)) - { - mParameters = params; - SPDLOG_LOGGER_INFO(mSLog, - "\nVSIController Type1 Parameters: " - "\nKpv: {:e}" - "\nKiv: {:e}" - "\nKpc: {:e}" - "\nKic: {:e}", - mParameters->Kpv, mParameters->Kiv, - mParameters->Kpc, mParameters->Kic); - mSLog->flush(); - } - else - { - std::cout << "Type of parameters class of " << this->name() << " has to be VSIControlType1!" << std::endl; - throw CPS::TypeException(); - } +void VSIControlType1::setParameters( + std::shared_ptr parameters) { + if (auto params = + std::dynamic_pointer_cast( + parameters)) { + mParameters = params; + SPDLOG_LOGGER_INFO(mSLog, + "\nVSIController Type1 Parameters: " + "\nKpv: {:e}" + "\nKiv: {:e}" + "\nKpc: {:e}" + "\nKic: {:e}", + mParameters->Kpv, mParameters->Kiv, mParameters->Kpc, + mParameters->Kic); + mSLog->flush(); + } else { + std::cout << "Type of parameters class of " << this->name() + << " has to be VSIControlType1!" << std::endl; + throw CPS::TypeException(); + } } -void VSIControlType1::initialize(const Complex &Vsref_dq, const Complex &Vcap_dq, - const Complex &Ifilter_dq, Real time_step, Bool modelAsCurrentSource) -{ - mTimeStep = time_step; - mModelAsCurrentSource = modelAsCurrentSource; - - mVcap_dq = Vcap_dq; - mIfilter_dq = Ifilter_dq; - - **mPhi_d = Ifilter_dq.real() / mParameters->Kiv; - **mPhi_q = Ifilter_dq.imag() / mParameters->Kiv; - **mGamma_d = (Vsref_dq.real() - Vcap_dq.real()) / mParameters->Kic; - **mGamma_q = (Vsref_dq.imag() - Vcap_dq.imag()) / mParameters->Kic; - - SPDLOG_LOGGER_INFO(mSLog, - "\nInitialize controller states:" - "\n\tPhi_d = {}" - "\n\tPhi_q = {}" - "\n\tGamma_d = {}" - "\n\tGamma_q = {}", - **mPhi_d, **mPhi_q, **mGamma_d, **mGamma_q); - mSLog->flush(); - - **mStateCurr << **mPhi_d, **mPhi_q, **mGamma_d, **mGamma_q; - - // initialize state matrix - // [x] = [phid, phiq, gammad, gammaq] - // [u] = [vdref, vqref, vdc, vqc, idc, idq] - // [y] = [vdout, vqout] - - mA << 0, 0, 0, 0, - 0, 0, 0, 0, - mParameters->Kiv, 0, 0, 0, - 0, mParameters->Kiv, 0, 0; - - mB << 1, 0, -1, 0, 0, 0, - 0, 1, 0, -1, 0, 0, - mParameters->Kpv, 0, -mParameters->Kpv, 0, -1, 0, - 0, mParameters->Kpv, 0, -mParameters->Kpv, 0, -1; - - mC << mParameters->Kpc * mParameters->Kiv, 0, mParameters->Kic, 0, - 0, mParameters->Kpc * mParameters->Kiv, 0, mParameters->Kic; - - mD << mParameters->Kpc * mParameters->Kpv, 0, -mParameters->Kpc * mParameters->Kpv + 1, 0, -mParameters->Kpc, 0, - 0, mParameters->Kpc * mParameters->Kpv, 0, -mParameters->Kpc * mParameters->Kpv + 1, 0, -mParameters->Kpc; - - Math::calculateStateSpaceTrapezoidalMatrices(mA, mB, Matrix::Zero(4, 1), mTimeStep, mATrapezoidal, mBTrapezoidal, mCTrapezoidal); - - **mInputCurr << mParameters->VdRef, mParameters->VqRef, Vcap_dq.real(), Vcap_dq.imag(), Ifilter_dq.real(), Ifilter_dq.imag(); - - // Log state-space matrices - SPDLOG_LOGGER_INFO(mSLog, - "\nState space matrices:" - "\nA = \n{}" - "\nB = \n{}" - "\nC = \n{}" - "\nD = \n{}", - mA, mB, mC, mD); - mSLog->flush(); +void VSIControlType1::initialize(const Complex &Vsref_dq, + const Complex &Vcap_dq, + const Complex &Ifilter_dq, Real time_step, + Bool modelAsCurrentSource) { + mTimeStep = time_step; + mModelAsCurrentSource = modelAsCurrentSource; + + mVcap_dq = Vcap_dq; + mIfilter_dq = Ifilter_dq; + **mPhi_d = Ifilter_dq.real() / mParameters->Kiv; + **mPhi_q = Ifilter_dq.imag() / mParameters->Kiv; + **mGamma_d = (Vsref_dq.real() - Vcap_dq.real()) / mParameters->Kic; + **mGamma_q = (Vsref_dq.imag() - Vcap_dq.imag()) / mParameters->Kic; + + SPDLOG_LOGGER_INFO(mSLog, + "\nInitialize controller states:" + "\n\tPhi_d = {}" + "\n\tPhi_q = {}" + "\n\tGamma_d = {}" + "\n\tGamma_q = {}", + **mPhi_d, **mPhi_q, **mGamma_d, **mGamma_q); + mSLog->flush(); + + **mStateCurr << **mPhi_d, **mPhi_q, **mGamma_d, **mGamma_q; + + // initialize state matrix + // [x] = [phid, phiq, gammad, gammaq] + // [u] = [vdref, vqref, vdc, vqc, idc, idq] + // [y] = [vdout, vqout] + + mA << 0, 0, 0, 0, 0, 0, 0, 0, mParameters->Kiv, 0, 0, 0, 0, mParameters->Kiv, + 0, 0; + + mB << 1, 0, -1, 0, 0, 0, 0, 1, 0, -1, 0, 0, mParameters->Kpv, 0, + -mParameters->Kpv, 0, -1, 0, 0, mParameters->Kpv, 0, -mParameters->Kpv, 0, + -1; + + mC << mParameters->Kpc * mParameters->Kiv, 0, mParameters->Kic, 0, 0, + mParameters->Kpc * mParameters->Kiv, 0, mParameters->Kic; + + mD << mParameters->Kpc * mParameters->Kpv, 0, + -mParameters->Kpc * mParameters->Kpv + 1, 0, -mParameters->Kpc, 0, 0, + mParameters->Kpc * mParameters->Kpv, 0, + -mParameters->Kpc * mParameters->Kpv + 1, 0, -mParameters->Kpc; + + Math::calculateStateSpaceTrapezoidalMatrices(mA, mB, Matrix::Zero(4, 1), + mTimeStep, mATrapezoidal, + mBTrapezoidal, mCTrapezoidal); + + **mInputCurr << mParameters->VdRef, mParameters->VqRef, Vcap_dq.real(), + Vcap_dq.imag(), Ifilter_dq.real(), Ifilter_dq.imag(); + + // Log state-space matrices + SPDLOG_LOGGER_INFO(mSLog, + "\nState space matrices:" + "\nA = \n{}" + "\nB = \n{}" + "\nC = \n{}" + "\nD = \n{}", + mA, mB, mC, mD); + mSLog->flush(); } -void VSIControlType1::calculateVBRconstants() -{ - mA_VBR = mParameters->Kpc + mParameters->Kic * mTimeStep / 2.; - mB_VBR = mParameters->Kpv + mParameters->Kiv * mTimeStep / 2.; - mB2_VBR = mParameters->Kpv + mParameters->Kiv * mTimeStep; - mC_VBR = mParameters->Kiv * mTimeStep / 2.; - mD_VBR = mParameters->Kic * mTimeStep / 2.; - mE_VBR = -mA_VBR * mB_VBR + 1.; - - // initialize Iref - mIref_d = mB2_VBR * mParameters->VdRef - mB_VBR * mVcap_dq.real() - mC_VBR * mVcap_dq.real() + **mPhi_d * mParameters->Kiv; - mIref_q = mB2_VBR * mParameters->VqRef - mB_VBR * mVcap_dq.imag() - mC_VBR * mVcap_dq.imag() + **mPhi_q * mParameters->Kiv; - - // calculate Vhist at t=k+1 - Real Vhist_d = -mA_VBR * mC_VBR * mVcap_dq.real() - mD_VBR * mIfilter_dq.real() + mD_VBR * mIref_d + mA_VBR * **mPhi_d * mParameters->Kiv + **mGamma_d * mParameters->Kic + mA_VBR * mB2_VBR * mParameters->VdRef; - Real Vhist_q = -mA_VBR * mC_VBR * mVcap_dq.imag() - mD_VBR * mIfilter_dq.imag() + mD_VBR * mIref_q + mA_VBR * **mPhi_q * mParameters->Kiv + **mGamma_q * mParameters->Kic + mA_VBR * mB2_VBR * mParameters->VqRef; - - // Test - Real Vsource_test_d = mVcap_dq.real() - mA_VBR * mB_VBR * mVcap_dq.real() - mA_VBR * mIfilter_dq.real() + Vhist_d; - Real Vsource_test_q = mVcap_dq.imag() - mA_VBR * mB_VBR * mVcap_dq.imag() - mA_VBR * mIfilter_dq.imag() + Vhist_q; - - // intial reference current - mIfilter_dq = Complex(mIref_d, mIref_q); - - SPDLOG_LOGGER_INFO(mSLog, - "\n--- VBR constants: ---" - "\nA_VBR = {}" - "\nB_VBR = {}" - "\nC_VBR = {}" - "\nD_VBR = {}" - "\nE_VBR = {}" - "\nIref_d = {}" - "\nIref_d = {}" - "\nInit Vhist_d = {}" - "\nInit Vhist_q = {}" - "\nInit Vsource_test_d = {}" - "\nInit Vsource_test_q = {}", - mA_VBR, mB_VBR, mC_VBR, mD_VBR, mE_VBR, - mIref_d, mIref_q, Vhist_d, Vhist_q, - Vsource_test_d, Vsource_test_q); - mSLog->flush(); +void VSIControlType1::calculateVBRconstants() { + mA_VBR = mParameters->Kpc + mParameters->Kic * mTimeStep / 2.; + mB_VBR = mParameters->Kpv + mParameters->Kiv * mTimeStep / 2.; + mB2_VBR = mParameters->Kpv + mParameters->Kiv * mTimeStep; + mC_VBR = mParameters->Kiv * mTimeStep / 2.; + mD_VBR = mParameters->Kic * mTimeStep / 2.; + mE_VBR = -mA_VBR * mB_VBR + 1.; + + // initialize Iref + mIref_d = mB2_VBR * mParameters->VdRef - mB_VBR * mVcap_dq.real() - + mC_VBR * mVcap_dq.real() + **mPhi_d * mParameters->Kiv; + mIref_q = mB2_VBR * mParameters->VqRef - mB_VBR * mVcap_dq.imag() - + mC_VBR * mVcap_dq.imag() + **mPhi_q * mParameters->Kiv; + + // calculate Vhist at t=k+1 + Real Vhist_d = + -mA_VBR * mC_VBR * mVcap_dq.real() - mD_VBR * mIfilter_dq.real() + + mD_VBR * mIref_d + mA_VBR * **mPhi_d * mParameters->Kiv + + **mGamma_d * mParameters->Kic + mA_VBR * mB2_VBR * mParameters->VdRef; + Real Vhist_q = + -mA_VBR * mC_VBR * mVcap_dq.imag() - mD_VBR * mIfilter_dq.imag() + + mD_VBR * mIref_q + mA_VBR * **mPhi_q * mParameters->Kiv + + **mGamma_q * mParameters->Kic + mA_VBR * mB2_VBR * mParameters->VqRef; + + // Test + Real Vsource_test_d = mVcap_dq.real() - mA_VBR * mB_VBR * mVcap_dq.real() - + mA_VBR * mIfilter_dq.real() + Vhist_d; + Real Vsource_test_q = mVcap_dq.imag() - mA_VBR * mB_VBR * mVcap_dq.imag() - + mA_VBR * mIfilter_dq.imag() + Vhist_q; + + // intial reference current + mIfilter_dq = Complex(mIref_d, mIref_q); + + SPDLOG_LOGGER_INFO(mSLog, + "\n--- VBR constants: ---" + "\nA_VBR = {}" + "\nB_VBR = {}" + "\nC_VBR = {}" + "\nD_VBR = {}" + "\nE_VBR = {}" + "\nIref_d = {}" + "\nIref_d = {}" + "\nInit Vhist_d = {}" + "\nInit Vhist_q = {}" + "\nInit Vsource_test_d = {}" + "\nInit Vsource_test_q = {}", + mA_VBR, mB_VBR, mC_VBR, mD_VBR, mE_VBR, mIref_d, mIref_q, + Vhist_d, Vhist_q, Vsource_test_d, Vsource_test_q); + mSLog->flush(); } -Complex VSIControlType1::step(const Complex &Vcap_dq, const Complex &Ifilter_dq) -{ - // - **mStatePrev = **mStateCurr; - - // get current inputs - **mInputPrev = **mInputCurr; - **mInputCurr << mParameters->VdRef, mParameters->VqRef, Vcap_dq.real(), Vcap_dq.imag(), Ifilter_dq.real(), Ifilter_dq.imag(); - - // calculate new states - //**mStateCurr = Math::StateSpaceTrapezoidal(**mStatePrev, mA, mB, mTimeStep, **mInputCurr, **mInputPrev); - **mStateCurr = Math::applyStateSpaceTrapezoidalMatrices(mATrapezoidal, mBTrapezoidal, mCTrapezoidal, **mStatePrev, **mInputCurr, **mInputPrev); - - // calculate controller outputs - if (mModelAsCurrentSource) - { - Real error_d = mParameters->VdRef - Vcap_dq.real(); - Real error_q = mParameters->VqRef - Vcap_dq.imag(); - (**mOutput)(0, 0) = ((**mStateCurr)(0, 0) * mParameters->Kiv + mParameters->Kpv * error_d) * (mTimeStep / mParameters->tau) + (1. - mTimeStep / mParameters->tau) * Ifilter_dq.real(); - (**mOutput)(1, 0) = ((**mStateCurr)(1, 0) * mParameters->Kiv + mParameters->Kpv * error_q) * (mTimeStep / mParameters->tau) + (1. - mTimeStep / mParameters->tau) * Ifilter_dq.imag(); - } - else - { - **mOutput = mC * **mStateCurr + mD * **mInputCurr; - } - - SPDLOG_LOGGER_DEBUG(mSLog, - "\n - InputCurr = \n{}" - "\n - InputPrev = \n{}" - "\n - StatePrev = \n{}" - "\n - StateCurr = \n{}" - "\n - Output values: \n{}", - **mInputCurr, **mInputPrev, **mStatePrev, **mStateCurr, **mOutput); - - // - return Complex((**mOutput)(0, 0), (**mOutput)(1, 0)); +Complex VSIControlType1::step(const Complex &Vcap_dq, + const Complex &Ifilter_dq) { + // + **mStatePrev = **mStateCurr; + + // get current inputs + **mInputPrev = **mInputCurr; + **mInputCurr << mParameters->VdRef, mParameters->VqRef, Vcap_dq.real(), + Vcap_dq.imag(), Ifilter_dq.real(), Ifilter_dq.imag(); + + // calculate new states + //**mStateCurr = Math::StateSpaceTrapezoidal(**mStatePrev, mA, mB, mTimeStep, **mInputCurr, **mInputPrev); + **mStateCurr = Math::applyStateSpaceTrapezoidalMatrices( + mATrapezoidal, mBTrapezoidal, mCTrapezoidal, **mStatePrev, **mInputCurr, + **mInputPrev); + + // calculate controller outputs + if (mModelAsCurrentSource) { + Real error_d = mParameters->VdRef - Vcap_dq.real(); + Real error_q = mParameters->VqRef - Vcap_dq.imag(); + (**mOutput)(0, 0) = + ((**mStateCurr)(0, 0) * mParameters->Kiv + mParameters->Kpv * error_d) * + (mTimeStep / mParameters->tau) + + (1. - mTimeStep / mParameters->tau) * Ifilter_dq.real(); + (**mOutput)(1, 0) = + ((**mStateCurr)(1, 0) * mParameters->Kiv + mParameters->Kpv * error_q) * + (mTimeStep / mParameters->tau) + + (1. - mTimeStep / mParameters->tau) * Ifilter_dq.imag(); + } else { + **mOutput = mC * **mStateCurr + mD * **mInputCurr; + } + + SPDLOG_LOGGER_DEBUG(mSLog, + "\n - InputCurr = \n{}" + "\n - InputPrev = \n{}" + "\n - StatePrev = \n{}" + "\n - StateCurr = \n{}" + "\n - Output values: \n{}", + **mInputCurr, **mInputPrev, **mStatePrev, **mStateCurr, + **mOutput); + + // + return Complex((**mOutput)(0, 0), (**mOutput)(1, 0)); } -Complex VSIControlType1::stepVBR(const Complex &Vcap_dq, const Complex &Ifilter_dq) -{ - // store previous values (t=k-1) - Complex Vcap_dq_prev = mVcap_dq; - Complex Ifilter_dq_prev = mIfilter_dq; - Real Iref_d_prev = mIref_d; - Real Iref_q_prev = mIref_q; - - // update variables at t=k - mVcap_dq = Vcap_dq; - mIfilter_dq = Ifilter_dq; - - // calculate reference current at time t=k - mIref_d = mB2_VBR * mParameters->VdRef - mB_VBR * mVcap_dq.real() - mC_VBR * Vcap_dq_prev.real() + **mPhi_d * mParameters->Kiv; - mIref_q = mB2_VBR * mParameters->VqRef - mB_VBR * mVcap_dq.imag() - mC_VBR * Vcap_dq_prev.imag() + **mPhi_q * mParameters->Kiv; - - // Update phi at time t=k - **mPhi_d = **mPhi_d + (mTimeStep / 2.) * (2 * mParameters->VdRef - mVcap_dq.real() - Vcap_dq_prev.real()); - **mPhi_q = **mPhi_q + (mTimeStep / 2.) * (2 * mParameters->VqRef - mVcap_dq.imag() - Vcap_dq_prev.imag()); - - // Update lambda at time t=k - **mGamma_d = **mGamma_d + (mTimeStep / 2.) * (mIref_d - mIfilter_dq.real() + Iref_d_prev - Ifilter_dq_prev.real()); - **mGamma_q = **mGamma_q + (mTimeStep / 2.) * (mIref_q - mIfilter_dq.imag() + Iref_q_prev - Ifilter_dq_prev.imag()); - - // calculate Vvbr at t=k+1 - Real Vhist_d = -mA_VBR * mC_VBR * mVcap_dq.real() - mD_VBR * mIfilter_dq.real() + mD_VBR * mIref_d + mA_VBR * **mPhi_d * mParameters->Kiv + **mGamma_d * mParameters->Kic + mA_VBR * mB2_VBR * mParameters->VdRef; - Real Vhist_q = -mA_VBR * mC_VBR * mVcap_dq.imag() - mD_VBR * mIfilter_dq.imag() + mD_VBR * mIref_q + mA_VBR * **mPhi_q * mParameters->Kiv + **mGamma_q * mParameters->Kic + mA_VBR * mB2_VBR * mParameters->VqRef; - - // Test - Real Vsource_test_d = mVcap_dq.real() - mA_VBR * mB_VBR * mVcap_dq.real() - mA_VBR * mIfilter_dq.real() + Vhist_d; - Real Vsource_test_q = mVcap_dq.imag() - mA_VBR * mB_VBR * mVcap_dq.imag() - mA_VBR * mIfilter_dq.imag() + Vhist_q; - // - - SPDLOG_LOGGER_DEBUG(mSLog, - "\n--- Test Step: ---" - "\nmVcap_dq_d = {}" - "\nmVcap_dq_q = {}" - "\nIfilter_dq_d = {}" - "\nIfilter_dq_q = {}" - "\nVsource_test_d = {}" - "\nVsource_test_q = {}" - "\nVhist_d = {}" - "\nVhist_q = {}" - "\nIref_d = {}" - "\nIref_d = {}" - "\nmPhi_d = {}" - "\nmPhi_q = {}" - "\nmGamma_d = {}" - "\nmGamma_q = {}", - mVcap_dq.real(), mVcap_dq.imag(), - Ifilter_dq.real(), Ifilter_dq.imag(), - Vsource_test_d, Vsource_test_q, - Vhist_d, Vhist_q, mIref_d, mIref_q, - **mPhi_d, **mPhi_q, **mGamma_d, **mGamma_q); - - return Complex(Vhist_d, Vhist_q); +Complex VSIControlType1::stepVBR(const Complex &Vcap_dq, + const Complex &Ifilter_dq) { + // store previous values (t=k-1) + Complex Vcap_dq_prev = mVcap_dq; + Complex Ifilter_dq_prev = mIfilter_dq; + Real Iref_d_prev = mIref_d; + Real Iref_q_prev = mIref_q; + + // update variables at t=k + mVcap_dq = Vcap_dq; + mIfilter_dq = Ifilter_dq; + + // calculate reference current at time t=k + mIref_d = mB2_VBR * mParameters->VdRef - mB_VBR * mVcap_dq.real() - + mC_VBR * Vcap_dq_prev.real() + **mPhi_d * mParameters->Kiv; + mIref_q = mB2_VBR * mParameters->VqRef - mB_VBR * mVcap_dq.imag() - + mC_VBR * Vcap_dq_prev.imag() + **mPhi_q * mParameters->Kiv; + + // Update phi at time t=k + **mPhi_d = + **mPhi_d + (mTimeStep / 2.) * (2 * mParameters->VdRef - mVcap_dq.real() - + Vcap_dq_prev.real()); + **mPhi_q = + **mPhi_q + (mTimeStep / 2.) * (2 * mParameters->VqRef - mVcap_dq.imag() - + Vcap_dq_prev.imag()); + + // Update lambda at time t=k + **mGamma_d = + **mGamma_d + (mTimeStep / 2.) * (mIref_d - mIfilter_dq.real() + + Iref_d_prev - Ifilter_dq_prev.real()); + **mGamma_q = + **mGamma_q + (mTimeStep / 2.) * (mIref_q - mIfilter_dq.imag() + + Iref_q_prev - Ifilter_dq_prev.imag()); + + // calculate Vvbr at t=k+1 + Real Vhist_d = + -mA_VBR * mC_VBR * mVcap_dq.real() - mD_VBR * mIfilter_dq.real() + + mD_VBR * mIref_d + mA_VBR * **mPhi_d * mParameters->Kiv + + **mGamma_d * mParameters->Kic + mA_VBR * mB2_VBR * mParameters->VdRef; + Real Vhist_q = + -mA_VBR * mC_VBR * mVcap_dq.imag() - mD_VBR * mIfilter_dq.imag() + + mD_VBR * mIref_q + mA_VBR * **mPhi_q * mParameters->Kiv + + **mGamma_q * mParameters->Kic + mA_VBR * mB2_VBR * mParameters->VqRef; + + // Test + Real Vsource_test_d = mVcap_dq.real() - mA_VBR * mB_VBR * mVcap_dq.real() - + mA_VBR * mIfilter_dq.real() + Vhist_d; + Real Vsource_test_q = mVcap_dq.imag() - mA_VBR * mB_VBR * mVcap_dq.imag() - + mA_VBR * mIfilter_dq.imag() + Vhist_q; + // + + SPDLOG_LOGGER_DEBUG(mSLog, + "\n--- Test Step: ---" + "\nmVcap_dq_d = {}" + "\nmVcap_dq_q = {}" + "\nIfilter_dq_d = {}" + "\nIfilter_dq_q = {}" + "\nVsource_test_d = {}" + "\nVsource_test_q = {}" + "\nVhist_d = {}" + "\nVhist_q = {}" + "\nIref_d = {}" + "\nIref_d = {}" + "\nmPhi_d = {}" + "\nmPhi_q = {}" + "\nmGamma_d = {}" + "\nmGamma_q = {}", + mVcap_dq.real(), mVcap_dq.imag(), Ifilter_dq.real(), + Ifilter_dq.imag(), Vsource_test_d, Vsource_test_q, + Vhist_d, Vhist_q, mIref_d, mIref_q, **mPhi_d, **mPhi_q, + **mGamma_d, **mGamma_q); + + return Complex(Vhist_d, Vhist_q); } \ No newline at end of file diff --git a/dpsim/examples/cxx/CMakeLists.txt b/dpsim/examples/cxx/CMakeLists.txt index 09dac0aec4..ba4c44026c 100644 --- a/dpsim/examples/cxx/CMakeLists.txt +++ b/dpsim/examples/cxx/CMakeLists.txt @@ -37,6 +37,7 @@ set(CIRCUIT_SOURCES Circuits/EMT_Slack_PiLine_PQLoad_with_PF_Init.cpp Circuits/EMT_Slack_PiLine_VSI_with_PF_Init.cpp Circuits/EMT_Slack_PiLine_VSI_Ramp_with_PF_Init.cpp + Circuits/EMT_Slack_PiLine_VSI_VoltageControlled_SteadyState_with_PF_Init.cpp # SP examples Circuits/SP_Circuits.cpp diff --git a/examples/Notebooks/Components/GridFormingInverter_ControllerType1.ipynb b/examples/Notebooks/Components/GridFormingInverter_ControllerType1.ipynb index f5c41e287a..73922d1398 100644 --- a/examples/Notebooks/Components/GridFormingInverter_ControllerType1.ipynb +++ b/examples/Notebooks/Components/GridFormingInverter_ControllerType1.ipynb @@ -27,7 +27,7 @@ "from villas.dataprocessing.readtools import *\n", "from villas.dataprocessing.timeseries import *\n", "\n", - "#%matplotlib widget" + "%matplotlib widget" ] }, { @@ -210,6 +210,7 @@ " vsi = dpsimpy_components.VSIVoltageControlDQ(uid=\"VSI\", name=\"VSI\", loglevel=log_level, \n", " with_interface_resistor=False,\n", " model_as_current_source=model_as_current_source)\n", + " \n", " vsi.set_parameters(sys_omega=nominal_omega, vdref=vsi_nominal_voltage, vqref=0)\n", " vsi.set_filter_parameters(Lf=Lf, Cf=Cf, Rf=Rf, Rc=Rc)\n", " \n", @@ -333,16 +334,18 @@ " \n", "def plot_emt_variable(varname_dpsim_emt, varname_dpsim_dp, varname_simulink, ts_dpsim_emt, ts_dpsim_dp, ts_dpsim_sp, ts_simulink, ylabels, y_lim, x_lim=[0.2, 0.8]):\n", " fig1 = plt.figure(figsize=(width, height))\n", - " dpsim_time = ts_dpsim_emt[varname_dpsim_emt].interpolate(timestep_common).time\n", - " dpsim_values_emt = ts_dpsim_emt[varname_dpsim_emt].interpolate(timestep_common).values\n", + " dpsim_time = ts_dpsim_dp[varname_dpsim_dp].interpolate(timestep_common).time\n", + " if (ts_dpsim_emt is not None):\n", + " dpsim_values_emt = ts_dpsim_emt[varname_dpsim_emt].interpolate(timestep_common).values\n", " dpsim_values_dp = np.sqrt(2/3) * ts_dpsim_dp[varname_dpsim_dp].interpolate(timestep_common).frequency_shift(60).values\n", " dpsim_values_sp = np.sqrt(2/3) * ts_dpsim_sp[varname_dpsim_dp].interpolate(timestep_common).frequency_shift(60).values\n", - " #time_simulink = ts_simulink[varname_simulink].interpolate(timestep_common).time\n", - " #simulink_values = ts_simulink[varname_simulink].interpolate(timestep_common).values\n", - " plt.plot(dpsim_time, dpsim_values_emt, label='DPSim - EMT')\n", + " time_simulink = ts_simulink[varname_simulink].interpolate(timestep_common).time\n", + " simulink_values = ts_simulink[varname_simulink].interpolate(timestep_common).values\n", + " if (ts_dpsim_emt is not None):\n", + " plt.plot(dpsim_time, dpsim_values_emt, label='DPSim - EMT')\n", " plt.plot(dpsim_time, dpsim_values_dp, '--', label='DPSim - DP')\n", " plt.plot(dpsim_time, dpsim_values_sp, '--', label='DPSim - SP')\n", - " #plt.plot(time_simulink, simulink_values, '-.', label='Simulink')\n", + " plt.plot(time_simulink, simulink_values, '-.', label='Simulink')\n", " plt.legend(loc='lower right')\n", " plt.xlabel('time (s)')\n", " plt.ylabel(ylabels)\n", @@ -374,9 +377,9 @@ "source": [ "# read Simulink log file\n", "\n", - "#simulink_log_name = \"Simulink_GridFormingInverter_Load_Fault.csv\"\n", - "#file_path_simulink = \"/home/mmo/git/inverter/dpsim/examples/Notebooks/Components/ReferenceResults/Simulink_GridFormingInverter_Load_Fault.csv\"\n", - "#ts_simulink = read_timeseries_dpsim(file_path_simulink)" + "simulink_log_name = \"Simulink_GridFormingInverter_Load_Fault.csv\"\n", + "file_path_simulink = \"/home/mmo/git/inverter/dpsim/examples/Notebooks/Components/ReferenceResults/Simulink_GridFormingInverter_Load_Fault.csv\"\n", + "ts_simulink = read_timeseries_dpsim(file_path_simulink)" ] }, { @@ -392,7 +395,7 @@ "metadata": {}, "outputs": [], "source": [ - "log_name = gridFormingInverter_Line_Load(domain=\"DP\", system_pf=system_pf, final_time=1, model_as_current_source=True, time_step=0.0001)\n", + "log_name = gridFormingInverter_Line_Load(domain=\"DP\", system_pf=system_pf, final_time=1, model_as_current_source=True, time_step=0.001)\n", "\n", "file_path = os.getcwd() + \"/logs/\" + log_name + \"/\" + log_name + \".csv\"\n", "ts_dpsim_dp = read_timeseries_dpsim(file_path)" @@ -411,9 +414,9 @@ "metadata": {}, "outputs": [], "source": [ - "log_name_emt = gridFormingInverter_Line_Load(domain=\"EMT\", system_pf=system_pf, final_time=1, model_as_current_source=True, time_step=0.0001)\n", - "file_path_emt = os.getcwd() + \"/logs/\" + log_name_emt + \"/\" + log_name_emt + \".csv\"\n", - "ts_dpsim_emt = read_timeseries_dpsim(file_path_emt)" + "#log_name_emt = gridFormingInverter_Line_Load(domain=\"EMT\", system_pf=system_pf, final_time=1, model_as_current_source=False, time_step=0.00001)\n", + "#file_path_emt = os.getcwd() + \"/logs/\" + log_name_emt + \"/\" + log_name_emt + \".csv\"\n", + "#ts_dpsim_emt = read_timeseries_dpsim(file_path_emt)" ] }, { @@ -429,7 +432,7 @@ "metadata": {}, "outputs": [], "source": [ - "log_name_sp = gridFormingInverter_Line_Load(domain=\"SP\", system_pf=system_pf, final_time=1, model_as_current_source=True, time_step=0.0001)\n", + "log_name_sp = gridFormingInverter_Line_Load(domain=\"SP\", system_pf=system_pf, final_time=1, model_as_current_source=True, time_step=0.001)\n", "file_path_sp = os.getcwd() + \"/logs/\" + log_name_sp + \"/\" + log_name_sp + \".csv\"\n", "ts_dpsim_sp = read_timeseries_dpsim(file_path_sp)" ] @@ -447,7 +450,9 @@ "metadata": {}, "outputs": [], "source": [ - "ts_simulink = None\n", + "#ts_simulink = None\n", + "ts_dpsim_emt = None\n", + "#ts_dpsim_sp = None\n", "plot_emt_variable('n2.v_0', 'n2.v', 'v2_1', ts_dpsim_emt, ts_dpsim_dp, ts_dpsim_sp, ts_simulink, 'V2', [-1200, 1200], [0.2, 0.6])" ] },