Годографы сопротивления
Введение
В цифровой релейной защите, чтобы посчитать действующее значение тока или угол между током и напряжением прибегают к векторному представлению синусоидальных величин. Для этого используют фильтр Фурье (пара КИХ фильтров длиной один период синусоидального сигнала со сдвигом в ФЧХ на 90 градусов).
Такие фильтры прексрасно работают в установившемся режиме. При изменении режима в ЭС (КЗ, отключение и т.д.) в них, как в любых фильтрах, возникает переходной процесс пока в линии задержки фильтра всё еще находятся отсчеты предаварийного режима.
Переходной процесс фильтра
Рассмотрим идеализированный (без апериодической составляющей) процесс увеличения тока. Вектор, характеризующий ток, при КЗ увеличивается и поворачивается. Длительность этого процесса - 1 период.
Вопрос: по какой траектории происходит движение вектора?

Движение по траектории 1 выглядит хорошо для защиты, которая срабатывает по величине тока, потому что на протяжении переходного процесса в фильтре значения не становятся больше большего из векторов. Траектория 2 в этом плане хуже. Защита без выдержки времени может пуститься и сработать на остриях траектории, так как значение тока (длины вектора) в этих точках переходной траектории могут оказаться выше уставки.
С другой стороны переходная траектория не отражает никакой режим в ЭС и в некоторых источниках предлагается всю переходную траекторию помечать как недостоверные измерения, которые не влият на пусковые органы защиты. Это тоже выход, но даёт стабильную задержку в реакции на повреждение.
Фильтр Фурье даёт переходную траекторию в виде циклоиды, как будто диск катится от одного вектора к другому и делает два оборота.

Дистанционная защита
В дистанционной защите, которая контролирует сопротивление (с учетом его активной и реактивной составляющих), понимание траекторий становися более важным. Годограф вектора сопротивления не должен попадать в область, ограниченную реле сопротивлением при внешних КЗ и должен там стабилизироваться при КЗ в зоне защиты.
Например, при КЗ в зоне защиты (спереди) желтая траектория входа нас полностью устроит, а при КЗ "за спиной" зелёная траектория нас не устроит никак. Первые ступени дистанционной защиты должны работать без выдержки времени, поэтому ждать прохода годографа через зону реле сопротивления не очень хорошо.

Построение годографов
При оценка комплекса сопротивления через деление вектора напряжения на вектор тока мы получим в комплексном виде деление одной траектории в виде циклоиды на другую. Причем вид циклоид (их начальная фаза) зависит от фазы аварии. Посмотрим как ведут себя траектории.
Примечание
Про дистанционную защиту и переходные годографы есть хорошие и объемные материалы у Шнеерсона Э.М. Данный демонстрационный проект основан на них.
using LinearAlgebra
# @markdown Функции и константы
const n = 80;
const T = 0.02 / n;
const w = 100pi;
function calc_ui(Z_load, Z_sc, Zsyst, direction, withAperiodic)
if direction
z_ratio = ( Z_sc + Zsyst ) / ( Z_load + Zsyst );
else
z_ratio = ( Z_sc - Zsyst ) / ( Z_load - Zsyst );
end
cols = Vector{Vector{Float64}}()
for phi = 0: 20*pi/180: pi
I_load = exp(im * phi);
I_sc = I_load / z_ratio;
U_load = I_load * Z_load;
U_sc = I_sc * Z_sc;
abs_U_load = abs(U_load);
arg_U_load = angle(U_load);
abs_I_sc = abs(I_sc);
arg_I_sc = angle(I_sc);
abs_U_sc = abs(U_sc);
arg_U_sc = angle(U_sc);
Nm = 6*n
i_vals = Vector{Float64}(undef, Nm)
u_vals = Vector{Float64}(undef, Nm)
wT = w * T
for k = 0: Nm-1
if k <= 2*n-1
i_vals[k+1] = sin(wT * k + phi)
u_vals[k+1] = abs_U_load * sin(wT * k + arg_U_load)
else
if withAperiodic
iexp = -(abs_I_sc * sin(wT * 2*n + arg_I_sc) - sin(wT * 2*n + phi)) * exp(-T*(k - 2*n)/0.02)
i_vals[k+1] = abs_I_sc * sin(wT * k + arg_I_sc) + iexp
else
i_vals[k+1] = abs_I_sc * sin(wT * k + arg_I_sc)
end
u_vals[k+1] = abs_U_sc * sin(wT * k + arg_U_sc)
end
end
push!(cols, u_vals)
push!(cols, i_vals)
end
return reduce(hcat, cols)
end
function calc_impedance(ui_data, expInvariant)
phi = 2*pi/n
coeff = 2/n .* exp.(-(0: n-1) .* phi .* im)
function justFourierFilter(vector)
return dot(vector, coeff)
end
function expInvariantFilter(vector)
minq = exp(-1/50/n/0.005)
constLevel = 2
s1 = sum(vector[end-n+1: end])
s2 = sum(vector[1: n])
if abs(s1) > constLevel
q = clamp(abs(s1 / s2)^0.05, minq, 1.0)
expv = [q^j for j in 0:n-1]
coeffMod = coeff .- sum(coeff .* expv) / sum(expv) # Формируем сдвинутые коэффициенты фильтров для подавления оцененной выше апериодики
return dot(vector[end-n+1: end], coeffMod)
else
return dot(vector[end-n+1: end], coeff)
end
end
s = size(ui_data)
cols = Vector{Vector{ComplexF64}}()
n4 = div(n, 4)
for v = 1: 2: s[2]
z = ComplexF64[]
for k = n + n4: s[1]
if expInvariant
v_i = expInvariantFilter(ui_data[k-n-n4+1: k, v+1])
v_u = expInvariantFilter(ui_data[k-n-n4+1: k, v])
else
v_i = justFourierFilter(ui_data[k-n+1: k, v+1])
v_u = justFourierFilter(ui_data[k-n+1: k, v])
end
push!(z, v_u/v_i)
end
push!(cols, z)
end
return reduce(hcat, cols)
end;
Сформируем значения токов и напряжений для различных начальных фаз аварии и посчитаем сопротивления.
# ИСХОДНЫЕ ДАННЫЕ
Z_load = 9 + 1im; #кажущееся сопротивление в режиме нагрузки
Z_sc = 0.8*exp(im*80*pi/180); #кажущееся сопротивление в режиме КЗ спереди
Z_sc_b = 0.5 * Z_sc; #кажущееся сопротивление в режиме КЗ за спиной
Z_back = 3im; #сопротивление системы за спиной
Z_forward = 2.1im; #сопротивление системы спереди
Опыт 1. КЗ в зоне и за спиной при синусоидальных воздействиях
data_f_01 = calc_ui(Z_load, Z_sc, Z_back, true, false)
data_b_01 = calc_ui(Z_load, -Z_sc_b, Z_forward, false, false)
z_f_01 = calc_impedance(data_f_01, false);
z_b_01 = calc_impedance(data_b_01, false);
gr()
anim = @animate for i = 60: 140
plt1 = plot(title="Переход замера Z от нагрузки к КЗ \n (в зоне защиты и за спиной)",
framestyle = :zerolines,
ylabel="jX",
yguidefontsize = 16,
xlabel="R",
xguidefontsize = 16,
xlims=(-5, 15),
ylims=(-10, 10),
linecolor = "#7E2F8E",
legend = false,
aspect_ratio = :equal);
plot!(plt1, real.(z_f_01[1: i, :]), imag.(z_f_01[1: i, :]), linecolor = "#D95319")
plot!(plt1, real.(z_b_01[1: i, :]), imag.(z_b_01[1: i, :]), linecolor = "#7E2F8E")
plot!(plt1,[-0.1, 1, 1.2, -0.15, -0.1], [-0.1, -0.15, 2, 2.1, -0.1], linewidth = 2, color = :black);
scatter!(plt1, [Z_load.re], [Z_load.im], markersize = 3, color = :red);
scatter!(plt1, [Z_sc.re], [Z_sc.im], markersize = 3, color = :red);
scatter!(plt1, [-Z_sc_b.re], [-Z_sc_b.im], markersize = 3, color = :red);
end
gif(anim, "v1.gif", fps=10)
Опасения были излишни, при КЗ за спиной (фиолетовый цвет) годографы подходят к точке НЕ через зону срабатывания. При КЗ спереди коричневые годографы заходят в зону срабатывания плотным пучком.
Опыт 2. КЗ в зоне и за спиной при учёте свободных составляющих в токах
Пересечение всех годографов в одной точке (Опыт 1) связано с чисто синусоидальными воздействиями и не имеет места при наличии в токах апериодических составляющих. Почситаем тот же процесс но с учетом свободных составляющих.
data_f_02 = calc_ui(Z_load, Z_sc, Z_back, true, true)
data_b_02 = calc_ui(Z_load, -Z_sc_b, Z_forward, false, true)
z_f_02 = calc_impedance(data_f_02, false);
z_b_02 = calc_impedance(data_b_02, false);
gr()
anim = @animate for i = 60: 180
plt2 = plot(title="Переход замера Z от нагрузки к КЗ (в зоне защиты\n и за спиной) с учетом свободных составляющих в токах",
framestyle = :zerolines,
ylabel="jX",
yguidefontsize = 16,
xlabel="R",
xguidefontsize = 16,
xlims=(-5, 15),
ylims=(-10, 10),
linecolor = "#7E2F8E",
legend = false,
aspect_ratio = :equal);
plot!(plt2, real.(z_f_02[1: i, :]), imag.(z_f_02[1: i, :]), linecolor = "#D95319")
plot!(plt2, real.(z_b_02[1: i, :]), imag.(z_b_02[1: i, :]), linecolor = "#7E2F8E")
plot!(plt2,[-0.1, 1, 1.2, -0.15, -0.1], [-0.1, -0.15, 2, 2.1, -0.1], linewidth = 2, color = :black);
scatter!(plt2, [Z_load.re], [Z_load.im], markersize = 3, color = :red);
scatter!(plt2, [Z_sc.re], [Z_sc.im], markersize = 3, color = :red);
scatter!(plt2, [-Z_sc_b.re], [-Z_sc_b.im], markersize = 3, color = :red);
end
gif(anim, "v2.gif", fps=10)
Траектории изменились довольно сильно, но не нарушили главного свойства: при КЗ "за спиной" пересечений с зоной срабатывания нет.
Опыт 3. КЗ спереди вне зоны, с учетом свободных составлющих (фильтр Фурье)
Посмотрим поближе, что происходит с траекториями если КЗ будет чуть дальше, чем мы ограничили нашу зону. То есть всё еще спереди, но реагировать на это КЗ уже не надо.
Z_sc = 2.4*exp(im*70*pi/180);
Z_back = 1.5im;
data_f_03 = calc_ui(Z_load, Z_sc, Z_back, true, true)
z_f_03 = calc_impedance(data_f_03, false);
plt3 = plot(title="Переход замера Z от нагрузки к КЗ \n вне зоны действия",
framestyle = :zerolines,
ylabel="jX",
yguidefontsize = 16,
xlabel="R",
xguidefontsize = 16,
xlims=(-1, 3),
ylims=(-1, 3),
legend = false,
aspect_ratio = :equal);
plot!(plt3, real.(z_f_03), imag.(z_f_03), linecolor = "#D95319")
plot!(plt3,[-0.1, 1, 1.2, -0.15, -0.1], [-0.1, -0.15, 2, 2.1, -0.1], linewidth = 2, color = :black);
scatter!(plt3, [Z_load.re], [Z_load.im], markersize = 3, color = :red);
scatter!(plt3, [Z_sc.re], [Z_sc.im], markersize = 3, color = :red);
scatter!(plt3, [-Z_sc_b.re], [-Z_sc_b.im], markersize = 3, color = :red);
display(plt3)
Видно, что годографы закручиваются улиткой вокруг истинного значения сопротивления. Эта улитка объясняется мешающим влиянием экспоненты на фильтр Фурье. Также видно, что эта улитка проходит через зону срабатывания дистанционной защиты, что плохо для защиты без выдержки времени.
Опыт 4. КЗ спереди вне зоны, с учетом свободных составлющих (модифицированный фильтр Фурье)
Как я писал в своём первом демонстрационном проекте коэффициенты фильтра Фурье можно всегда сместить на некоторую константу, чтобы снизить влияние экспоненты с конкретным значением затухания. Если это затухание оценить независимо, то в реальном времени можно подстраиваться под любую апериодику. Посмотрим на результаты. Для этого используем флаг expInvariant функции calc_impedance.
data_f_04 = calc_ui(Z_load, Z_sc, Z_back, true, true)
z_f_04 = calc_impedance(data_f_04, false);
z_f_05 = calc_impedance(data_f_04, true);
plt4 = plot(title="Переход замера Z от нагрузки к КЗ \n вне зоны действия (Фурье и модификация)",
framestyle = :zerolines,
ylabel="jX",
yguidefontsize = 16,
xlabel="R",
xguidefontsize = 16,
xlims=(0, 3),
ylims=(1, 4),
legend = false,
aspect_ratio = :equal);
plot!(plt4, real.(z_f_04), imag.(z_f_04), linecolor = "#D95319")
plot!(plt4, real.(z_f_05), imag.(z_f_05), linecolor = :green)
plot!(plt4,[-0.1, 1, 1.2, -0.15, -0.1], [-0.1, -0.15, 2, 2.1, -0.1], linewidth = 2, color = :black);
scatter!(plt4, [Z_load.re], [Z_load.im], markersize = 3, color = :red);
scatter!(plt4, [Z_sc.re], [Z_sc.im], markersize = 3, color = :red);
scatter!(plt4, [-Z_sc_b.re], [-Z_sc_b.im], markersize = 3, color = :red);
display(plt4)
Я бы сказал, что стало лучше в части того, что более нет улитки, но опасная близость с зоной срабатывания всё-таки остается. Поэтому такое решение нельзя считать универсальным.
Вместо выводов
- Анимации переходных годографов выглядят красиво;
- На графиках и анимациях можно увидеть, опасные сближения годографов с зоной срабатывания и принять дополнительные меры по улучшению селективных свойств защиты (в идеале без ущерба для быстродействия);
- Можно попробовать при разработке новых способов расчета интегральных параметров сразу учитыавть требования к переходным годографам.



