Реализация цифровых фильтров‌‌ в микропроцессорных устройствах релейной защиты

DOI: 10.21122/1029-7448-2016-59-5-397-417 УДК 621.316.925

Реализация цифровых фильтров в микропроцессорных устройствах релейной защиты

Ю. В. Румянцев1*, Ф. А. Романюк1*, В. Ю. Румянцев1*, И. В. Новаш1*

^Белорусский национальный технический университет

© Белорусский национальный технический университет, 2016 Бе1аги81ап National Technical University, 2016

Реферат. Рассмотрена реализация цифровых фильтров, используемых в измерительных органах тока микропроцессорных устройств релейной защиты. Приведены математические описания данных фильтров, а также программы для расчета их коэффициентов. Показано, что достоверно оценить функционирование цифровых фильтров можно только при входных воздействиях, приближенных по виду к реальным формам вторичного тока трансформатора тока, к которому подключена защита, использующая оцениваемый фильтр. Для этих целей в среде динамического моделирования MatLab-Simulink воспроизведена структурная модель энергосистемы и измерительного органа тока микропроцессорной (цифровой) защиты. Расчеты позволили установить, что при наличии в первичном сигнале тока экспоненциально затухающей апериодической составляющей магнитопровод трансформатора тока входит в насыщение даже при допустимой кратности первичного тока и номинальной вторичной нагрузке. Следствием этого являются искажение вторичного тока трансформатора тока и появление в нем высших гармоник, а также гармоник, не кратных основной. Сравнение реакций фильтров на указанные воздействия позволило выявить такие реализации фильтров, которые дают возможность определять амплитуду входного сигнала с наименьшими погрешностями. Предложены пути повышения качества фильтрации, связанные с применением оконных функций. Так, совместное использование цифровых фильтров и окна Хэм-минга позволяет достичь близкого к нулю значения коэффициента усиления амплитуды сигнала в области высоких частот и практически подавить все спектральные составляющие, начиная с частоты 100 Гц. Происходящее при этом увеличение времени установления амплитуды сигнала может быть уменьшено за счет выбора наиболее оптимальных параметров всех элементов измерительного органа тока микропроцессорной защиты.

Ключевые слова: цифровой фильтр, дискретное преобразование Фурье, метод наименьших квадратов, ортогональные составляющие, окно Хэмминга, цифровая обработка сигналов, трансформатор тока, насыщение, микропроцессорная релейная защита, MatLab, Simulink

Для цитирования: Реализация цифровых фильтров в микропроцессорных устройствах релейной защиты / Ю. В. Румянцев [и др.] // Энергетика. Изв. высш. учеб. заведений и энерг. объединений СНГ. 2016. Т. 59, № 5. C. 397-417

Адрес для переписки

Романюк Федор Алексеевич Белорусский национальный технический просп. Независимости, 65/2, 220013, г. Минск, Республика Беларусь Тел.: +375 17 331-00-51 faromanuk@bntu.by

Address for correspondence

Romaniuk Fiodar A.

Belarusian National Technical University 65/2 Nezavisimosty Ave., 220013, Minsk, Republic of Belarus Tel.: +375 17 331-00-51 faromanuk@bntu.by

Digital Filters Implementation in Microprocessor-Based Relay Protection

Yu. V. Rumiantsev1), F. A. Romaniuk1), V. Yu. Rumiantsev1), I. V. Novash1)

'-Belarusian National Technical University

Abstract. This article presents the implementation of digital filters used in digital relay protection current measuring elements. Mathematical descriptions of the considered digital filters as well as the computer programs for their coefficients calculation are described. It has been shown that in order to reliable estimate the digital filter performance its input signals waveforms must be close to the actual secondary current waveform of the current transformer to which the digital protection with the estimated digital filter is connected. For these purposes in MatLab-Simulink dynamic simulation environment the power system and the current measuring element models were developed. Performed calculations allowed to reveal that the exponentially decaying DC component which in some cases contains in primary fault current drives the current transformer core into saturation even when its nominal parameters are not exceeded. This results in distortion of the current transformer secondary current which in this case contains higher and inter-harmonics. Moreover, such harmonic content is not completely taking into account during coefficients calculation of the considered digital filters what results in signal magnitude estimation inaccuracy. Comparison of the digital filters response to the above-mentioned input signals allowed to find out such digital filter implementations which enable signal magnitude estimation with a minimum error. Ways of filtering quality improvement concerned with the window functions are proposed. Thus, the joint usage of digital filter and Hamming window allows to achieve the zero value of the signal magnitude gain factor in high-frequency range and substantially suppress all spectral components above 100 Hz. The increasing of the signal magnitude settling time in this case can be reduced by choosing the most optimal parameters of the all components of the current measuring element.

Keywords: digital filter, discrete Fourier transform, least error squares technique, orthogonal component, Hamming window, digital signal processing, current transformer, saturation, digital relay protection, MatLab, Simulink

For citation: Rumiantsev Yu. V., Romaniuk F. A., Rumiantsev V. Yu., Novash I. V.

Введение

Постоянное совершенствование цифровых устройств релейной защиты (РЗ) позволяет реализовывать все более сложные алгоритмы обработки контролируемых сигналов [1, 2]. При этом увеличивается частота дискретизации входных сигналов, достигая 24—64 выборок на период промышленной частоты, что, однако, не приводит к существенному повышению качества определения амплитуды сигнала цифровыми фильтрами (ЦФ). Это связано с тем, что параметры ЦФ рассчитываются применительно к простейшим входным воздействиям, форма которых существенно отличается от формы реальных аварийных сигналов.

Основные повреждения в энергосистеме - короткие замыкания (КЗ), токи которых не являются периодическими сигналами из-за наличия в них экспоненциально затухающей апериодической составляющей. В свою очередь, информация о токах повреждения передается к устройствам РЗ через электромагнитные измерительные трансформаторы тока (ТТ), магнитопро-

воды которых входят в насыщение от данной составляющей тока КЗ, в результате чего происходит искажение формы вторичного тока ТТ. Поэтому для проверки работоспособности ЦФ в условиях, приближенных к реальным, в качестве тестового входного воздействия целесообразно использовать вторичный ток ТТ [3].

Целью настоящей работы являлось исследование работоспособности различных типов ЦФ и выбор наиболее приемлемой их реализации. Для этого разработали алгоритмы расчета коэффициентов ЦФ, а в среде MatLab-Simulink [4] создали имитационную модель для возможности получения искаженных форм сигнала вторичного тока ТТ.

Основная часть

В отличие от аналоговых, в цифровых измерительных органах устройств РЗ все операции производятся над цифровыми отсчетами

T = где f - частота основной гармоники, Гц.

Тогда аналоговому сигналу u (t) = sin(0t) соответствует цифровой сигнал u(nT) = sin| 2nn|, где n = 1, 2, ..., N - целое число, соответствую-

l N )

щее номеру выборки. Кроме того, цифровой сигнал принято представлять в виде

U (nT) = Um Sin С2Nn + ф] = Um Sin (^Nn) C0S9 ++ Um C0S ^) ^ =

T TC . , 2nn l TTS

= UC sin I-1 + Ul cosí

N ) У N

где иьт = ит 8Шф, ист = ит ео8ф - амплитуда синусной и косинусной ортогональных составляющих сигнала (1).

Входной сигнал цифровой защиты у

у

Назначение ЦФ заключается в выделении полезного сигнала и(пТ) (обычно это основная гармоника, реже - высшая, кратная основной) и в максимально возможном подавлении помехи е(пТ). В общем случае сигнал помехи содержит: высшие гармоники, кратные основной частоте; затухающие апериодические составляющие; затухающие и незатухающие гармонические составляющие с частотой, отличной от основной [5]. Поэтому точный учет всех составляющих помехи принципиально невозможен,

а следовательно, ЦФ всегда будет выделять полезный сигнал с ошибкой тем большей, чем больше значения выборок входного сигнала y(nT) будут отличаться от полезного сигнала u(nT).

Для выделения полезного сигнала широко используется метод наименьших квадратов (МНК), в котором минимизируется сумма квадратов ошибок в соответствии с выражением

N

S = ^ (y(nT) - u(nT))2 ^ min.

На практике сигнал помехи обычно представляют упрощенным и учитывают в нем только основные составляющие: экспоненциально затухающую апериодическую составляющую и несколько высших гармоник, кратных основной. В зависимости от набора компонент в помехе формируются различные алгоритмы построения ЦФ, такие как: ЦФ на основе МНК [6], на основе дискретного преобразования Фурье (ДПФ) [7], а также формирователи ортогональных составляющих (ФОС) [8] и др. Рассмотрим более подробно построение перечисленных выше ЦФ.

Метод наименьших квадратов

Идентификация входного сигнала РЗ на основе МНК впервые предложена в [6]. Сигнал представляется состоящим из полезного сигнала основной гармоники и помехи в виде экспоненциально затухающей апериодической составляющей с постоянной времени т и 3-й гармоники

y

Выражение (3) отражает в первом приближении сигнал тока при КЗ в энергосистеме. Состав компонентов в нем объясняется тем, что четные гармоники не содержатся в токе повреждения, а гармоники, начиная с 5-й, существенно ослабляются аналоговым фильтром нижних частот (ФНЧ). Необходимо отметить, что еще на стадии обсуждения данной публикации авторам было указано на некоторую идеализацию модели сигнала, однако сама идея использования МНК в дальнейшем получила широкое распространение.

Экспоненциально затухающий сигнал с постоянной времени т разлагают в ряд Тейлора

nT rri 2^2 Зт-гЗ

— , nT n T nT

e T = 1--^---

2t2 6T3

и используют первые три-четыре члена ряда, а гармонические сигналы представляют согласно (1). В результате для любого номера выборки п = 1, 2, ..., N выражение (3) принимает вид

У(пТ) = Ua

Гл nT n T2 !

1--+-2

т 2t2

um isin

2nn

N

U'm 1 Cos

2nn

N

N

TTC . I 6nn ) TTS (6nn um asm I — 1 + um 3Cos

N

Обозначив X, = Ua; X2 =-Ua-; X3 = X4 = U!; X = US,; X = UC

т 2т

ml ' 5 ml

Х7 = иьт3 как неизвестные, можно записать (4) для всех номеров выборок п = 1, 2, ..., N в матричном виде У = АХ:

У

У

2п

N

2п

1 Т Т sin I — I cos I — I sin I — I cos I —

N

6n

N

1 TT /1Т2 • I 4"! ( 4n | . (12n

1 2T 4T sin I — I cos I — I sin I-I cos,

N J I N J I N J I N

i it qt,2 • I (8n! . (18n! (18n

1 3T 9T sin — cos — sin - cos

N

N

N

6n

N 12n

N

2^2 • I 2nn ! ( 2nn! . (2nn ! ( 2nn

1 nT n T sin - cos - sin - cos

N

N

N

N

X1 X 2 X3 X4 X5

X6 X

где А - матрица коэффициентов (главная матрица системы); X - вектор-столбец неизвестных; У - вектор-столбец свободных членов.

Квадратную систему линейных уравнений в матричном виде можно решить, умножив вектор-столбец свободных членов У на обратную матрицу А-1

X = А-1У.

На практике число выборок на период основной частоты для современных защит N = 24, 32, ..., т. е. число уравнений N оказывается большим числа неизвестных. Такие системы уравнений называются переопределенными, так как в них матрица коэффициентов А - прямоугольная с числом строк, большим числа столбцов. Для решения таких систем используют псевдообращенную матрицу А+ [9]

A + = (Aт A)-1 Ат,

где Ат - транспонированная матрица.

Соответственно приближенное решение системы уравнений

X = А+Y.

Получение псевдообращенной матрицы является стандартной процедурой для большинства математических пакетов, например в системе MatLab - это функция p=pinv (a) [4]. Для рассматриваемого примера файл-программа на языке программирования MatLab для вычисления псевдообращенной матрицы:

f=50;% Промышленная частота, Гц N=2 4;% Число выборок за период T=1/(N*f);% Период дискретизации k=3;% Номер высшей гармонипки for n=0:N-1; a(n+1,1)=1; a(n+1,2)=n*T; a(n+1,3)=TA2*nA2; a(n+1,4)=sin(2*pi*n/N); a(n+1,5)=cos(2*pi*n/N); a(n+1,6)=sin(2*pi*k*n/N); a(n+1,7)=cos(2*pi*k*n/N); end;

x=pinv(a);% Получение псевдообращенной матрицы x cos=x(4,:);% Коэффициенты acn x sin=x(5,:);% Коэффициенты asn.

В результате матрица коэффициентов А[24х7] размером 24 строки на 7 столбцов преобразуется в матрицу А+[7х24], причем 4-я и 5-я строки этой матрицы представляют собой коэффициенты acn и asn цифровых фильтров для выделения косинусной и синусной ортогональных составляющих основной гармоники по выражениям:

= Z acny(nT);

n=1

N

,=Ё asny

sn^

n=1

Амплитуда сигнала основной гармоники для произвольной выборки п определяется как

U , =Ju2 + u2

mln \ sn cn

u

Дискретное преобразование Фурье

Невозможность точного учета всех аварийных составляющих в сигнале повреждения позволяет предложить минимально возможную по составу элементов модель сигнала, в которой содержится только основная гармоника

у(пт)=ит ^ [ N ]+и: ^ [ ^

В этом случае рассчитанные по аналогии с моделью (4) коэффициенты ЦФ полностью совпадают с коэффициентами при ортогональных составляющих дискретного преобразования Фурье

2 N—

i

N n=0

2nknЛ . . (2nkn

cos I-1 + j sin

N ) I N

где к - номер выделяемой гармоники.

Выражение (9) удобно представить в виде

N-1 N-1

и

п=0 п=0

где Ас(пТ), Ах(пТ) - коэффициент косинусного и синусного фильтров Фурье.

В практике построения алгоритмов функционирования зарубежных микропроцессорных токовых защит широко используются косинусные фильтры

Ас

Наследуя все преимущества ДПФ, косинусный фильтр является практически нечувствительным к апериодическим составляющим сигнала.

Формирователи ортогональных составляющих

В рассмотренных выше моделях сигналов количество уравнений превышало число неизвестных. Однако если зафиксировать на конечном интервале наблюдения (окно данных) длительностью порядка одного периода промышленной частоты четное число выборок сигнала Ь, то можно сформировать квадратную матрицу размером Ь х Ь и находить коэффициенты ЦФ по (6) путем определения обратной матрицы А- . В этом случае модель сигнала строится с учетом того, что экспоненциально затухающая апериодическая составляющая представляется двумя членами ряда Тейлора, и тогда на основную и высшие гармоники в каждой строке матрицы коэффициентов приходится М = Ь - 2 элемента [8]. Поскольку каждая гармоника

404 Реализация цифровых фильтров в микропроцессорных устройствах релейной защиты

сигнала записывается в виде двух ортогональных составляющих, для четного Ь максимальная гармоника в сигнале имеет номер М/2

Файл-программа для расчета коэффициентов ФОС на языке программирования Ма1;ЬаЬ с выводом результатов в рабочую область показан ниже. В отличие от вышеприведенной программы для определения обратной матрицы используется оператор x=inv(a) [4].

L=24;% Число выборок - четное

T=0.02/2 4,•% Время между двумя соседними выборками

N=0.02/^% Число выборок за период основной частоты

for n=0:L-1;% Формирование матрицы коэффициентов a(n+1,1)=1;% Апериодическая a(n+1,L)=n*T;% составляющая for k=2:2:L-2;

a(n+1,k)=sin(pi*n*k/N);% Гармоники от 1-й a(n+1,k+1)=cos(pi*n*k/N);% до (L-2)/2 end; end;

x=inv(a);% Вычисление обратной матрицы x cos=x(2,:);% Вывод коэффициентов фильтра x sin=x(3,:);% В рабочую область MATLAB.

Особенностью расчета коэффициентов ФОС является то, что число выборок за период основной частоты N может быть как целым (если NT = = 0,02 с), так и действительным (если NT Ф 0,02 с). Такая возможность выбора времени T между соседними отсчетами сигнала позволяет варьировать коэффициенты ЦФ для получения наиболее оптимальной амплитудно-частотной (АЧХ) и переходной характеристик [11]. В фильтрах, основанных на ДПФ, число выборок L всегда равно числу выборок за период основной частоты N и принимается целым и четным из ряда 16, 20, 24 и т. д.

Если не учитывать экспоненциально затухающую апериодическую составляющую в сигнале, т. е. положить в вышеприведенной программе члены a(n+1,1)=0 и a(n+1,L)=0, то полученные коэффициенты ФОС x cos и x sin будут совпадать с коэффициентами синусного и косинусного фильтров ДПФ.

Совершенствование элементной базы микропроцессорных устройств РЗ позволяет учитывать в модели сигнала все большее число высших гармоник, даже тех, влияние которых существенно ослабляется аналоговым ФНЧ. Все это отрицательно сказывается на точности вычисления коэффициентов ЦФ и, в конечном счете, на качестве фильтрации. Поэтому целе-

5 0 Гц

сообразно включать в модель сигнала гармоники не выше 7-8-й, а остальные принимать равными нулю.

Тестовые воздействия

Основные критерии, по которым оценивается качество фильтрации, -время установления истинного значения выходного сигнала ЦФ и его максимальная колебательность. Все рассмотренные выше ЦФ представлялись парой ортогональных составляющих. Из каждой пары выбирали один фильтр с наилучшими показателями качества фильтрации основной гармоники (со значениями коэффициентов для N = 24), приведенными в табл. 1.

Таблица 1

Коэффициенты цифровых фильтров Digital filters coefficients

№ выборки КФ ФОС МНК Окно Хэмминга

1 0,0833 0,0000 -0,1407 0,0800

2 0,0805 0,1638 -0,0690 0,0971

3 0,0722 -0,0112 -0,0129 0,1470

4 0,0589 0,1423 0,0146 0,2260

5 0,0417 -0,0417 0,0200 0,3284

6 0,0216 0,1049 0,0229 0,4464

7 0,0000 -0,0833 0,0420 0,5714

8 -0,0216 0,0618 0,0823 0,6940

9 -0,0417 -0,1250 0,1301 0,8053

10 -0,0589 0,0244 0,1603 0,8968

11 -0,0722 -0,1555 0,1511 0,9619

12 -0,0805 0,0028 0,0969 0,9957

13 -0,0833 -0,1667 0,0127 0,9957

14 -0,0805 0,0028 -0,0727 0,9619

15 -0,0722 -0,1555 -0,1315 0,8968

16 -0,0589 0,0244 -0,1505 0,8053

17 -0,0417 -0,1250 -0,1363 0,6940

18 -0,0216 0,0618 -0,1085 0,5714

19 0,0000 -0,0833 -0,0859 0,4464

20 0,0216 0,1049 -0,0733 0,3284

21 0,0417 -0,0417 -0,0571 0,2260

22 0,0589 0,1423 -0,0123 0,1470

23 0,0722 -0,0112 0,0830 0,0971

24 0,0805 0,1638 0,2346 0,0800

Вторую ортогональную составляющую определяли путем явного или вычисляемого сдвига выходного сигнала ЦФ на угол л/2. Явный сдвиг реализуется путем задержки сигнала основной гармоники на N/4 отсчетов N — четное число) и всегда приводит к увеличению времени установления выходного сигнала ЦФ.

Для вычисления второй ортогональной составляющей гармонического сигнала по известной первой можно использовать две выборки - текущую ип и предыдущую ип-1, зафиксированные через период дискретизации Т [8]. Тогда амплитуда сигнала определяется по выражению

k2 + u«-i- 2unun-1cos

Umkn

'2nkn^

I N )

sin

2nkn

N

где к - номер выделяемой гармоники.

К достоинствам использования (13) следует отнести минимальную задержку в один период дискретизации Т при определении амплитуды полезного сигнала, а также принципиальную невозможность обращения знаменателя в нуль. Однако при отклонении частоты контролируемого сигнала от номинальной определение его амплитуды при помощи ЦФ сопровождается погрешностями. Результат определения амплитуды синусоидального сигнала и (/) = Бт(2л/?) с частотой /= 48 Гц приведен на рис. 1. Амплитуда сигнала изменяется во времени с колебаниями двойной частоты и наибольшей погрешностью, не превышающей величины

Епах - |150 - 3/|,%.

Если по условиям функционирования ЦФ такие колебания амплитуды являются недопустимыми, то можно использовать большее число выборок сигнала, но при этом усложняется расчетное выражение и увеличивается время определения амплитуды. Поэтому в дальнейшем (13) будет применяться для всех фильтров, а сами ЦФ станут именоваться по названию алгоритма их реализации - соответственно МНК, КФ и ФОС.

1,0

U, В 0

0,6 0,4 0,2 0

____________________ ...........3l.........

6 % " /

^max 150 3f, %

1

1

0,01

0,02

0,03

0,04

0,05 t, с 0,06

Рис. 1. Определение амплитуды синусоидального сигнала частотой 48 Гц цифровым фильтром, настроенным на частоту 50 Гц

Fig. 1. 48 Hz sine signal magnitude estimation by digital filter tuned to 50 Hz

Распространенным тестовым воздействием на ЦФ является сигнал основной гармоники, наложенный на затухающую экспоненту. Такой сигнал

в первом приближении имитирует простейший переходный процесс в энергосистеме.

На рис. 2 приведена реакция трех фильтров (МНК, КФ и ФОС) на воздействие вида

nT

y(nT) = e т + sin

2nn

V ~N~ ,

где N = 24; т = 0,05 с; Т = 0,02/Ж

По истечении переходного процесса, который длится один период промышленной частоты, все фильтры правильно определяют амплитуду основной гармоники.

и, В 1,00,5"

0

0,005

0,010

0,015

0,020 0,025 0,030 0,035 t, с 0,040

Рис. 2. Определение амплитуды основной гармоники сигнала цифровыми фильтрами: МНК - метод наименьших квадратов; КФ - косинусный фильтр;

ФОС - формирователь ортогональных составляющих

Fig. 2. Signal fundamental component magnitude estimation by digital filter: МНК - least error squares; КФ - cosine filter; ФОС - orthogonal components former

В действительности информация о параметрах функционирования энергосистемы передается в устройства релейной защиты через измерительные трансформаторы тока. Поэтому проверку работоспособности ЦФ целесообразно проводить сигналом, приближенным по своей форме к вторичному току ТТ, который можно получить, используя методы математического моделирования.

С использованием блоков системы Simulink-SimPowerSystems [4] разработана модель фрагмента электрической сети и основных элементов измерительного органа тока цифровой защиты [12], представленная на рис. 3.

Приведем описание структурных элементов модели на рис. 3.

1. Трехфазная группа трансформаторов тока (ТТ). Математическая модель трехфазной группы ТТ с соединением вторичных обмоток и нагрузок в «звезду с нулевым проводом» содержит уравнения трех однофазных ТТ [13]. После выражения неизвестных геометрических параметров магнитной системы ТТ через его паспортные данные и принятую характеристику намагничивания [14] система уравнений, описывающая упрощенную модель трехфазной группы ТТ, в общем случае примет следующий вид (для ТТ с /2ном = 5 А, обмотка класса 10P) [15, 16]:

^ 83,34 {(0,0004/1ноМ + Д„ Д ] + 11' ] авс Л " Ком (0,0004/!„ом + й„ом;) ' ] " ' ' '

Ч ]=Р~-2,25 -10-4 К„ом Ну,

1ном

Вщ= I (Н]);

' =12Л +12Б +12€,

где /1ном - номинальный первичный ток ТТ; К„ом - номинальная предельная кратность ТТ; Лном - номинальное активное сопротивление вторичной нагрузки ТТ; Ян - действительное активное сопротивление вторичной нагрузки ТТ; Я0 - активное сопротивление нулевого провода; /0 - мгновенное значение тока, протекающего в нулевом проводе ТТ; Бт = 1(Н) - средняя характеристика намагничивания электротехнической стали; '2 -мгновенное значение первичного и вторичного токов ТТ соответственно.

Рис. 3. Структура модели Fig. 3. Simulation model structure

2. Энергосистема, в которой задаются: номинальное напряжение ином, ток трехфазного КЗ 1кз и эквивалентная постоянная времени энергосистемы т (по умолчанию т = 0,05 с).

3. Блок задания КЗ (КЗ) используется для моделирования времени наступления и вида КЗ. Во всех расчетах рассматривается трехфазное КЗ.

4. Входной преобразователь (I ^ U), который обеспечивает преобразование входного вторичного тока ТТ в пропорциональное напряжение заданного уровня.

5. Антиалайзинговый ФНЧ второго порядка (ФНЧ-2). В соответствии с теоремой Котельникова - Шеннона сигнал, подаваемый на вход аналого-цифрового преобразователя (АЦП), должен предварительно проходить через ФНЧ для подавления спектральных компонентов, частота которых превышает половину частоты дискретизации. Таким образом, для основной гармоники f = 50 Гц) граничная частота fg = 25N должна быть подавлена как минимум на порядок, т. е. коэффициент передачи на этой частоте должен быть K < 0,1. Данное условие позволяет из выражения для АЧХ звена

второго порядка выбрать частоту среза f (Гц) ФНЧ: fc = fg4K = 25N\/K.

Длительность переходного процесса для этого звена составляет tn ~ 377£.

При аппроксимации частотной характеристики фильтра по Баттерворту

1

постоянная времени принимает значение T =-, а параметр затуха-

2 п fc

ния Z = 0, 5 V2. Если принять коэффициент передачи на граничной частоте К = 0,09, то длительность установления амплитуды сигнала можно приблизительно оценить как tn ~ 5/N периодов основной частоты. При этом коэффициент передачи на частоте 50 Гц равен единице, т. е. основная гармоника сигнала не ослабляется фильтром.

6. АЦП с квантованием по времени T = 0,02/N, равным периоду дискретизации сигнала.

7. ЦФ, в котором осуществляется определение основной гармоники сигнала. Параметры блока - предварительно рассчитанные коэффициенты фильтров, например из табл. 1.

8. Блок определения амплитуды (Амплитуда) основной гармоники, реализованный в соответствии с (13).

Параметры элементов модели

Для получения входных сигналов ЦФ использовали модель ТТ с номинальным первичным 11ном = 300 A и вторичным 12ном = 5 A токами, номинальной предельной кратностью Кном = 20, действительное активное сопротивление нагрузки Rh принимали меньше номинального активного сопротивления вторичной нагрузки ТТ Rhmm. Ток трехфазного КЗ при эквивалентной постоянной времени энергосистемы т = 0,05 c принимали 1кз = 5000 А, т. е. заведомо меньше максимально допустимого тока по условиям работы данного ТТ в заданном классе точности, определяемом Кном.

Зависимости функции времени вторичного тока ТТ (сплошная линия) и приведенного к вторичной стороне первичного тока (штриховая линия) представлены на рис. 4.

250

I, А

150 100 50 0

-50 -100 -150

0

0,02 0,04 0,06 0,08 0,10 0,12 0,14 0,16 0,18 t, с 0,20

Рис. 4. Вторичный ток трансформатора тока

Fig. 4. Current transformer secondary current waveform

Наличие экспоненциально затухающей апериодической составляющей во входном сигнале тока, даже при допустимых кратностях первичного тока ТТ, приводит к насыщению его магнитопровода и искажению формы вторичного тока. Гармонический состав токов на входе и выходе ТТ за пять периодов промышленной частоты представлен на рис. 5.

b

1

Ii., .....

5 40 <

1

,lli. i.i Mllllll..

50 100 150 200 f, Гц 250

100 150 200 f, Гц 250

Рис. 5. Спектры сигнала на входе (а) и выходе (b) трансформатора тока Fig. 5. Harmonic content of the current transformer input (a) and output (b) signal

Спектр экспоненциального сигнала вида у(I) = иае т находится с помощью одностороннего преобразования Фурье [17]

<х t

Т,

F(j) = JUae т¿"dt = т

U„

т

+ j

где ю - угловая частота.

Модуль этой функции дает спектр амплитуд (спектральную плотность) сигнала

A() = Ua Т

V2 т2 +1

a

Амплитудный спектр имеет вид убывающей функции, которая содержит непрерывный спектр частот ю = 0 ... да. Приведенный на рис. 5 дискретный спектр сигнала с шагом 10 Гц объясняется тем, что гармонический анализ проводился за время в пять периодов промышленной частоты, за которое из сигнала можно выделить минимальную частоту = = 1/

N

Т an = 0.

/ * п п=1

По-иному экспоненциально затухающая апериодическая составляющая влияет на вторичный ток ТТ. Вследствие насыщения стали магнитопро-вода ТТ, уровни некоторых спектральных составляющих вторичного тока в разы превышают уровни соответствующих составляющих первичного тока. В конечном счете это приводит к тому, что форма кривой вторичного тока существенно отличается от синусоидальной. Причем в первые несколько периодов сигнал настолько искажен, что величина амплитуды вторичного тока (рис. 4) становится значительно меньшей приведенного к вторичной стороне первичного тока, что недопустимо затягивает время установления выходного сигнала. Характерной особенностью искаженного вторичного тока является наличие в нем интергармоник (гармоник, не кратных основной), которые не учитываются ни в одной из рассмотренных моделей сигналов.

Изменение во времени амплитуды основной гармоники (50 Гц) и двух составляющих входного сигнала с наибольшими амплитудными спектрами из рис. 5Ь (2-й гармоники с частотой 100 Гц и интергармоники с частотой 60 Гц) показано на рис. 6. По мере затухания апериодической составляющей, которое длится t ~ 3т с, сигнал основной гармоники стремится к установившемуся значению, равному амплитуде входного сигнала, а сигналы с частотами 60 и 100 Гц уменьшаются до нуля.

Рис. 6. Определение цифровыми фильтрами амплитуд гармонических компонентов

искаженного сигнала

Fig. 6. Distorted signal harmonic components magnitude estimation by digital filter

Однако влияние этих сигналов на конечный результат фильтрации будет различным. Если допустить, что сигналы проходят через ЦФ раздельно, то 2-я гармоника полностью подавляется фильтром, а интергармоника, хотя и ослабленная, будет участвовать в формировании выходного сигнала.

Необходимо отметить, что путем математического и физического моделирования одиночных электромагнитных ТТ установлено, что их преоб-

разовательные свойства в линейных режимах работы достаточно точно отображаются передаточной функцией полосового частотного фильтра [17]. Поэтому спектр вторичного тока ТТ вблизи нулевой частоты ослабляется самим ТТ (рис. 5Ь), что при прочих равных условиях облегчает работу ЦФ.

Таким образом, ТТ обеспечивает передачу неискаженного вторичного тока к цепям РЗ только в случае гармонического характера тока повреждения, величина которого !кз < KномI1, и при величине вторичной нагрузки, равной или менее номинальной. Поэтому ЦФ должен не столько отстраиваться от экспоненциально затухающей апериодической составляющей тока повреждения, сколько адекватно функционировать при искаженном вторичном токе ТТ.

Переходный процесс в рассматриваемых фильтрах МНК, КФ и ФОС при определении амплитуды основной гармоники сигнала, пропорционального вторичному току ТТ, при исключении из структуры модели (рис. 3) ФНЧ-2 представлен на рис. 7а. ЦФ на основе МНК демонстрирует наихудшие показатели качества и не может быть рекомендован к использованию в цифровых устройствах РЗ. Повышенные биения выходного сигнала ФОС существенно ослабляются, если сигнал предварительно проходит через аналоговый ФНЧ-2 (рис. 7Ь). Наиболее достоверно амплитуда основной гармоники вторичного тока ТТ определяется КФ.

b

20 U, В 15

10

20 U, В 15

10

0 02 0 04 0 06 0 08 t с 0 10

0

0 02 0 04 0 06 0 08 t с 0 10

Рис. 7. Переходный процесс в цифровом фильтре без

Fig. 7. Transient response of digital filter without

а

5

5

0

В модели сигналов, использующихся для расчета коэффициентов МНК и ФОС, входит экспоненциально затухающая апериодическая составляющая, которая учитывается несколькими членами ряда Тейлора. В результате проведенных расчетов установлено, что чем больше членов ряда учитывается в модели сигнала, тем хуже качество фильтрации. Так, в МНК учитываются три члена ряда Тейлора, в ФОС - два, а в КФ - апериодическая составляющая отсутствует. Это объясняется тем, что в псевдообращенной или обратной матрице члены в строках, соответствующих коэффициентам

вида nT и , получаются на несколько порядков больше коэффициент 2т2

тов ЦФ, что вносит погрешность в их расчет. Причем это свойство не проявляется на идеализированных сигналах вида (14), но оказывает существенное влияние на результаты фильтрации при наличии в сигналах интергармоник.

Существенно повысить качество фильтрации искаженных сигналов помогают оконные функции. Любой ЦФ основной гармоники имеет АЧХ, подобную изображенной на рис. 8а, АЧХ ФОС. Такая АЧХ содержит главный лепесток, в центре которого располагаются выделяемая частота (50 Гц), и боковые лепестки, обращающиеся в нуль только на частотах, кратных основной. При этом интергармоники, которые располагаются между подавляемыми частотами, будут значительно усиливаться. Обработка сигнала с помощью оконных функций заключается в уменьшении (в пределах до 0) амплитуд всех боковых лепестков.

2,0-

1,5-

1,0-

0,5-

а ^ 0,6 is

0,4 -

0,2

100

200

300

400

500 F, Гц

Рис. 8. Амплитудно-частотные характеристики: а - формирователя ортогональных составляющих; b - то же совместно с окном Хэмминга

Fig. 8. Frequency response: а - оrthogonal components former; b - оrthogonal components former with Hamming window

a

1

Для этих целей широко используются косинусные окна, представляющие собой сумму нескольких гармонических составляющих, периоды которых укладываются на длине окна целое число раз [18]. Наиболее простым и эффективным является окно Хэмминга, коэффициенты которого вычисляются как:

м>(п) = 0,54 - 0,46ОО8

Окно Хэмминга представляет собой ЦФ, включающийся последовательно с основным ЦФ измерительного органа тока представленной на рис. 3 модели. Поэтому контролируемый сигнал претерпевает двойную фильтрацию - оконной функцией и основным ЦФ. Математически такая операция называется дискретной линейной сверткой s(n) двух векторов дискретных сигналов: коэффициентов ЦФ а(п) (п = 0, 1, ..., N - 1) и коэффициентов окна Хэмминга ^(п) (п = 0, 1, ..., М- 1), причем в общем случае N Ф М. Для вычисления линейной свертки сигналы а(п) и ^(п) сдвигают друг относительно друга, почленно перемножают и складывают

В Ма1;ЬаЬ дискретная свертка вычисляется с помощью функции оопу, для чего предварительно формируются два вектора коэффициентов. Например, для 24 коэффициентов ФОС и окна Хэмминга (табл. 1) она выглядит так:

а=[0.0 0.1638 ... 0.1638]% 24 коэффициента ФОС w=[0.08 0.0971 ... 0.08]% 24 коэффициента окна Хэмминга s=conv(a,w).

Длина выходного вектора на единицу меньше суммы элементов обоих векторов, т. е. в приведенном примере 47 коэффициентов. В результате все боковые лепестки оказываются подавленными и АЧХ ФОС совместно с оконным преобразованием принимает вид, как на рис. 8Ь. Рис. 9 позволяет оценить эффект от оконного преобразования: форма выходного сигнала комбинированного фильтра (ФОС совместно с окном Хэмминга) оказывается более гладкой (кривая 2), но время его установления существенно дольше, чем у ФОС (кривая 1).

Поэтому если для целей РЗ допустимо иметь время обработки контролируемого сигнала более одного периода промышленной частоты, то применение оконных функций оказывается эффективным. Кроме того, если использовать меньшее число коэффициентов оконной функции по сравнению с числом выборок, т. е. принять М< N, то быстродействие комбинированного фильтра увеличится. На рис. 9 этому случаю соответствует кривая 3, для которой используются 24 коэффициента ФОС и восемь коэффициентов окна Хэмминга.

В общем случае на динамику измерительного органа тока оказывают влияние три фильтра: ФНЧ-2, выбранный цифровой фильтр и оконная функция. Причем если время установления сигнала в любом из выбранных ЦФ фиксировано и примерно равняется одному периоду промышленной частоты, то длительность переходных процессов в аналоговом фильтре и

п

s(n) = X а(п)м>(п - т), п = 0,1, ..., ^ + М - 2).

окне Хэмминга напрямую связана с качеством фильтрации: чем дольше преобразование, тем лучше фильтрация. Поэтому, варьируя частоту среза ФНЧ и размер оконной функции, можно добиться оптимальных динамических характеристик измерительного органа тока цифровой защиты.

Рис. 9. Определение цифровыми фильтрами амплитуды основной гармоники искаженного сигнала: 1 - формирователь ортогональных составляющих без аналоговой фильтрации; 2 - дискретная свертка формирователя ортогональных составляющих и окна Хэмминга длиной в 24 отсчета; 3 - то же длиной в 8 отсчетов; 4 - искаженный аварийный сигнал, пропорциональный вторичному току трансформатора тока

Fig. 9. Distorted signal fundamental component magnitude estimation by digital filters: 1 - orthogonal components former without analog pre-filtering; 2 - discrete linear convolution of orthogonal components former and Hamming window (24 samples); 3 - discrete linear

convolution of orthogonal components former and Hamming window (8 samples); 4 - fault signal which is proportional to the current transformer distorted secondary current

ВЫВОДЫ

1. Разработаны алгоритмы расчета коэффициентов цифровых фильтров для микропроцессорных устройств релейной защиты.

2. В среде имитационного динамического моделирования MatLab-Simulink разработана структура модели энергосистемы и измерительного органа тока для проверки работоспособности цифровых защит.

3. Установлено, что из рассмотренных цифровых фильтров наилучшие показатели качества в условиях воздействий, приближенных к реальным, демонстрируют косинусный фильтр и фильтр на основе формирователей ортогональных составляющих. Для минимизации времени установления выходного сигнала рекомендуется вычислять его амплитуду по двум выборкам в соответствии с выражением (13).

4. Рассмотрены пути улучшения качества фильтрации сигнала, в частности использование оконных функций.

ЛИТЕРАТУРА

1. Rebizant, W. Digital Signal Processing in Power System Protection and Control / W. Rebizant, J. Szafran, A. Wiszniewski // Springer. London, 2011. 316 с.

2. Phadke, A. G. Computer Relaying for Power Systems / A. G. Phadke, J. S. Thorp. 2 edition. Chichester, England: John Wiley & Sons, 2009. 325 p.

3. Влияние насыщения трансформаторов тока на работу токовых защит / Ф. А. Романюк [и др.] // Энергетика. Изв. высш. учеб. заведений и энерг. объединений СНГ. 2010. № 1. С. 5-9.

4. Дьяконов, В. П. MatLab и Simulink для радиоинженеров / В. П. Дьяконов. М.: ДМК-Пресс, 2011. 975 с.

5. Шнеерсон, Э. М. Цифровая релейная защита / Э. М. Шнеерсон. М.: Энергоатомиздат, 2007. 549 с.

6. Sachdev, M. S. A New Algorithm for Digital Impedance Relays / M. S. Sachdev, M. A. Ba-ribeau // IEEE Transactions on Power Apparatus and Systems. 1979. No 6. Р. 2232-2240.

7. Phadke, A. G. Synchronized Phasor Measurements and their Applications / A. G. Phadke, J. S. Thorp // Springer Science & Business Media. 2008. 248 p.

8. Романюк, Ф. А. Информационное обеспечение микропроцессорных защит электроустановок / Ф. А. Романюк. Минск: Технопринт, 2001. 133 с.

9. Лоусон, Ч. Численное решение задач метода наименьших квадратов / Ч. Лоусон, P. Хен-сон: пер. с англ. М.: Наука, Гл. ред. физ.-мат. лит., 1986. 232 с.

10. Schweitzer III, E. O. Filtering for Protective Relays / E. O. Schweitzer III, D. Hou // WESCA-NEX 93. "Communications, Computers and Power in the Modern Environment". Conference Proceedings. IEEE, 1993. P. 15-23.

11. Романюк, Ф. А. Аналого-цифровая частотная фильтрация входных сигналов релейной защиты / Ф. А. Романюк, В. Ю. Румянцев // Энергетика. Изв. высш. учеб. заведений и энерг. объединений СНГ. 1994. № 3-4. С. 7-12.

12. Перспективные технологии реализации микропроцессорных защит линий распределительных сетей / Ф. А. Романюк [и др.] // Перспективные материалы и технологии: в 2 т. / под ред. В. В. Клубовича. Витебск: Витебский государственный технологический университет, 2015. T. 1. С. 115-139.

13. Романюк, Ф. А. Информационное обеспечение вычислительного эксперимента в релейной защите и автоматике энергосистем / Ф. А. Романюк, В. И. Новаш. Минск: ВУЗ-ЮНИТИ, 1998. 174 с.

14. Королев, Е. П. Расчеты допустимых нагрузок в токовых цепях релейной защиты /

E. П. Королев, Э. М. Либерзон. М.: Энергия, 1980. 210 с.

15. Новаш, И. В. Упрощенная модель трехфазной группы трансформаторов тока в системе динамического моделирования / И. В. Новаш, Ю. В. Румянцев // Энергетика. Изв. высш. учеб. заведений и энерг. объединений СНГ. 2015. № 5. С. 23-38.

16. Wye-Connected Current Transformers Simplified Model Validation in MatLab-Simulink /

F. Romanyuk [et al.] // Przeglqd Electrotechniczny. 2015. Vol. 91, No 11. P. 292-295.

17. Ванин, В. К. Релейная защита на элементах вычислительной техники / В. К. Ванин, Г. М. Павлов. Л.: Энергоатомиздат, 1991. 336 с.

18. Сергиенко, А. Б. Цифровая обработка сигналов / А. Б. Сергиенко. СПб.: БХВ-Петербург, 2011. 768 с.

Поступила 31.05.2016 Подписана в печать 10.08.2016 Опубликована онлайн 26.09.2016

REFERENCES

1. Rebizant W., Szafran J., Wiszniewski A. (2011) Digital Signal Processing in Power System Protection and Control. London, Springer. 316.

2. Phadke A. G., Thorp J. S. (2009) Computer Relaying for Power Systems. 2 edition. Chichester, England, John Wiley & Sons. 325 p.

3. Romaniouk F. A., Tishechkin A. A., Roumiantsev V. Yu., Novash I. V., Bobko N. N., Glinsky E. V. (2010) Influence of Current Transformer Saturation on Operation of Current Protection. Ener-getika. Izvestiya Vysshikh Uchebnykh Zavedenii i Energeticheskikh Obedinenii SNG [Ener-getika. Proceedings of CIS Higher Education Institutions and Power Engineering Associations], (1), 5-9 (in Russian).

4. Dyakonov V. P. (2011) MatLab and Simulink for Radio Engineers. Moscow, DMK-Press. 975 (in Russian).

5. Sheyerson E. M. (2007) Digital Relay Protection. Moscow, Energoatomizdat Publ. 549 (in Russian).

6. Sachdev M. S., Baribeau M. A. (1979) A New Algorithm for Digital Impedance Relays. IEEE Transactions on Power Apparatus and Systems, (6), 2232-2240. DOI: 10.1109/TPAS. 1979.319422.

7. Phadke A. G., Thorp J. S. (2008) Synchronized Phasor Measurements and their Applications. Springer Science & Business Media. 248.

8. Romaniuk F. A. (2001) Information Support of Microprocessor-Based Protection of Electrical Installations. Minsk, Tekhnoprint. 133 (in Russian).

9. Lawson Ch., Henson R. (1986) Numerical Solution of Problems of the Method of Least Squares. Moscow, Nauka. 232 (in Russian).

10. Schweitzer III E. O., Hou D. (1993) Filtering for Protective Relays. WESCANEX 93. "Communications, Computers and Power in the Modern Environment". Conference Proceedings. IEEE, 15-23. DOI: 10.1109/WESCAN.1993.270548.

11. Romaniuk F. A., Rumiantsev V. Yu. (1994) Analog-to-Digital Frequency Filtering of the Input Signals of Relay Protection. Energetika. Izvestiya Vysshikh Uchebnykh Zavedenii i Energeticheskikh Obedinenii SNG [Energetika. Proceedings of CIS Higher Education Institutions and Power Engineering Associations], (3-4), 7-12 (in Russian).

12. Romanyuk F., Novash I., Rumiantsev V. Yu., Rumiantsev Yu.V. (2015) Promising Implementation Technologies of Microprocessor-Based Protections of the Lines of Distribution Grids. V. V. Klubovich Prospective Materials and Technology. Vol. 1. Vitebsk, Vitebsk State Technological University, 115-139 (in Russian).

13. Romanyuk F. A., Novash V. I. (1998) Information Support of Computing Experiment in Relay Protection and Automation of Power Systems. Minsk, VUZ-YUNITI. 174 (in Russian).

14. Korolyev Ye. P., Liberzon E. M. (1980) Calculations of Permissible Loads in the Current Circuits of Relay Protection. Moscow, Energiya. 210 (in Russian).

15. Novash I. V., Rumiantsev Y. V. (2015) A Simplified Model of Three-Phase Bank of Current Transformers in the Dynamic Simulation System. Energetika. Izvestiya Vysshikh Uchebnykh Zavedenii i Energeticheskikh Obedinenii SNG [EnergetikA. Proceedings of CIS Higher Education Institutions and Power Engineering Associations], (5), 23-38 (in Russian).

16. Romanyuk F., Novash I., Rumiantsev Y., Wçgierek P. (2015). Wye-Connected Current Transformers Simplified Model Validation in MatLab-Simulink. Przeglqd Elektrotechniczny, 91 (11), 292-295. DOI: 10.15199/48.2015.11.67.

17. Vanin V. K., Pavlov G. M. (1991) Relay Protection on the Elements of Computers. Leningrad, Energoatomisdat. 336 (in Russian).

18. Sergiyenko A. B. (2011) Digital Signal Processing. St.-Petersburg, BHV-Petersburg. 768 (in Russian).

Received: 31 May 2016 Accepted: 10 August 2016 Published online: 26 September 2016