Engee documentation
Notebook

Implicit generator without regulators

Model description

This example considers the operation of a turbo-generator-transformer unit through a two-circuit overhead line (OHL) to an infinite power busbar (IPB) without excitation and speed controllers. At a time instant of 5 s. in the model, a three-phase short circuit (SC) occurs at the beginning of one of the overhead line circuits. Analytical calculation of the short-circuit tripping time limit using the method of areas and comparison with the results of mathematical modelling [1-3] is carried out.

The process of starting and setting up the model from the script development environment using command control, processing of simulation results, visualisation of simulation results and suggested scripts for independent work with the model are shown. Logging of signals is carried out, their time plots are shown. External view of the model:

synchronous_generator_round_rotor_1737723205823.png

The main blocks used in this example and their parameters:

  1. Synchronous Machine Round Rotor - implicitly pole synchronous generator ($P_{ном} = 500 МВт$,$U_{ном} = 24 кВ$).
  2. Voltage Source (Three-Phase) - a three-phase voltage source, used to simulate the SWBM, the steady state is set by specifying the effective line voltage and phase shift ($E = 500\angle0° кВ$,$S_{кз} = 10 ГВт$).
  3. Fault (Three-Phase) - short circuit.
  4. From Library of Passive Elements blocks of three-phase transmission line Three-Phase PI Section Line ($L = 200км$, AC 500/64, 2 circuits) and two-winding three-phase transformer Two-Winding Transformer (Three-Phase) (TD-630000/500).

Calculation of the tripping angle limit by the area method

Import the necessary modules for working with graphs:

In [ ]: 
using Plots;
gr();

GMS Parameters:

In [ ]: 
# номинальное напряжение, В
Us = 500e3
# мощность КЗ, ВА
Skz = 10e9
# полное сопротивление системы, Ом
Zkz = Us ^ 2 / Skz
# соотношение активного и индуктивного сопротивления
XR = 15
# индуктивное сопротивление системы, Ом
Xkz = Zkz * sin(atan(XR))
# активное сопротивление системы, Ом
Rkz = Zkz * cos(atan(XR))
# комплексное сопротивление системы, Ом
Zs = Rkz + 1im * Xkz;

Transmission Line Parameters:

In [ ]: 
# удельное сопротивление прямой последовательности, Ом/км
Z0 = 0.01967 + 0.2899im
# количество цепей
n = 2
# длина ВЛ, км
L = 200
# удельная емкостная проводимость, См/км
b0 = 3.908e-6
bc = b0 * L * n
# полное сопротивление ВЛ, Ом
Zl = Z0 * L / n
# емкостное сопротивление половины ВЛ, Ом
Zc = -2im / bc;

Transformer parameters:

In [ ]: 
# активное сопротивление ТР, Ом
Rt = 0.4514
# индуктивное сопротивление ТР, Ом
Xt = 30.62
Zt = Rt + Xt * 1im;

Generator parameters:

In [ ]: 
# номинальное напряжение, В
Ug = 24e3
# номинальная мощность, ВА
Sg = 500e6
# активное сопротивление, о.е.
Ra = 0.003
# синхронное индуктивное сопротивление по оси d, о.е.
Xd = 1.81
# переходное индуктивное сопротивление по оси d, о.е.
Xdd = 0.3
# комплексное переходное сопротивление генератора, Ом
Zg = (Ra + Xdd * 1im) * Us ^ 2 / Sg;

Load Parameters:

In [ ]: 
# активная мощность нагрузки, Вт
Pn = 20e6
# активное сопротивление нагрузки, Ом
Zn = Us ^ 2 / Pn;

Calculation of intrinsic and mutual resistances $Z_{11}$ and $Z_{12}$ in normal mode by the unit current method:

In [ ]: 
I2 = 1
Uc = I2 * Zs
Ic0 = Uc / Zc
Ibc = I2 + Ic0
Ub = Uc + Zl * Ibc
Ib0 = Ub / Zc
Iab = Ib0 + Ibc
Ua = Ub + Iab * Zt
Ia0 = Ua / Zn
I1 = Ia0 + Iab
U1 = Ua + I1 * (Xdd * 1im * Us ^ 2 / Sg)
Z12d = (Ua + I1 * (Xd * 1im * Us ^ 2 / Sg)) / I2
Z12 = U1 / I2
Z11 = U1 / I1
alpha11 = pi / 2 - angle(Z11)
alpha12 = pi / 2 - angle(Z12);

Calculation of intrinsic and mutual resistances $Z_{11}$ and $Z_{12}$ in post-emergency mode by the single current method:

In [ ]: 
I2 = 1
Uc = I2 * Zs
Ic0 = Uc / (2 * Zc)
Ibc = I2 + Ic0
Ub = Uc + 2 * Zl * Ibc
Ib0 = Ub / (2 * Zc)
Iab = Ib0 + Ibc
Ua = Ub + Iab * Zt
Ia0 = Ua / Zn
I1 = Ia0 + Iab
U1 = Ua + I1 * (Xdd * 1im * Us ^ 2 / Sg)
Z12_pav = U1 / I2
Z11_pav = U1 / I1
alpha11_pav = pi / 2 - angle(Z11_pav)
alpha12_pav = pi / 2 - angle(Z12_pav);

The angular characteristics are calculated using the transient EMF $E`_q$, which remains unchanged at the initial moment of short-circuit time and after its tripping:

In [ ]: 
# вырабатываемая генератором мощность в нормальном режиме
P0 = 400e6;
Q0 = 0;
# синхронная ЭДС по оси q
Eq = sqrt((Us + Q0 * abs(Z12d) / Us) ^ 2 + (P0 * abs(Z12d) / Us) ^ 2)
# угол ротора в исходном режиме
d0 = atan(P0 * abs(Z12d) / (Us ^ 2 + Q0 * abs(Z12d)))
# переходная ЭДС
Ep = sqrt((Us + Q0 * abs(Z12) / Us) ^ 2 + (P0 * abs(Z12) / Us) ^ 2)
# угол переходной ЭДС
d0p = atan(P0 * abs(Z12) / (Us ^ 2 + Q0 * abs(Z12)))
# переходная ЭДС по оси q в нормальном режиме
Eqq = Ep * cos(d0p - d0)
# переходная ЭДС по оси q в аварийном и послеаварийном режимах
Eqq_pav = Eqq
Eqq_av = Eqq
# максимальная активная мощность в нормальном режиме
Pmax = (Eqq ^ 2 / abs(Z11) * sin(alpha11) + Eqq * Us / abs(Z12)) / 1e6
# угол переходной ЭДС по оси q
d0p = asin(P0 / Pmax / 1e6)
# максимальная мощность в аварийном режиме
Pmax_av = 0
# максимальная мощность в послеаварийном режиме
Pmax_pav = (Eqq_pav ^ 2 / abs(Z11_pav) * sin(alpha11_pav) + Eqq_pav * Us / abs(Z12_pav)) / 1e6
# угол ротора при работе на послеаварийной характеристике
d0_pav = asin(P0 / Pmax_pav / 1e6)
# критический угол
dkr = pi - d0_pav
# предельный угол отключения
d_pr = acos((P0 * 1e-6 * (pi - d0_pav - d0p) + Pmax_pav * cos(pi - d0_pav) - Pmax_av * cos(d0p)) / (Pmax_pav - Pmax_av));

Visualisation of the area method:

In [ ]: 
# функции угловых характеристик
delta = 0:0.1:180
Pmax = x -> (Eqq ^ 2 / abs(Z11) * sin(alpha11) .+ (Eqq * Us / abs(Z12)) * sin.(x .- alpha12)) / 1e6
Pmax_av = x -> fill(0, size(x))
Pmax_pav = x -> (Eqq_pav ^ 2 / abs(Z11_pav) * sin(alpha11_pav) .+ Eqq_pav * Us / abs(Z12_pav) * sin.(x .- alpha12_pav)) / 1e6 
Pt = x -> fill(P0 * 1e-6, size(x))
# визуализация метода площадей
plot(delta, [Pmax(delta .* pi ./ 180) Pmax_av(delta .* pi ./ 180) Pmax_pav(delta .* pi ./ 180) Pt(delta)],
label = [L"Норм.реж." L"Авар.реж." L"П/АВ.реж." L"Мощн.турбины"], ylabel = L"P, МВт", xlabel = L"\delta\prime, град", legend=false,
xaxis = (0:30:180), size = (900,440), left_margin=5Plots.mm, bottom_margin=5Plots.mm)
plot!((d0p:0.01:d_pr) .* 180 / pi, Pmax_av((d0p:0.01:d_pr)), fillrange = Pt((d0p:0.01:d_pr)), fillalpha = 0.2, c = 2, label = "Fторм")
plot!((d_pr:0.01:dkr) .* 180 / pi, Pmax_pav((d_pr:0.01:dkr)), fillrange = Pt((d_pr:0.01:dkr)), fillalpha = 0.2, c = 1, label = "Fторм")
vline!([d0p,d_pr,dkr] .* 180 / pi, c = 5)
annotate!([(90, Pmax(pi / 2) - 50, (L"Норм.реж.", 10, :black)),
        (50, Pmax_av(pi / 2) .+ 30, (L"Авар.реж.", 10, :black)),
        (90, Pmax_pav(pi / 2) - 50, (L"ПА/В.реж.", 10, :black)),
        (d0p * 180 / pi + 20, P0 * 1e-6 - 200, (L"F_{уск}", 15, :black, :left)),
        (d_pr * 180 / pi + 15, P0 * 1e-6 + 200, (L"F_{торм}", 15, :black, :left)),
        (0, P0 * 1e-6 + 50, (L"P_{турб}", 10, :black, :left)),
        (d0p * 180 / pi, 30, (L"\delta\prime_0  = %$(round(d0p * 180 / pi, digits=2))", 12, :black, :right)),
        (d_pr * 180 / pi, 30, (L"\delta\prime_{пред} = %$(round(d_pr * 180 / pi, digits=2))", 12, :black, :left)),
        (dkr * 180 / pi, 30, (L"\delta\prime_{крит} = %$(round(dkr * 180 / pi, digits=2))", 12, :black, :right))])
Out[0]:

Limit time of short-circuit tripping:

In [ ]: 
Tj = 7.05
tpr = sqrt((d_pr - d0p) * 2 * Tj / (0.8 * 2 * pi * 50))
print("t_пред = $(round(tpr, digits = 3))")

t_пред = 0.235

Modelling

Loading the model:

In [ ]: 
model_name = "synchronous_generator_round_rotor";
model_name in [m.name for m in engee.get_all_models()] ? engee.open(model_name) : engee.load( "$(@__DIR__)/$(model_name).engee");

Setting the model by command control for a line trip and short-circuit time equal to 0.24 c at which dynamic stability will be preserved:

In [ ]: 
step = "Отключение КЗ";
# отключение цепи ВЛ
engee.set_param!(model_name * "/" * step, "Time" => 5.24);

Running the loaded model:

In [ ]: 
results = engee.run(model_name);

To import the simulation results, logging of the required signals has been enabled in advance and their names have been set. Let's convert the simulation results from the results variable:

In [ ]: 
# вектор времени симуляции
sim_time = results["w"].time;
# вектор скорости вращения ротора
w1 = results["w"].value;
# вектор тока генератора
P1 = results["P"].value;

Simulation when increasing the short-circuit duration up to 0.2402 s:

In [ ]: 
# отключение цепи ВЛ
engee.set_param!(model_name*"/"*step, "Time" => 5.2402);
results = engee.run(model_name);
# вектор скорости вращения ротора
w2 = results["w"].value;
# вектор тока генератора
P2 = results["P"].value;

Generator active power plots:

In [ ]: 
plot(sim_time, P1, title = L"Мощность\; генератора", ylabel = L"P, о.е.", xlabel="Время, c",  label = L"t_{КЗ} = 0.24 с", legendfontsize = 10)
plot!(sim_time, P2,  label = L"t_{КЗ} = 0.2402 с", ylims = (-1, 2),size = (900,440), left_margin=5Plots.mm, bottom_margin=5Plots.mm)
Out[0]:

Rotor speed graphs:

In [ ]: 
plot(sim_time, w1, title = L"Скорость\; ротора", ylabel = L"\omega, о.е.", xlabel=L"Время, c", label = L"t_{КЗ} = 0.24 с", legendfontsize = 10)
plot!(sim_time, w2, label = L"t_{КЗ} = 0.2402 с", ylims = (0.95, 1.1),size = (900,440), left_margin=5Plots.mm, bottom_margin=5Plots.mm)
Out[0]:

In the second experiment, the dynamic stability of the generator was violated, thus the limit time of short-circuit tripping in the model was found empirically. The result converges with the theoretical calculation with a small error.

In [ ]: 
println("Предельное время отключения аналитически = "*string(round(tpr, digits=3)) * "с")
println("Предельное время отключения практически = "*string(0.24) * "с")

Предельное время отключения аналитически = 0.235с
Предельное время отключения практически = 0.24с

Addendum

Try to change the following model parameters yourself and investigate how this affects the dynamic stability:

  1. increase the overhead line length by 200 km;
  2. the type of short circuit in the Fault (Three-Phase) block;
  3. reduce the active power of the generator in the initial mode by 200 MW.

Conclusions

In this example the tools for command control of the Engee model and uploading of the simulation results were used and the work with the Plots module was shown. The operation of an implicitly pole generator through a transformer and a double-circuit overhead line on a SWBM without excitation and speed controllers has been considered. The analytical calculation of dynamic stability has been carried out and the results have been compared with the simulation data.

References

  1. Venikov V. A. Transient electromechanical processes in electrical systems: Textbook for electric power engineering specialities of universities. - 4th ed., revision and supplement. - M.: Vysh. shk., 1985.
  2. Transients of electrical systems in examples and illustrations: Textbook for universities (V. V. Ezhkov, N. I. Zelenohat, I. V. Litkens, etc.; Edited by V. A. Stroev). - Moscow: Znak, 1996. - 224 pp., ill.
  3. Dolgov, A.P. Stability of electrical systems : a textbook / Dolgov A.P.. - Novosibirsk : Novosibirsk State Technical University, 2010. - 177 c. - ISBN 978-5-7782-1320-3.

Blocks used in example