Commit 09739a05 authored by Viviane Sapucaia's avatar Viviane Sapucaia
Browse files

Added Turbine Governor model


Former-commit-id: 5c56e4b1
parent ce026ccd
/** Synchron Generator Tests
*
* @file
* @author Markus Mirz <mmirz@eonerc.rwth-aachen.de>
* @copyright 2017, Institute for Automation of Complex Power Systems, EONERC
* @license GNU General Public License (version 3)
*
* DPsim
*
* This program is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* any later version.
*
* This program is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this program. If not, see <http://www.gnu.org/licenses/>.
*********************************************************************************/
#include "Simulation.h"
using namespace DPsim;
int main() {
// Define Object for saving data on a file
Logger log("log.txt"),
vtLog("data_vt.csv"),
jLog("data_j.csv");
// Define machine parameters in per unit
Real nomPower = 555e6;
Real nomPhPhVoltRMS = 24e3;
Real nomFreq = 60;
Real nomFieldCurr = 1300;
Int poleNum = 2;
Real J = 2.8898e+04;
Real H = 3.7;
//Exciter
Real Ka = 20;
Real Ta = 0.2;
Real Ke = 1;
Real Te = 0.314;
Real Kf = 0.063;
Real Tf = 0.35;
Real Tr = 0.02;
// Turbine
Real Ta_t = 0.3;
Real Fa = 0.3;
Real Tb = 7;
Real Fb = 0.3;
Real Tc = 0.2;
Real Fc = 0.4;
Real Tsr = 0.1;
Real Tsm = 0.3;
Real Kg = 20;
Real Rs = 0.003;
Real Ll = 0.15;
Real Lmd = 1.6599;
Real Lmd0 = 1.6599;
Real Lmq = 1.61;
Real Lmq0 = 1.61;
Real Rfd = 0.0006;
Real Llfd = 0.1648;
Real Rkd = 0.0284;
Real Llkd = 0.1713;
Real Rkq1 = 0.0062;
Real Llkq1 = 0.7252;
Real Rkq2 = 0.0237;
Real Llkq2 = 0.125;
//Real Rkq2 = 0;
//Real Llkq2 = 0;
// Declare circuit components
ElementPtr gen = std::make_shared<VoltageBehindReactanceEMT>("gen", 1, 2, 3,
nomPower, nomPhPhVoltRMS, nomFreq, poleNum, nomFieldCurr,
Rs, Ll, Lmd, Lmd0, Lmq, Lmq0, Rfd, Llfd, Rkd, Llkd, Rkq1, Llkq1, Rkq2, Llkq2, H, true);
Real loadRes = 1037.8378;
ElementPtr r1 = std::make_shared<ResistorEMT>("r1", 1, 0, loadRes);
ElementPtr r2 = std::make_shared<ResistorEMT>("r2", 2, 0, loadRes);
ElementPtr r3 = std::make_shared<ResistorEMT>("r3", 3, 0, loadRes);
ElementList circElements;
circElements.push_back(gen);
circElements.push_back(r1);
circElements.push_back(r2);
circElements.push_back(r3);
// Declare circuit components for resistance change
Real breakerRes = 1037.8378;
ElementPtr rBreaker1 = std::make_shared<ResistorEMT>("rbreak1", 1, 0, breakerRes);
ElementPtr rBreaker2 = std::make_shared<ResistorEMT>("rbreak2", 2, 0, breakerRes);
ElementPtr rBreaker3 = std::make_shared<ResistorEMT>("rbreak3", 3, 0, breakerRes);
ElementList circElementsBreakerOn;
circElementsBreakerOn.push_back(gen);
circElementsBreakerOn.push_back(rBreaker1);
circElementsBreakerOn.push_back(rBreaker2);
circElementsBreakerOn.push_back(rBreaker3);
circElementsBreakerOn.push_back(r1);
circElementsBreakerOn.push_back(r2);
circElementsBreakerOn.push_back(r3);
// Set up simulation
Real tf, dt, t;
Real om = 2.0*M_PI*60.0;
tf = 10; dt = 0.0001; t = 0;
Int downSampling = 1;
Simulation newSim(circElements, om, dt, tf, log, SimulationType::EMT, downSampling);
newSim.setNumericalMethod(NumericalMethod::Trapezoidal_flux);
newSim.addSystemTopology(circElementsBreakerOn);
newSim.switchSystemMatrix(0);
// Initialize generator
Real initActivePower = 555e3;
Real initReactivePower = 0;
Real initTerminalVolt = 24000 / sqrt(3) * sqrt(2);
Real initVoltAngle = -DPS_PI / 2;
Real fieldVoltage = 7.0821;
Real mechPower = 5.5558e5;
shared_ptr<VoltageBehindReactanceEMT> genPtr = std::dynamic_pointer_cast<VoltageBehindReactanceEMT>(gen);
genPtr->init(om, dt, initActivePower, initReactivePower, initTerminalVolt, initVoltAngle, fieldVoltage, mechPower);
genPtr->AddExciter(Ta, Ka, Te, Ke, Tf, Kf, Tr, Lmd, Rfd);
genPtr->AddGovernor(Ta_t, Tb, Tc, Fa, Fb, Fc, Kg, Tsr, Tsm, initActivePower/ nomPower);
// Calculate initial values for circuit at generator connection point
Real initApparentPower = sqrt(pow(initActivePower, 2) + pow(initReactivePower, 2));
Real initTerminalCurr = initApparentPower / (3 * initTerminalVolt)* sqrt(2);
Real initPowerFactor = acos(initActivePower / initApparentPower);
std::cout << "A matrix:" << std::endl;
std::cout << newSim.getSystemMatrix() << std::endl;
std::cout << "vt vector:" << std::endl;
std::cout << newSim.getLeftSideVector() << std::endl;
std::cout << "j vector:" << std::endl;
std::cout << newSim.getRightSideVector() << std::endl;
Real lastLogTime = 0;
Real logTimeStep = 0.0001;
newSim.setSwitchTime(1, 1);
//newSim.setSwitchTime(6, 0);
// Main Simulation Loop
while (newSim.getTime() < tf) {
std::cout << newSim.getTime() << std::endl;
newSim.stepGeneratorTest(vtLog, jLog, gen, newSim.getTime());
newSim.increaseByTimeStep();
}
std::cout << "Simulation finished." << std::endl;
return 0;
}
......@@ -66,13 +66,13 @@ int main() {
Real Llkq2 = 0;
// Declare circuit components
ElementPtr gen = std::make_shared<VoltageBehindReactanceEMT>("gen", 1, 2, 3,
ElementPtr gen = std::make_shared<SynchronGeneratorDP>("gen", 1, 2, 3,
nomPower, nomPhPhVoltRMS, nomFreq, poleNum, nomFieldCurr,
Rs, Ll, Lmd, Lmd0, Lmq, Lmq0, Rfd, Llfd, Rkd, Llkd, Rkq1, Llkq1, Rkq2, Llkq2, H, true);
Real loadRes = 1037.8378;
ElementPtr r1 = std::make_shared<ResistorEMT>("r1", 1, 0, loadRes);
ElementPtr r2 = std::make_shared<ResistorEMT>("r2", 2, 0, loadRes);
ElementPtr r3 = std::make_shared<ResistorEMT>("r3", 3, 0, loadRes);
ElementPtr r1 = std::make_shared<ResistorDP>("r1", 1, 0, loadRes);
ElementPtr r2 = std::make_shared<ResistorDP>("r2", 2, 0, loadRes);
ElementPtr r3 = std::make_shared<ResistorDP>("r3", 3, 0, loadRes);
ElementList circElements;
circElements.push_back(gen);
......@@ -81,25 +81,25 @@ int main() {
circElements.push_back(r3);
// Declare circuit components for resistance change
Real breakerRes = 103.78378;
ElementPtr rBreaker1 = std::make_shared<ResistorEMT>("rbreak1", 1, 0, breakerRes);
ElementPtr rBreaker2 = std::make_shared<ResistorEMT>("rbreak2", 2, 0, breakerRes);
ElementPtr rBreaker3 = std::make_shared<ResistorEMT>("rbreak3", 3, 0, breakerRes);
Real breakerRes = 0.001;
ElementPtr rBreaker1 = std::make_shared<ResistorDP>("rbreak1", 1, 0, breakerRes);
ElementPtr rBreaker2 = std::make_shared<ResistorDP>("rbreak2", 2, 0, breakerRes);
ElementPtr rBreaker3 = std::make_shared<ResistorDP>("rbreak3", 3, 0, breakerRes);
ElementList circElementsBreakerOn;
circElementsBreakerOn.push_back(gen);
circElementsBreakerOn.push_back(rBreaker1);
circElementsBreakerOn.push_back(rBreaker2);
circElementsBreakerOn.push_back(rBreaker3);
//circElementsBreakerOn.push_back(r1);
//circElementsBreakerOn.push_back(r2);
//circElementsBreakerOn.push_back(r3);
circElementsBreakerOn.push_back(r1);
circElementsBreakerOn.push_back(r2);
circElementsBreakerOn.push_back(r3);
// Set up simulation
Real tf, dt, t;
Real om = 2.0*M_PI*60.0;
tf = 10; dt = 0.0001; t = 0;
Int downSampling = 1;
Simulation newSim(circElements, om, dt, tf, log, SimulationType::EMT, downSampling);
tf = 0.3; dt = 0.000001; t = 0;
Int downSampling = 50;
Simulation newSim(circElements, om, dt, tf, log, SimulationType::DynPhasor, downSampling);
newSim.setNumericalMethod(NumericalMethod::Trapezoidal_flux);
newSim.addSystemTopology(circElementsBreakerOn);
newSim.switchSystemMatrix(0);
......@@ -111,9 +111,9 @@ int main() {
Real initVoltAngle = -DPS_PI / 2;
Real fieldVoltage = 7.0821;
Real mechPower = 5.5558e5;
shared_ptr<VoltageBehindReactanceEMT> genPtr = std::dynamic_pointer_cast<VoltageBehindReactanceEMT>(gen);
shared_ptr<SynchronGeneratorDP> genPtr = std::dynamic_pointer_cast<SynchronGeneratorDP>(gen);
genPtr->init(om, dt, initActivePower, initReactivePower, initTerminalVolt, initVoltAngle, fieldVoltage, mechPower);
genPtr->AddExciter(Ta, Ka, Te, Ke, Tf, Kf, Tr, Lmd, Rfd);
//genPtr->AddExciter(Ta, Ka, Te, Ke, Tf, Kf, Tr, Lmd, Rfd);
// Calculate initial values for circuit at generator connection point
Real initApparentPower = sqrt(pow(initActivePower, 2) + pow(initReactivePower, 2));
......@@ -129,8 +129,8 @@ int main() {
Real lastLogTime = 0;
Real logTimeStep = 0.0001;
newSim.setSwitchTime(2, 1);
newSim.setSwitchTime(6, 0);
newSim.setSwitchTime(0.1, 1);
newSim.setSwitchTime(0.2, 0);
// Main Simulation Loop
while (newSim.getTime() < tf) {
......
......@@ -3,14 +3,14 @@
%% read PLECS results
Results_PLECS = csvread('../../../vsa/Results/SimpSynGenDqEmtExciter_LoadChange_Simulink/Voltages_and_currents.csv');
Te_PLECS = csvread('../../../vsa/Results/SimpSynGenDqEmtExciter_LoadChange_Simulink/electrical_torque.csv');
omega_PLECS = csvread('../../../vsa/Results/SimpSynGenDqEmtExciter_LoadChange_Simulink/omega.csv');
Results_PLECS = csvread('../../../vsa/Results/TestExciterAndTurbine/SynGenDqEMTExciter_loadChange_Simulink/Voltages_and_currents.csv');
%Te_PLECS = csvread('../../../vsa/Results/SynGenDqEmt_ABCFault_PLECS/electrical_torque.csv');
omega_PLECS = csvread('../../../vsa/Results/TestExciterAndTurbine/SynGenDqEMTExciter_loadChange_Simulink/omega.csv');
%theta_PLECS = csvread('../../vsa/Results/SynGenVBREmt_ABCFault_PLECS/theta.csv');
%% read results from c++ simulation
VoltageVector = csvread('../../../vsa/Results/Testing/data_vt.csv');
CurrentVector = csvread('../../../vsa/Results/Testing/data_j.csv');
Log_SynGen = csvread('../../../vsa/Results/Testing/SynGen_gen.csv');
VoltageVector = csvread('../../../vsa/Results/TestExciterAndTurbine/SynGenVBREmtExciter_LoadChange_DPsim/data_vt.csv');
CurrentVector = csvread('../../../vsa/Results/TestExciterAndTurbine/SynGenVBREmtExciter_LoadChange_DPsim/data_j.csv');
Log_SynGen = csvread('../../../vsa/Results/TestExciterAndTurbine/SynGenVBREmtExciter_LoadChange_DPsim/SynGen_gen.csv');
%% Plot
figure(1)
hold off
......@@ -71,21 +71,21 @@ legend('ia DPSim','ia Simulink');
figure(3)
hold off
plot(Log_SynGen(:,1),Log_SynGen(:,8));
plot(Log_SynGen(:,1),Log_SynGen(:,9));
hold on
plot(Results_PLECS(:,1),omega_PLECS);
title('Rotor speed');
legend('\omega DPSim','\omega PLECS');
legend('\omega DPSim','\omega Simulink');
% %
% figure(4)
% hold off
% plot(Log_SynGen(:,1),Log_SynGen(:,20));
% hold on
% plot(Results_PLECS(:,1),-Te_PLECS);
%
figure(4)
hold off
plot(Log_SynGen(:,1),Log_SynGen(:,7));
hold on
plot(Results_PLECS(:,1),-Te_PLECS);
title('Electrical Torque');
legend('Te DPSim','Te PLECS');
% title('Electrical Torque');
% legend('Te DPSim','Te PLECS');
%
% figure(5)
% hold off
......
......@@ -2,7 +2,7 @@ clc
clear
%% read PLECS results
Results_PLECS = csvread('../../vsa/Results/SynGenDqEmt_ABCFault_Simulink/Voltages_and_currents.csv');
Results_PLECS = csvread('../../../vsa/Results/SynGenDqEmt_ABCFault_PLECS/Voltages_and_currents.csv');
%Te_PLECS = csvread('../../vsa/Results/SynGenVBREmt_ABCFault_PLECS/electrical_torque.csv');
%omega_PLECS = csvread('../../vsa/Results/SynGenVBREmt_ABCFault_PLECS/omega.csv');
%theta_PLECS = csvread('../../vsa/Results/SynGenVBREmt_ABCFault_PLECS/theta.csv');
......@@ -10,8 +10,8 @@ Results_PLECS = csvread('../../vsa/Results/SynGenDqEmt_ABCFault_Simulink/Voltage
%% Read data from DP simulation and calculate absolute value and phase
% Read values from CSV files
voltageDP = csvread('../../vsa/Results/SynGenVBRDynPh_ABCFault_DPsim/data_vt.csv');
currentDP = csvread('../../vsa/Results/SynGenVBRDynPh_ABCFault_DPsim/data_j.csv');
voltageDP = csvread('../../../vsa/Results/Testing/data_vt.csv');
currentDP = csvread('../../../vsa/Results/Testing/data_j.csv');
%Log_SynGen = csvread('../../vsa/Results//SynGenDqDynPh_ABCFault_DPsim/Euler/SynGen_gen.csv');
compOffsetDP = (size(voltageDP,2) - 1) / 2;
......@@ -97,7 +97,7 @@ figure(4)
hold off
PLECSplotc = plot(Results_PLECS(:,1), Results_PLECS(:,5), '--');
hold on
DPplotc = plot(currentShiftDP(:,1),currentShiftDP(:,2));
DPplotc = plot(currentShiftDP(:,1),-currentShiftDP(:,2));
DPabsPlotc = plot(currentAbsDP(:,1),currentAbsDP(:,2));
title('Current phase A');
legend('Current Phase a Simulink', 'DP shift a', 'DP abs a')
......@@ -109,7 +109,7 @@ figure(5)
hold off
PLECSplot2c = plot(Results_PLECS(:,1), Results_PLECS(:,6), '--');
hold on
DPplot2c = plot(currentShiftDP(:,1),currentShiftDP(:,3));
DPplot2c = plot(currentShiftDP(:,1),-currentShiftDP(:,3));
DPabsPlot2c = plot(currentAbsDP(:,1),currentAbsDP(:,3));
title('Current phase B');
legend('Current Phase b Simulink', 'DP shift b', 'DP abs b')
......@@ -121,7 +121,7 @@ figure(6)
hold off
PLECSplot3c = plot(Results_PLECS(:,1), Results_PLECS(:,7), '--');
hold on
DPplot3c = plot(currentShiftDP(:,1),currentShiftDP(:,4));
DPplot3c = plot(currentShiftDP(:,1),-currentShiftDP(:,4));
DPabsPlot3c = plot(currentAbsDP(:,1),currentAbsDP(:,4));
title('Currents phase C');
legend('Current Phase c Simulink', 'DP shift c', 'DP abs c')
......
[1207/155452:WARNING:resource_bundle.cc(291)] locale_file_path.empty() for locale
[1207/155452:WARNING:resource_bundle.cc(291)] locale_file_path.empty() for locale
......@@ -38,6 +38,7 @@ Exciter::Exciter(Real Ta, Real Ka, Real Te, Real Ke, Real Tf, Real Kf, Real Tr,
mLad = Lad;
mRfd = Rfd;
mVf = Vf_init*(mLad / mRfd);
mVf = 0;
}
......@@ -47,15 +48,15 @@ Real Exciter::step(Real mVd, Real mVq, Real Vref, Real dt) {
mVh = sqrt(pow(mVd, 2.) + pow(mVq, 2.));
// Voltage Transducer equation
mVm = Trapezoidal(mVm, -1, 1, dt / mTr, mVh);
mVm = Euler(mVm, -1, 1, dt / mTr, mVh);
// Stabilizing feedback equation
mVfl = mVfl + dt*mVis / mTf;
mVis = (mKf / mTf)*mVf - mVfl;
// Amplifier equation
mVr = Trapezoidal(mVr, -1, mKa, dt / mTa, Vref - mVm - mVis);
mVr = Euler(mVr, -1, mKa, dt / mTa, Vref - mVm - mVis);
// Exciter
mVse = 0.0039*exp(mVf*1.555);
mVf = Trapezoidal(mVf, -mKe, 1, dt / mTe, mVr - mVse);
mVf = Euler(mVf, -mKe, 1, dt / mTe, mVr - mVse);
return (mRfd / mLad)*mVf;
......
/** Turbine Governor
*
* @author Markus Mirz <mmirz@eonerc.rwth-aachen.de>
* @copyright 2017, Institute for Automation of Complex Power Systems, EONERC
* @license GNU General Public License (version 3)
*
* DPsim
*
* This program is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* any later version.
*
* This program is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this program. If not, see <http://www.gnu.org/licenses/>.
*********************************************************************************/
#include "TurbineGovernor.h"
#include "../IntegrationMethod.h"
using namespace DPsim;
TurbineGovernor::TurbineGovernor(Real Ta, Real Tb, Real Tc, Real Fa, Real Fb, Real Fc, Real K, Real Tsr, Real Tsm, Real Tm_init)
{
mTa = Ta;
mTb = Tb;
mTc = Tc;
mFa = Fa;
mFb = Fb;
mFc = Fc;
mK = K;
mTsr = Tsr;
mTsm = Tsm;
mTm = Tm_init;
}
Real TurbineGovernor::step(Real Om, Real OmRef, Real PmRef, Real dt) {
// ### Governing ###
// Input of speed relay
Psr_in = PmRef + (OmRef - Om)*mK;
// Input of servor motor
Psm_in = Euler(Psm_in, -1, 1, dt / mTsr, Psr_in);
// rate of change of valve
mpVcv = (Psm_in - mVcv) / mTsm;
if (mpVcv >= 0.1)
mpVcv = 0.1;
else if (mpVcv <= -1)
mpVcv = -1;
//Valve position
mVcv = mVcv + dt*mpVcv;
if (mVcv >= 1)
mVcv = 1;
else if (mVcv <= 0)
mVcv = 0;
//### Turbine ###
// Simplified equation
mTm = Euler(mTm, -1 / mTb, 1 / mTb, mpVcv*mFa, dt, mVcv);
return mTm;
}
/** Turbine Governor
*
* @file
* @author Markus Mirz <mmirz@eonerc.rwth-aachen.de>
* @copyright 2017, Institute for Automation of Complex Power Systems, EONERC
* @license GNU General Public License (version 3)
*
* DPsim
*
* This program is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* any later version.
*
* This program is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this program. If not, see <http://www.gnu.org/licenses/>.
*********************************************************************************/
#pragma once
#include "BaseComponent.h"
namespace DPsim {
/// Synchronous generator model
/// If parInPerUnit is not set, the parameters have to be given with their respective stator or rotor
/// referred values. The calculation to per unit is performed in the initialization.
/// The case where parInPerUnit is not set will be implemented later.
/// parameter names include underscores and typical variables names found in literature instead of
/// descriptive names in order to shorten formulas and increase the readability
class TurbineGovernor {
protected:
// ### Steam Turbine Parameters ####
/// Time constant of main inlet volume and steam chest
Real mTa;
/// Time constant of reheater
Real mTb;
/// Time constant of cross over piping and LP inlet volumes
Real mTc;
/// Fraction of total turbine power generated by HP section
Real mFa;
/// Fraction of total turbine power generated by IP section
Real mFb;
/// Fraction of total turbine power generated by LP section
Real mFc;
// ### Regulator ####
/// Main droop
Real mK;
// ### Speed Relay and servo motor ####
/// Time constant of speed relay
Real mTsr;
/// Time constant of servo motor
Real mTsm;
/// Opening rate limit
Real mLc1 = 0.1;
/// Closing rate limit
Real mLc2 = -1;
// ### Vaariables ###
/// Mechanical Torque in pu (output of steam turbine)
Real mTm = 0;
/// Valve position
Real mVcv = 0;
/// Valve position changing rate
Real mpVcv = 0;
/// Boiler pressure
Real mPb = 1;
/// Speed reference
Real mOmRef = 1;
/// Speed
Real mOm = 0;
/// Power Reference
Real mPmRef;
/// Input of speed realy
Real Psr_in = 0;
/// Input of servor motor
Real Psm_in = 0;
bool mLogActive;
Logger* mLog;
public:
TurbineGovernor() { };
~TurbineGovernor() {};
/// Initializes exciter parameters
TurbineGovernor(Real Ta, Real Tb, Real Tc, Real Fa, Real Fb, Real Fc, Real K, Real Tsr, Real Tsm, Real Tm_init);
/// Performs an step to update field voltage value
Real step(Real mOm, Real mOmRef, Real PmRef, Real dt);
};
}
......@@ -48,6 +48,11 @@ void VoltageBehindReactanceEMT::AddExciter(Real Ta, Real Ka, Real Te, Real Ke, R
WithExciter = true;
}
void VoltageBehindReactanceEMT::AddGovernor(Real Ta, Real Tb, Real Tc, Real Fa, Real Fb, Real Fc, Real K, Real Tsr, Real Tsm, Real Tm_init){
mTurbineGovernor = TurbineGovernor(Ta, Tb, Tc, Fa, Fb, Fc, K, Tsr, Tsm, Tm_init);
WithTurbineGovernor = true;
}
void VoltageBehindReactanceEMT::init(Real om, Real dt,
......@@ -110,6 +115,7 @@ void VoltageBehindReactanceEMT::init(Real om, Real dt,
// steady state per unit initial value
initStatesInPerUnit(initActivePower, initReactivePower, initTerminalVolt, initVoltAngle, initFieldVoltage, initMechPower);
// Correcting variables
mThetaMech = mThetaMech + PI/2;
mMechTorque = -mMechTorque;
......@@ -191,6 +197,11 @@ void VoltageBehindReactanceEMT::stepInPerUnit(Real om, Real dt, Real time, Numer
mIc;
// Calculate mechanical variables with euler
if (WithTurbineGovernor == true)
{
mMechTorque = -mTurbineGovernor.step(mOmMech, 1, 0.001, dt);
}
mElecTorque = (mPsimd*mIq - mPsimq*mId);
mOmMech = mOmMech + dt * (1. / (2. * mH) * (mElecTorque - mMechTorque));
mThetaMech = mThetaMech + dt * (mOmMech* mBase_OmMech);
......
......@@ -25,6 +25,7 @@
#include "SynchGenBase.h"
#include "Exciter.h"
#include "TurbineGovernor.h"
namespace DPsim {
......@@ -42,6 +43,10 @@ namespace DPsim {
Exciter mExciter;