From fca578d4813e4a4e3c826d734725035aedf591b3 Mon Sep 17 00:00:00 2001 From: "martin.moraga" Date: Fri, 22 Mar 2024 09:20:41 +0100 Subject: [PATCH] add VBR GridFormingInverter SP domain Signed-off-by: martin.moraga --- .../SP/SP_Ph1_VSIVoltageControlVBR.h | 100 +++++++ .../src/SP/SP_Ph1_VSIVoltageControlVBR.cpp | 281 ++++++++++++++++++ 2 files changed, 381 insertions(+) create mode 100644 dpsim-models/include/dpsim-models/SP/SP_Ph1_VSIVoltageControlVBR.h create mode 100644 dpsim-models/src/SP/SP_Ph1_VSIVoltageControlVBR.cpp diff --git a/dpsim-models/include/dpsim-models/SP/SP_Ph1_VSIVoltageControlVBR.h b/dpsim-models/include/dpsim-models/SP/SP_Ph1_VSIVoltageControlVBR.h new file mode 100644 index 0000000000..4d680c8809 --- /dev/null +++ b/dpsim-models/include/dpsim-models/SP/SP_Ph1_VSIVoltageControlVBR.h @@ -0,0 +1,100 @@ +/* Copyright 2017-2021 Institute for Automation of Complex Power Systems, + * EONERC, RWTH Aachen University + * + * This Source Code Form is subject to the terms of the Mozilla Public + * License, v. 2.0. If a copy of the MPL was not distributed with this + * file, You can obtain one at https://mozilla.org/MPL/2.0/. + *********************************************************************************/ +#pragma once + +#include +#include +#include +#include + +namespace CPS +{ + namespace SP + { + namespace Ph1 + { + class VSIVoltageControlVBR : public CompositePowerComp, + public Base::VSIVoltageSourceInverterDQ, + public MNAVariableCompInterface, + public SharedFactory + { + public: + /// Defines name amd logging level + VSIVoltageControlVBR(String name, Logger::Level logLevel = Logger::Level::off) + : VSIVoltageControlVBR(name, name, logLevel) {} + /// Defines UID, name, logging level and connection trafo existence + VSIVoltageControlVBR(String uid, String name, Logger::Level logLevel = Logger::Level::off, + Bool modelAsCurrentSource = true); + + // #### General #### + /// Initializes component from power flow data + void initializeFromNodesAndTerminals(Real frequency) final; + + // #### MNA section #### + /// Initializes internal variables of the component + void mnaParentInitialize(Real omega, Real timeStep, Attribute::Ptr leftVector) final; + /// Stamps system matrix + void mnaParentApplySystemMatrixStamp(SparseMatrixRow &systemMatrix) final; + /// Add MNA pre step dependencies + void mnaParentAddPreStepDependencies(AttributeBase::List &prevStepDependencies, AttributeBase::List &attributeDependencies, AttributeBase::List &modifiedAttributes) final; + /// MNA pre step operations + void mnaParentPreStep(Real time, Int timeStepCount) final; + /// Add MNA post step dependencies + void mnaParentAddPostStepDependencies(AttributeBase::List &prevStepDependencies, AttributeBase::List &attributeDependencies, AttributeBase::List &modifiedAttributes, Attribute::Ptr &leftVector) final; + /// Updates current through the component + void mnaCompUpdateCurrent(const Matrix &leftVector) final; + /// Updates voltage across component + void mnaCompUpdateVoltage(const Matrix &leftVector) final; + /// MNA post step operations + void mnaParentPostStep(Real time, Int timeStepCount, Attribute::Ptr &leftVector) final; + /// + void mnaParentApplyRightSideVectorStamp(Matrix &rightVector) final; + + /// Mark that parameter changes so that system matrix is updated + Bool hasParameterChanged() final { return true; }; + + protected: + // ### Electrical Subcomponents ### + /// Capacitor Cf as part of LC filter + std::shared_ptr mSubCapacitorF; + + private: + /// + void createSubComponents() final; + /// + void updatePower(); + + /// + void calculateResistanceMatrix(); + /// + void update_DqToComplexATransformMatrix(); + + private: + // ### VBR Model specific variables + /// History voltage in dq domain (output of VSI_Controller) + Complex mVhist; + /// Current flowing throw the filter inductor in dp domain + MatrixComp mFilterCurrent; + + /// Auxiliar Matrix + Matrix mA; + Matrix mE; + Matrix mY; + Matrix mAMatrix; + Matrix mEMatrix; + Complex mAdmitance; + + /// Park Transformation + /// + Matrix mDqToComplexA; + /// + Matrix mComplexAToDq; + }; + } + } +} diff --git a/dpsim-models/src/SP/SP_Ph1_VSIVoltageControlVBR.cpp b/dpsim-models/src/SP/SP_Ph1_VSIVoltageControlVBR.cpp new file mode 100644 index 0000000000..06b44237ed --- /dev/null +++ b/dpsim-models/src/SP/SP_Ph1_VSIVoltageControlVBR.cpp @@ -0,0 +1,281 @@ +/* Copyright 2017-2021 Institute for Automation of Complex Power Systems, + * EONERC, RWTH Aachen University + * + * This Source Code Form is subject to the terms of the Mozilla Public + * License, v. 2.0. If a copy of the MPL was not distributed with this + * file, You can obtain one at https://mozilla.org/MPL/2.0/. + *********************************************************************************/ + +#include +#include + +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); +} + +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(); +} + +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(); +} + +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); +} + +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); +} + +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); + + // TODO: droop + // if (mWithDroop) + // mDroop->signalStep(time, timeStepCount); + + // VCO Step + **mThetaInv = **mThetaInv + mTimeStep * **mOmega; + + // Update nominal system angle + **mThetaSys = **mThetaSys + mTimeStep * mOmegaNom; + + // calculat history term + if (mWithControl) + mVhist = mVSIController->stepVBR(**mVcap_dq, **mIfilter_dq); + + update_DqToComplexATransformMatrix(); + mComplexAToDq = mDqToComplexA.transpose(); + + // calculate resistance matrix at t=k+1 + calculateResistanceMatrix(); + + 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); +} + +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); +} + +void SP::Ph1::VSIVoltageControlVBR::mnaParentPostStep( + 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; +} + +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()); +} + +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); +}