# 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/.
# SPDX-License-Identifier: MPL-2.0
from typing import List
import numpy as np
import math
from VeraGridEngine.enumerations import DeviceType
from VeraGridEngine.Devices.Dynamic.emt_template import EmtModelTemplate
from VeraGridEngine.Utils.Symbolic.block import (Block, Var, VarPowerFlowReferenceType)
from VeraGridEngine.Utils.Symbolic.block_helpers import tf_to_block
from VeraGridEngine.Devices.Dynamic.var_factory import VarFactory
import VeraGridEngine.Utils.Symbolic.symbolic as sym
from VeraGridEngine.enumerations import ConverterControlType
[docs]
def build_gfl_converter_model_emt(vfactory: VarFactory, inputs,
control1: ConverterControlType = ConverterControlType.Pac,
control2: ConverterControlType = ConverterControlType.Qac,
multilinear:bool = False):
"""
Build power control loop model for Grid Following Converter for EMT simulation.
Supports multiple control modes via ConverterControlType.
Args:
inputs: [Va, Vb, Vc, Vdc] - AC instantaneous three-phase voltages, DC voltage
control1: First control mode (typically active power or DC voltage related)
control2: Second control mode (typically reactive power or AC voltage related)
multilinear: Use multilinearization for Park transforms
Returns:
gfl_block: The complete converter block
i_a, i_b, i_c: Three-phase current outputs
P, Q: Active and reactive power measurements
"""
v_a = inputs[0]
v_b = inputs[1]
v_c = inputs[2]
v_dc = inputs[3] # DC voltage for DC voltage control mode
Pt_vsc = inputs[4] # Power flow initial values
Qt_vsc = inputs[5]
algebraic_eqs = []
algebraic_vars = []
gfl_block = Block()
control_blocks = []
# ==============================
# Inputs and variables
# ==============================
v_abc = [v_a, v_b, v_c]
pll_block, v_dq, omega, theta, aux_vars = pll_transform(vfactory, v_abc, multilinear=multilinear, name = 'vg')
v_d_g = v_dq[0]
v_q_g = v_dq[1]
# Parameters
Kp_icl = vfactory.add_var('Kp_icl') # proportional gain for inner current loop
Ki_icl = vfactory.add_var('Ki_icl') # integral gain for inner current loop
Kp_pol = vfactory.add_var('Kp_pol') # proportional gain for outer power loop
Ki_pol = vfactory.add_var('Ki_pol') # integral gain for outer power loop
Kp_vdc = vfactory.add_var('Kp_vdc') # proportional gain for DC voltage control
Ki_vdc = vfactory.add_var('Ki_vdc') # integral gain for DC voltage control
Kp_vac = vfactory.add_var('Kp_vac') # proportional gain for AC voltage control
Ki_vac = vfactory.add_var('Ki_vac') # integral gain for AC voltage control
R = vfactory.add_var('R')
L = vfactory.add_var('L')
# VARS
i_d = vfactory.add_var('i_d')
i_q = vfactory.add_var('i_q')
i_d_ref = vfactory.add_var('i_d_ref')
i_q_ref = vfactory.add_var('i_q_ref')
# POWER LOOP
P = vfactory.add_var('P')
Q = vfactory.add_var('Q')
P_ref = vfactory.add_var('P_ref')
Q_ref = vfactory.add_var('Q_ref')
# Voltage references for control modes
Vdc_ref = vfactory.add_var('Vdc_ref')
Vm_ac_ref = vfactory.add_var('Vm_ac_ref')
v_d_c = vfactory.add_var('v_d_c')
v_q_c = vfactory.add_var('v_q_c')
# Three-phase currents
i_a = vfactory.add_var('i_a')
i_b = vfactory.add_var('i_b')
i_c = vfactory.add_var('i_c')
event_dict = {
Kp_icl: vfactory.add_const(0.05),
Ki_icl: vfactory.add_const(1.0),
Kp_pol: vfactory.add_const(0.05),
Ki_pol: vfactory.add_const(1.0),
Kp_vdc: vfactory.add_const(0.1),
Ki_vdc: vfactory.add_const(2.0),
Kp_vac: vfactory.add_const(0.1),
Ki_vac: vfactory.add_const(2.0),
R: vfactory.add_const(0.1),
L: vfactory.add_const(0.1),
P_ref: vfactory.add_const(0.0),
Q_ref: vfactory.add_const(0.0),
Vdc_ref: vfactory.add_const(1.0),
Vm_ac_ref: vfactory.add_const(1.0),
}
# P and Q at the point of common coupling
algebraic_eqs.append(P - vfactory.add_const(3/2)*(v_q_g*i_q + v_d_g*i_d))
algebraic_eqs.append(Q - vfactory.add_const(3/2)*(v_q_g*i_d - v_d_g*i_q))
algebraic_vars.append(P)
algebraic_vars.append(Q)
# ==============================
# CONTROL 1: Active Power Axis (i_q_ref)
# ==============================
# Control1 maps to i_q_ref (q-axis current controls active power)
if control1 == ConverterControlType.Pac:
# Active power control at AC side
control_block_1, _ = tf_to_block(vfactory,
num=[Ki_pol, Kp_pol],
den=[0, 1],
x= P_ref - P,
y = i_q_ref,
name='Pac_ctrl'
)
control_blocks.append(control_block_1)
elif control1 == ConverterControlType.Pdc:
# Active power control at DC side (Pdc ~ Pac in steady state)
# For simplicity, we control Pac which is approximately Pdc
control_block_1, _ = tf_to_block(vfactory,
num=[Ki_pol, Kp_pol],
den=[0, 1],
x= P_ref - P, # P_ref represents desired Pdc
y = i_q_ref,
name='Pdc_ctrl'
)
control_blocks.append(control_block_1)
elif control1 == ConverterControlType.Vm_dc:
# DC voltage control - regulates v_dc to Vdc_ref
# Positive i_q_ref increases power flow from DC to AC, reducing v_dc
control_block_1, _ = tf_to_block(vfactory,
num=[Ki_vdc, Kp_vdc],
den=[0, 1],
x= Vdc_ref - v_dc,
y = i_q_ref,
name='Vdc_ctrl'
)
control_blocks.append(control_block_1)
else:
raise ValueError(f"Control1 type {control1} not supported for GFL converter. "
f"Supported: Pac, Pdc, Vm_dc")
# ==============================
# CONTROL 2: Reactive Power Axis (i_d_ref)
# ==============================
# Control2 maps to i_d_ref (d-axis current controls reactive power)
if control2 == ConverterControlType.Qac:
# Reactive power control
control_block_2, _ = tf_to_block(vfactory,
num=[Ki_pol, Kp_pol],
den=[0, 1],
x= Q_ref - Q,
y = i_d_ref,
name='Qac_ctrl'
)
control_blocks.append(control_block_2)
elif control2 == ConverterControlType.Vm_ac:
# AC voltage magnitude control
# Use v_q_g as the measured voltage magnitude in dq frame (vd ~ 0)
control_block_2, _ = tf_to_block(vfactory,
num=[Ki_vac, Kp_vac],
den=[0, 1],
x= Vm_ac_ref - v_q_g, # v_q_g is approximately Vm in dq frame
y = i_d_ref,
name='Vac_ctrl'
)
control_blocks.append(control_block_2)
else:
raise ValueError(f"Control2 type {control2} not supported for GFL converter. "
f"Supported: Qac, Vm_ac")
# Physical Current Limits, TODO add AntiWindup
I_max = vfactory.add_const(1.2)
operation = 'normal'
if operation == 'normal':
id_max = sym.sqrt(sym.max(I_max**2 - sym.max(i_q, i_q_ref)**2, vfactory.add_const(1e-5)))
i_d_ref_sat = sym.hard_sat(i_d_ref, -id_max, id_max)
i_q_ref_sat = sym.hard_sat(i_q_ref, -I_max, I_max)
# Voltage Control Loop (Inner Current Loop)
control_block_iq , vq_hat = tf_to_block(vfactory,
num=[Ki_icl, Kp_icl],
den=[0, 1],
x= i_q_ref_sat - i_q,
name='vq_hat'
)
control_block_id , vd_hat = tf_to_block(vfactory,
num=[Ki_icl, Kp_icl],
den=[0, 1],
x= i_d_ref_sat - i_d,
name='vd_hat'
)
# EMT block (differential equations for currents)
dt_id = vfactory.add_diff_var('dt_id', base_var= i_d)
dt_iq = vfactory.add_diff_var('dt_iq', base_var= i_q)
EMT_block=Block(
algebraic_eqs = [
v_d_c - v_d_g - (R*i_d + L*dt_id - omega*L*i_q),
v_q_c - v_q_g - (R*i_q + L*dt_iq - omega*L*i_d),
],
diff_vars= [dt_id, dt_iq]
)
algebraic_eqs.append(v_d_c - (vd_hat + v_d_g + L*(omega)*i_q))
algebraic_eqs.append(v_q_c - (vq_hat + v_q_g - L*(omega)*i_d))
algebraic_vars.extend([v_q_c, v_d_c, i_d, i_q])
# Inverse Park transform to get three-phase currents
inv_park_block, i_abc, _ = inverse_park_transform_block(vfactory, [i_d, i_q], theta, name='i')
i_a_out, i_b_out, i_c_out = i_abc
algebraic_eqs.extend([
i_a - i_a_out,
i_b - i_b_out,
i_c - i_c_out,
])
algebraic_vars.extend([i_a, i_b, i_c])
# Build initialization equations based on control modes
# P and Q are initialized from power flow results via external mapping
init_eqs = {
P: -Pt_vsc,
Q: -Qt_vsc,
theta: vfactory.add_const(0), # EMT: initialize theta to 0 (will track actual angle)
omega: vfactory.add_const(1),
v_d_g: vfactory.add_const(0),
v_q_g: sym.sqrt(v_a**2 + v_b**2 + v_c**2) / vfactory.add_const(np.sqrt(3/2)), # Approximate Vm from abc
i_q: (2 / 3) * P / v_q_g,
i_d: (2 / 3) * Q / v_q_g,
i_q_ref: i_q,
i_d_ref: i_d,
v_d_c: v_d_g + (R*i_d - omega*L*i_q),
v_q_c: v_q_g + (R*i_q - omega*L*i_d),
vd_hat: v_d_c - (v_d_g + L*(omega)*i_q),
vq_hat: v_q_c - (v_q_g - L*(omega)*i_d),
}
# Add control-specific initialization
if control1 in [ConverterControlType.Pac, ConverterControlType.Pdc]:
init_eqs[P_ref] = P
elif control1 == ConverterControlType.Vm_dc:
init_eqs[Vdc_ref] = v_dc
# For DC voltage control, initialize i_q_ref based on initial power flow
init_eqs[i_q_ref] = (2 / 3) * P / v_q_g
if control2 == ConverterControlType.Qac:
init_eqs[Q_ref] = Q
elif control2 == ConverterControlType.Vm_ac:
init_eqs[Vm_ac_ref] = v_q_g
# For AC voltage control, initialize i_d_ref based on initial reactive power
init_eqs[i_d_ref] = (2 / 3) * Q / v_q_g
gfl_block_aux = Block(
algebraic_eqs=algebraic_eqs,
algebraic_vars=algebraic_vars,
event_dict= event_dict,
init_eqs=init_eqs,
external_mapping={
VarPowerFlowReferenceType.P:P,
VarPowerFlowReferenceType.Q:Q,
}
)
gfl_block.add(gfl_block_aux)
# Add all control blocks
for ctrl_block in control_blocks:
gfl_block.add(ctrl_block)
gfl_block.add(control_block_id)
gfl_block.add(control_block_iq)
gfl_block.add(pll_block)
gfl_block.add(EMT_block)
gfl_block.add(inv_park_block)
gfl_block.unify_blocks()
return gfl_block, i_a, i_b, i_c, P, Q
[docs]
def VscGflEmtBuild(vfactory: VarFactory, name: str = "",
control1: ConverterControlType = ConverterControlType.Pac,
control2: ConverterControlType = ConverterControlType.Qac) -> EmtModelTemplate:
"""
VSC GFL (Grid Following) EMT model
with from side the DC bus and to side the AC bus
Args:
name: Model name
control1: First control mode (Pac, Pdc, or Vm_dc)
control2: Second control mode (Qac or Vm_ac)
Supported control combinations:
- Pac + Qac: Active and reactive power control
- Pac + Vm_ac: Active power and AC voltage control
- Pdc + Qac: DC power and reactive power control
- Pdc + Vm_ac: DC power and AC voltage control
- Vm_dc + Qac: DC voltage and reactive power control
- Vm_dc + Vm_ac: DC voltage and AC voltage control
"""
templ = EmtModelTemplate()
templ.tpe = DeviceType.VscDevice
# Inputs: Va, Vb, Vc, Vdc (Pt_vsc and Qt_vsc come from power flow initialization)
inputs = [vfactory.add_var("Va_"), vfactory.add_var("Vb_"), vfactory.add_var("Vc_"), vfactory.add_var("Vdc_")]
v_a = inputs[0]
v_b = inputs[1]
v_c = inputs[2]
v_dc = inputs[3]
# Vars:
Pt_vsc = vfactory.add_var('Pt_vsc')
Qt_vsc = vfactory.add_var('Qt_vsc')
# Parameters:
bt = vfactory.add_var('bt')
gt = vfactory.add_var('gt')
Qf = vfactory.add_var('Qf')
a0 = vfactory.add_var('a0')
a1 = vfactory.add_var('a1')
a2 = vfactory.add_var('a2')
# Build the converter model with specified control modes
gfl_block, i_a, i_b, i_c, P, Q = build_gfl_converter_model_emt(
vfactory=vfactory,
inputs=[v_a, v_b, v_c, v_dc, Pt_vsc, Qt_vsc],
control1=control1,
control2=control2
)
event_dict = {
Qf: vfactory.add_const(0.0),
bt: vfactory.add_const(0.0),
gt: vfactory.add_const(0.1),
a0: vfactory.add_const(0.0),
a1: vfactory.add_const(0.0),
a2: vfactory.add_const(0.0),
}
# EMT model outputs three-phase currents
vsc_block = Block(
algebraic_eqs=[
Pt_vsc + P,
Qt_vsc + Q,
],
algebraic_vars=[Pt_vsc, Qt_vsc],
event_dict= event_dict,
external_mapping={
VarPowerFlowReferenceType.P: P,
VarPowerFlowReferenceType.Q: Q,
},
in_vars= inputs,
out_vars = [i_a, i_b, i_c] # EMT: output three-phase currents
)
vsc_block.external_mapping = {
VarPowerFlowReferenceType.v_A: inputs[0],
VarPowerFlowReferenceType.v_B: inputs[1],
VarPowerFlowReferenceType.v_C: inputs[2],
VarPowerFlowReferenceType.Vdc: inputs[3],
VarPowerFlowReferenceType.Pt: Pt_vsc,
VarPowerFlowReferenceType.Qt: Qt_vsc,
VarPowerFlowReferenceType.Pf: P,
VarPowerFlowReferenceType.Qf: Qf,
VarPowerFlowReferenceType.i_A: i_a,
VarPowerFlowReferenceType.i_B: i_b,
VarPowerFlowReferenceType.i_C: i_c,
}
vsc_block.add(gfl_block)
vsc_block.name = 'gfl_block_emt'
templ.block = vsc_block
return templ