Методические вопросы расчета распространения грозовых волн в коронирующей линии электропередачи методом бегущих волн Текст научной статьи по специальности «Физика»

УДК 621.311

Б.В.Ефимов, Н.И.Гумерова

МЕТОДИЧЕСКИЕ ВОПРОСЫ РАСЧЕТА РАСПРОСТРАНЕНИЯ ГРОЗОВЫХ ВОЛН В КОРОНИРУЮЩЕЙ ЛИНИИ ЭЛЕКТРОПЕРЕДАЧИ МЕТОДОМ БЕГУЩИХ ВОЛН*

Аннотация

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

Ключевые слова:

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

B.V.Efimov, N.I.Gumerova

METHODOLOGICAL ASPECTS OF CALCULATION THE SURGE WAVES PROPAGATION ON TRANSMISSION LINE WITH CORONA USING TRAVELING WAVES METHOD

Abstract

Problems of calculation of microseconds wave processes at power transmission lines by moving wave technique, allowing for corona influence to shunt admittance matrix, are considered in detail. Corona was simulated by discretely situated nodes with multiport circuits, consisting by nonlinear additions to coefficients of electrostatic induction. For the corona simulating inaccuracies analysis was carried out. Mode and wave calculation method results coincide for the 110 kV transmission line model.

Keywords:

surge waves, distorting nodes, corona, modeling of nonlinear circuits with lumped parameters.

В основе всех известных методов расчета деформации фронтов волн при напряжении выше начала короны лежат нелинейные вольт-кулоновые характеристики коронирующих проводов, полученные на высоковольтных экспериментальных установках. В частности, такие эксперименты, выполненные на участках реальных проводов ВЛ, подробно описаны в статьях [1, 2]. В этих же работах даны аппроксимирующие выражения для динамического потенциального коэффициента коронирующего провода. Показано, что с большой степенью точности в многопроводной постановке задачи можно считать зависящими от зарядов на проводе и вокруг него только собственные потенциальные коэффициенты. Все взаимные потенциальные коэффициенты можно полагать постоянными, вычисленными по обычным формулам в электростатическом приближении. Исходя из этих соображений разработан метод расчета волновых процессов в коронирующих ВЛ, основанный на разложении приращений напряжений и зарядов на проводах на модальные

* Работа выполнена при финансовой поддержке РФФИ (проект № 11-08-00690).

составляющие, каждая из которых распространяется по линии со своей скоростью. Будем называть этот метод модальным (ММ).

Физика процесса и практические вопросы численной реализации этого метода для двухпроводной и многопроводной линий подробно рассмотрены в работах [3, 4]. Предполагается, что каждая точка фронта исходной волны (фиксированного значения напряжения) при х=0 распространяется по линии с эквивалентной скоростью, несколько меньшей скорости света в воздухе (с), то есть при наличии короны происходит дополнительное запаздывание точек фронта, нелинейно нарастающее (для фиксированного х) по мере увеличения заряда на проводе и вокруг него. С увеличением х дополнительное запаздывание растет линейно. Эти два обстоятельства и формируют фронт волны. Такое запаздывание достаточно просто вычисляется для однопроводной линии [5]. Более сложно оно определяется в многопроводной постановке задачи [3, 4], но суть модели от этого не меняется. Метод пригоден только для однородной линии без потерь, которая начинается при х=0 и уходит в бесконечность (х ^ да). Тем не менее, при решении телеграфных уравнений однородной линии (при заданной модели короны) не вводится никаких новых допущений. Решение для скоростей является точным как в однопроводной, так и в многопроводной постановках задачи. Метод хорошо сходится при увеличении числа точек на фронте волны и на современных ПЭВМ является достаточно быстродействующим. Поэтому такой способ расчета деформации фронтов волн можно принять в качестве эталонного.

Принципиально другой подход используется при расчете волнового процесса в линии методом бегущих волн с дискретно расставленными искажающими узлами, моделирующими корону. Для той же однородной линии считается, что волны по всем проводам распространяются только со скоростью с. Потери в линии отсутствуют. Линия начинается при х=0, и ее длина настолько велика, что отражения от дальнего конца не успевают прийти в заданную точку х>0 во всем рассматриваемом диапазоне времен. При напряжении ниже начала короны волны, заданные при х=0, распространяются без искажений. Алгоритм расчета в этом случае обычный, описанный во многих работах. Исходя из некоторых неформализованных соображений о необходимой точности описания фронта волны, задается шаг по времени А?. Например, для волны с экспоненциальным нарастанием напряжения на фронте с постоянной времени гФ=0.2 мкс задается А? =0.01 мкс, или примерно 50 точек на весь фронт. Он определяет шаг по длине линии Ах=А? с (в примере - 3 м). Для каждого провода отводится 2 массива памяти для хранения значений падающих (прямых) и отраженных (обратных) волн. Прямые волны на каждому-м шаге перемещаются на один элемент в сторону возрастания номеров элементов массива, что моделирует продвижение волн на Ах и одновременно переход от к ^+1=^+А?. Пока отраженных волн нет, но можно считать, что во втором массиве в обратную сторону, то есть уменьшения номеров элементов массива, на том же у-м шаге по длине и времени перемещаются нули. Фактически элементы обоих массивов с одинаковыми номерами моделируют промежуточный узел схемы замещения линии, который волны проходят без искажений. В соответствии с работой [6] назовем такие узлы «пустыми». Продолжим пример. Для линии длиной 3 км и шага 3 м число таких узлов составит 1000. Время пробега волны

по этой линии без потерь - 10 мкс. Предположим, что будут рассчитываться формы волн в конце линии также до времени 10 мкс. Тогда общее время счета составит 20 мкс или 2000 шагов по длине и времени.

При превышении коронного порога хотя бы на одном из проводов, включаются искажающие узлы, расставленные через расстояние , кратное Ах. Пусть 1Т =10 Ах. Тогда искажающие узлы вставляются в схему замещения вместо каждого десятого узла, или через 30 м. Отсюда корона моделируется участками 1у, в пределах каждого из которых параметры коронирующей линии считаются постоянными по длине. По времени эти параметры изменяются также ступенчато, но значительно чаще - с шагом А?. Волны в каждом искажающем узле преломляются и отражаются. При этом многократно отраженные и преломленные волны в искажающих узлах играют очень существенную роль для формирования фронта волны на заданном расстоянии от начала линии. В частности, именно они формируют повышенные наведенные напряжения на изолированных от земли проводах по сравнению с наводками, рассчитанными в электростатическом приближении. Если в расчетах рассматривается промежуточная точка линии, то нужно учитывать отраженные волны, идущие с участка линии, расположенного за рассматриваемой точкой. Точное согласование линии в конце невозможно из-за переменных параметров коронирующих проводов. В этом состояла одна из ошибок экспериментов на опытных пролетах воздушных линий [7]. Единственный выход в сопоставительных расчетах модальным и волновым методом состоит в том, что в волновом методе линия задается настолько большой длины, чтобы отражения от конца не влияли на процессы в рассматриваемой точке. В примере нужно удлинить линию на 1.5 км. Тогда первое отражение от конца линии придет в точку х=3 км через 10 мкс. Таким образом, в примере для расчета формы волны при х=3 км нужно иметь схему замещения линии длиной не менее 4.5 км, состоящую из 1500 узлов, из которых 50 будут являться искажающими. Далее в методических целях шаг по длине уменьшается до 0.5 м, а число искажающих узлов увеличивается до 300. Набрать такую схему в комплексе АТР трудно, поэтому в данной работе проводится сопоставление результатов расчетов только по двум алгоритмам - модальному и волновому.

Принципиальное различие двух рассматриваемых подходов к расчету искажения волны короной демонстрирует рис.1.

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

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

Рис.1. Способы расчета искажения фронта волны:

а - путем дополнительного запаздывания мгновенных значений волны в модальном методе (ММ); б - путем снижения напряжений вследствие заряда емкости коронирующего провода в волновом методе (ВМ)

В методе бегущих волн волна моделируется столбиками (при А?^0 также бесконечно малой толщины). Столбики друг за другом распространяются по линии со скоростью с независимо от наличия короны. Коронный искажающий узел, в любом случае имеющий емкостной характер сопротивления, постепенно уменьшает высоту столбика за счет заряда емкости через волновое сопротивление линии. Энергия волны «уходит» с провода. Все точки исходной волны бегут по линии синхронно, но по мере продвижения по ней опускаются вниз. Поэтому одинаковые значения на фронте деформируемой волны в двух методах получаются из различных точек волны напряжения при х=0. На рис.1, для примера, отмечена точка, полученная сдвигом по горизонтали напряжения, ненамного превышающего коронный порог, и опусканием напряжения, близкого к амплитудному значению волны.

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

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

Поскольку оба подхода дают очень близкие результаты расчетов и в целом соответствуют экспериментальным данным, можно считать, что и модальный, и волновой метод при всем различии алгоритмов эквивалентны.

Далее принята следующая модель короны [2]. Рассматривается одиночный цилиндрический провод с гладкой поверхностью, подвешенный на постоянной высоте над землей. Поскольку основная часть дальнейшего материала посвящена многопроводным системам, то для упрощения перехода к /-му проводу и в однопроводном случае будем индексировать все геометрические и электрические параметры, присвоив ему номер 1.

Независимо от скорости нарастания напряжения, приложенного к проводу, начиная с некоторого критического значения напряженности поля на поверхности провода ЕКР1 (соответственно критического значения заряда на

единицу длины провода - ЯКР1), вокруг провода образуется объемный заряд той же полярности, что и приложенное напряжение.

Для критических значений напряженности электрического поля, заряда и напряжения на единицу длины провода можно использовать формулы [2]:

(

ЕКР; = 24.6 • 0.82 • 105

1 +

0.0301

\

(В/м);

ЯКР1 = 2^^ 0 Г1 Е КР1 (Кл/м);

икр1 = ГЕкр1 1п—1- (В),

(1)

(2)

(3)

где Т\ и Н\ - радиус и средняя высота подвеса одиночного провода (м); 0.82 -коэффициент, учитывающий местное повышение напряженности поля на поверхности из-за витой структуры проводов.

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

Для отрицательных разрядов молнии в работе [8] были получены такие зависимости в виде:

Г1 при Я ^ Якр1

2Н-

= а Г

Я:

1 -±Щ ^КР^

-1 + к

(4)

при Я > ЯкР1 ,

где п11 = 1п—1 - логарифмический множитель собственного потенциального

коэффициента провода;

1Г1 - , "11

2ле„

геометрический потенциальный

коэффициент на единицу длины провода для одиночного провода; к = [0.7(пп - 4)2 + 4 I-1 - эмпирический коэффициент, полученный в результате обработки данных ряда экспериментов.

к

п

1

Можно получить и интегральную зависимость для отрицательной

коРоны и!-} = /(д,) при д, > дт :

и[ -1 = а Г

41

і'^

дКР 1

ді

і __ЦП Зт.

-1 + к

ёд1

(5)

Взяв интеграл, получим:

,(_)

ді

і

и к

дКР

ді

дКР

' і

1 + • \п

к V

дКР

- + к -1

/у

дКР

■ + к -1

(6)

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

Для приближенного учета характеристик короны положительного знака в статье [2] получена следующая зависимость:

і du1(

(+)

а

д і

dді

і

ді

- і + к • т2

і - т^\п д крі

к • т2

пРи ді ^ дкрі

пРи ді > дкрі =

(7)

где «1=0.77 и т 2=0.082 - эмпирические коэффициенты.

Соответственно интегральное выражение и1 = /(д1) будет:

ді

т,

и

д,

дКР

+ к • т2 - і

кр і (

дКР

• \п і

к

д,

дКР

+ к • т2 - і

/у

д,

дКР

+ к • т2 - і

(8)

В дальнейшем все рассуждения ведутся для отрицательной короны и верхний индекс (-) опущен.

В статье [2] также показано, что при наличии короны на любом числе проводов изменяются только их собственные потенциальные коэффициенты, которые зависят только от зарядов на этом проводе и вокруг него.

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

• собственный потенциальный коэффициент /-го провода (на единицу длины):

1 , 2Нг 1

а г гг = -1п----= -----пп ; (9)

2лв0 гг 2лв0 ^ ’

к

п

д

д

п

п

X

п

X

взаимный потенциальный коэффициент между >м и к-м проводами:

1 , Ягк 1

“ Гг* =2;т1п и=п‘ • о»)

0 г* 0

где - кратчайшее расстояние между г-м и *-м проводами; Б* кратчайшее

расстояние между г-м проводом и отражением в земле *-го провода.

Далее введем собственные динамические потенциальные коэффициенты:

( п \

- - 1 + *г

йиг 1 , * Ч 1 I

-(пгг -А г ) =---- Пгг - 1п

Дгг 7 4 гг г * г\

апг 2т 0 2т 0

п

ПКРг

*г

(11)

Величину

А г = 1п

ПКРг

■ - 1 + *,

*г

> 0

(12)

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

Все изложенное относится к фронтовой части волн. Экспериментальные данные говорят о независимости процессов образования короны от скорости изменения зарядов на проводах в пределах изменения фронтов волн от 0.1 до 10 мкс. Поэтому для расчета грозовых перенапряжений запаздывание изменения потенциального коэффициента даже при самых крутых волнах напряжения на проводе можно не учитывать.

В момент достижения максимума заряда на проводе создание нового объемного заряда прекращается. В начале спадающей части импульса объемный заряд остается в окружающем провод пространстве. Это видно из вольт-кулоновых характеристик [2]. Заряд на самом проводе уменьшается пропорционально спаду напряжения. Из этого следует, что собственный потенциальный коэффициент провода становится равным геометрическому.

Матрицу погонных потенциальных коэффициентов, умноженных на коэффициент 2лв 0, обычно называют матрицей логарифмов многопроводной линии и обозначают через N. Эта матрица (и только она) характеризует геометрию однородной некоронирующей линии без потерь и определяет ее поперечные проводимости, то есть потенциальные коэффициенты:

N =

21

21

2п

(13)

, 2й. , Б*

где Пг, = 1п—^ и П* = Пк1 = 1П

Гг йгк

Для коронирующей линии для потенциальных коэффициентов необходимо вводить матрицу с переменными коэффициентами:

п

п

п

п

11

12

пп1 пп2

п

пп

п11 - А1 п

п21

12

- А

21 2

п

п

п2

ппп - Ап

Каждая из поправок Аг является функцией пг, то есть зависит от заряда на г-м проводе и для волн апериодической формы определяется по следующему алгоритму:

• при заряде меньше критического она равна нулю;

• при превышении зарядом критического значения она определяется по формулам, приведенным выше;

• в момент достижения зарядом амплитудного значения поправка скачком снова становится равной нулю.

При идеально проводящей земле и пренебрежении продольными составляющими тока короны все собственные и взаимные индуктивные параметры считаются постоянными. Их нужно определять по формулам:

Ь

^11

м

М

12

^22

м

м

1п

2п

мп1 м

п2

= г°- N

2л

(15)

индуктивностей проводов и взаимных

где Ь матрица продольных индуктивностей между ними.

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

Корона влияет на потенциальные коэффициенты, коэффициенты электростатической индукции и частичные емкости линии, то есть на поперечные проводимости. Поэтому будем называть искажающие узлы -7-узел [6].

Обозначим матрицу потенциальных коэффициентов, рассчитанных на длину участка 17, через Ау. В этой матрице, как и ранее, все внедиагональные коэффициенты постоянны, а диагональные могут изменяться в зависимости от зарядов на проводах, причем а7 гг = а д гг /17 , где а д гг - погонный собственный

потенциальный коэффициент коронирующего г-го провода, вычисленный по формуле (11).

Обращая А7 , получим матрицу собственных и взаимных коэффициентов электростатической индукции, в которой будут зависеть от зарядов на проводах уже все коэффициенты. Обозначим эту матрицу В7 = А-1. При отсутствии короны В7 переходит в обычную матрицу коэффициентов электростатической индукции для многопроводной линии без потерь длиной 17. Обозначим эту матрицу, коэффициенты которой определены на длине /7, через ВГ. Соответственно, теперь исходная матрица геометрических потенциальных коэффициентов (А Г) тоже считается заданной на длине 17 (ВГ = АГ1).

п

п

Введем новую матрицу добавок к коэффициентам электростатической индукции - АБ, равную поэлементной разности Ву и ВГ, то есть АВ = Ву - ВГ, где собственные и взаимные коэффициенты АВ зависят от зарядов на всех проводах.

Теперь будем считать, что участок длиной їу является многопроводной линией без потерь. Матрица ВГ определяет параметры этой линии. Остается АВ. Влияние короны учитывается многополюсником с сосредоточенными параметрами, включенным в конце участка. Пример схемы для трехпроводной линии дан на рис.2а.

Рис. 2. Схема замещения коронирующего участка линии:

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

Параметры многополюсника - это собственные и взаимные коэффициенты матрицы АВ. Изобразить такой многополюсник

с использованием обычных элементов электрических схем нельзя. Для наглядности от АВ можно перейти к частичным емкостям по формулам:

АСік - -АР ,к и АС и - ^АР ,к (рис.2б).

' гк

к=1

Неискажающий участок линии имеет матрицу волновых сопротивлений Z^ , которая вычисляется по формуле: Zw = 60 • N.

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

схему, состоящую только из сосредоточенных элементов (рис.2в).

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

2Zw + и7 (^ = и Л . (16)

При начальных условиях иу = иу (tj) = ы¥(j). Аналогичное уравнение можно записать для волн, набегающих справа (иП). Решение уравнения (16) на j-м шаге при времени Аt после его начала даст искомый вектор иу(j+1).

Отметим, что порядок сомножителей Zw • АВ фиксирован. Матрица АВ имеет специфическую структуру. Ясно, что при напряжении ниже начала короны В7 = ВГ и АВ = 0, то есть равны нулю все ее коэффициенты.

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

Далее все численные данные приведены для примера трехпроводной линии с горизонтальным расположением проводов, подвешенных на высоте к 1=10 м и расстояниями между ними Ь 12=Ь2з=4 м; 613=8 м. Радиус всех витых проводов Г1=1 см. Критическая напряженность поля - £’КР1=31.9кВ/см. Критическое напряжение на одиночном проводе - иКР1=242.3 кВ. Погонный критический заряд на проводе пКР1=1.77 мкКл/м.

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

Матрицы потенциальных коэффициентов рассматриваемой линии будут:

А г =

1

2л8д 1Г

" 7.6 1.63" • А = 1 ' 3.8 1.63"

1.63 7.6 2 Л80 /у 1.63 7.6

Соответственно обратные им матрицы:

В г — 2^80/7

" 0.1379 - 0.0296" " 0.2898 -0.0621"

; у 2л80/у *

-0.0296 0.1379 -0.0621 0.1449

Частичные емкости (Ф/м):

Сг — 2^80/7

Теперь перейдем к добавкам, выносимым в искажающий узел:

"0.1083 0.0296" "0.2277 0.0621"

; у — 2л80/у *

0.0296 0.1083 0.0621 0.0828

АВ = 2л80/7

0.1519

-0.0325

-0.0325

0.00698

; АС = 2тс80/.

0* у

0.1193 0.0325

0.0325 - 0.0256

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

Для численной реализации алгоритма знаки частичных емкостей не имеют значения. Однако это принципиально важно при физическом моделировании короны в многопроводной постановке задачи.

Для трехпроводной линии все соображения остаются такими же. Добавляется взаимный потенциальный коэффициент между первым и третьим проводом равный аГ13 = 0.9905 / 2л80/у. Матрица АВ при двукратном изменении потенциального коэффициента первого провода теперь имеет вид:

Теперь уже и для взаимных коэффициентов индукции (и, следовательно, взаимных частичных емкостей) происходит изменение знака. Это связано с ослаблением непосредственных связей между вторым и третьим проводом при увеличении частичных емкостей обоих проводов на первый. Можно повторить, что для численной реализации метода знаки коэффициентов в АВ не имеют практического значения.

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

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

Характеристическое матричное уравненй для выражения (16) будет:

АС2 = АС22 +

АС11 +АС12

0.1556 - 0.0304 - 0.01376

АВ = 2л80/у - 0.0304 0.00594 0.00269

- 0.01376 0.0269 0.00122

С

V

где Е - единичная диагональная матрица; W = - собственные векторы СВ

матрицы Н; А = -А' - собственные значения СЗ матрицы Н; 1=1..п (п - число проводов линии).

Из физических соображений ясно, что собственные значения для матрицы, характеризующей схему, состоящую только из емкостей и активных сопротивлений, будут вещественными и положительными числами. Они равны обратной величине постоянных времени схемы (= 1/Аг-). Естественно при реализации метода матрицу Z№ -ДБ/2 обращать не нужно. Ее собственные

значения сразу дают постоянные времени.

Искомое решение уравнения (16) в конце у-го шага по времени будет:

где и7;- - вектор напряжений в узле, найденный на предыдущем шаге по времени; W - квадратная матрица, составленная из собственных векторов матрицы Н;

е г - диагональная матрица, являющаяся функцией собственных значений матрицы Н.

Отраженные и преломленные волны, уходящие вправо и влево от Г-узла, определятся как разности найденных напряжений и соответствующих падающих волн.

Теперь можно сформулировать алгоритм расчета напряжения в Г-узле. Перед началом цикла по времени и длине вычисляются все матрицы с постоянными коэффициентами.

1. Находятся суммы правых и левых волн для всех проводов.

2. Определяется вектор зарядов на проводах по выражению q. = А,

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

заряд, включается алгоритм расчета деформации волн.

3. Определяется соответствующий динамический потенциальный коэффициент участка линии, относящийся к данному Г-узлу, и находится матрица АГ. .

4. Обращается матрица потенциальных коэффициентов и находится БГ..

5. Находится ДБ .

6. Находится — Z № - ДБ

иТи+1) = К,- + «Л, ) - Wе г W+ «Л, ) - Ь

(19)

7. Находятся собственные значения и собственные векторы этой матрицы.

8. Матрица собственных векторов обращается.

ы

9. Находится диагональная матрица е т .

10. Вычисляется вектор [(иП;. + иЛ;.) - иГ;.].

11. Полученный вектор в соответствии с выражением (19) последовательно умножается на квадратную матрицу, затем на диагональную и снова на квадратную, после чего вычисляются разности векторов по формуле (19) и тем самым находится иу (.+1).

12. Определяются новые значения зарядов по выражению q.+1 = А у,иГ (.+1). Если новое значение какого либо из элементов q(.+1) становится меньше соответствующего значения q ., то на данном проводе

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

Несмотря на кажущуюся огромную трудоемкость алгоритма (нахождение СЗ и СВ нужно выполнять в каждом Г-узле на каждом шаге по времени), он оказывается вполне реализуемым на современных ПЭВМ. В программе, результаты расчетов по которой приведены ниже, использована стандартная процедура ЕУСЯО() из библиотеки математических программ 1М8Ь.

Сопоставление модального и волнового метода расчета проводилось для описанной выше линии при подаче волны на первый провод. Амплитуда волны была задана равной 600 кВ. Это близко к горизонтальной части вольт-секундной характеристики гирлянд изоляторов линий класса напряжения 110 кВ. Форма волны - разность двух экспонент с постоянными времени 1000 и 0.2 мкс.

Длина линии - 4.5 км. Шаг по длине варьировался от 6 до 0.5 м. Искажающие узлы расставлялись через 2, 5, 10, 20 и 30 шагов. Оптимальным оказался шаг по длине 1 м. Расстояние между Г-узлами - 30 м. При этом общее число узлов составило 4500. Число искажающих узлов 150. Выводились напряжения в диапазоне времен 1=(0^10) мкс на расстояниях 1.5 и 3 км от начала, то есть в 1500-м и 3000-м узлах. Время счета на ПЭВМ с тактовой частотой 2.5 ГГц составило 3.5 сек.

Результаты расчетов по модальному и волновому методам приведены на рис.3. Волна подавалась на первый провод. Сопротивления, подключенные к параллельным проводам при х=0, полагались или много больше волновых сопротивлений линии (рис.3а), или много меньше их (рис.3б). В первом случае волны, приходящие к началу линии от искажающих узлов. полностью отражались с тем же знаком. Во втором - с обратным знаком, при этом заряд на втором проводе оказывался настолько большим, что на нем начиналась корона обратного знака, как на грозозащитном тросе. Для исключения этого явления (только для этого расчетного варианта) критический заряд второго провода был искусственно увеличен в 2 раза. На рис.3 данные по обоим алгоритмам приведены для сопоставимых условий коронирования только первого провода.

Рис. 3. Сопоставление деформации волн при расчетах модальным (ММ) и волновым (ВМ) методами:

а - при изолированных параллельных проводах; б - для параллельных проводов, заземленных при х=0

Из рис.3 видно, что получено практически полное наложение результатов расчетов обоими методами. Действительно, при правильно подобранных шагах по длине, времени и расстановке Г-узлов эти методы эквивалентны. Ранее приводились данные со значительно большими расхождениями [9]. Основным фактом, влияющим на точность счета волновым методом при заданном 1Г, является необходимость введения большого числа пустых узлов, то есть нужно выбирать Ах<< 1Г. Предельная величина определяется крутизной фронта волны. В рассмотренном примере увеличение более 60 м приводило к появлению заметных колебаний на графиках напряжений деформированных волн. Эти колебания вызываются резкими скачками параметров линии и развитием волновых процессов между искажающими узлами. При 1Г, равном 30 м и менее, амплитуда этих колебаний становится меньше толщины линий на графиках. Казалось бы, можно выбрать шаг по длине линии Ах такой же, как и 1Г. Однако практика использования

волнового метода показала, что влияние короны при этом будет значительно слабее, чем в модальном методе. Одним из факторов возникающих погрешностей является недостаточная точность совпадения на одном из шагов текущих и критических значений зарядов, то есть начало коронирования может запаздывать на шаг. Но даже при очень малом Ах, например 0.3 м (соответственно А^1 нс), требование Ах << 1Г сохраняется. Вообще, серия расчетов показывает, что выбор >10Ах и более необходим для получения высокой точности счета волновым методом. Это требование при заданном 1Г приводит к необходимости значительного уменьшения шага по длине и времени. Причины этого явления в настоящее время остались невыясненными.

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

0.01 мкс и более. Соответственно нижний предел шага по длине можно принять 3 м, а это уже величины Ах, которые нужны для точного описания схем замещения самой подстанции. Для современных схем грозозащиты характерными длинами подходов ВЛ к подстанции являются расстояния не более 1^1.5 км (а не 3 км, как это принято с запасом в настоящей работе). Более дальние удары молнии всегда безопасны для оборудования подстанций. При этих условиях можно отметить, что полученная высокая скорость и точность счета волновым методом снимает практически все ограничения на подробность описания реальных характеристик подхода ВЛ к подстанции.

Выводы

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

2. Получено хорошее совпадение форм кривых при временах от 0.01 до 10 мкс при расчетах по двум независимым методам, алгоритмам и программам: по методу, использующему многоскоростную модель распространения электромагнитных волн вдоль коронирующей линии, и методу бегущих волн с моделированием короны дискретно расставленными искажающими узлами. Каждый из этих алгоритмов имеет свои области применения. В частности, моделирование короны искажающими узлами наилучшим образом вписывается в общий алгоритм расчета волновых процессов в линиях с потерями и неоднородностями по длине.

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

Литература

1. Вольт-кулоновые характеристики короны на расщепленных проводах при импульсном напряжении / И.Н.Богатенков, Н.И.Гумерова, М.В.Костенко и др. // Тр. ЛПИ. 1974. № 340. С. 8-13.

2. Погонные параметры коронирующей многопроводной линии электропередачи, расположенной над идеально проводящей землей / НИ.Гумерова, Б.В.Ефимов // Моделирование переходных процессов и установившихся режимов высоковольтной сети. Апатиты: Изд. КНЦ РАН, 2008. С.7-16.

3. Анализ влияния короны в двухпроводной линии, подвешенной над идеально проводящей землей / Н.И.Гумерова, Б.В.Ефимов // Моделирование переходных процессов и установившихся режимов высоковольтной сети. Апатиты: Изд. КНЦ РАН, 2008. С.16-38.

4. Распространение грозовых волн в многопроводной коронирующеей линии, подвешенной над идеально проводящей землей / Н.И.Гумерова, Б.В.Ефимов // Труды КНЦ РАН. Энергетика. 2011. Вып. 2. С.66-78.

5. Техника высоких напряжений / под ред. М.В.Костенко. М.: Высш. шк., 1973. 528 с.

6. Анализ надежности грозозащиты подстанций / М.В.Костенко, Б.В.Ефимов, И.М.Зархи, Н.И Гумерова. Л.: Наука, 1981. 127 с.

7. Effekt of corona in traveling waves / Wagner C.F., Lloyd B.L. // Trans AIEE. 1955. Vol. 74. Pt 3.

8. Перенапряжения и защита от них в воздушных и кабельных электропередачах высокого напряжения / М.В.Костенко, К.П.Кадомская, М.Л.Левинштейн, И.А.Ефремов Л.: Наука, 1988. 302 с.

9. Ефимов Б.В. Грозовые волны в воздушных линиях. Апатиты: Изд-во КНЦ РАН, 2000. 134 с.

Сведения об авторах

Ефимов Борис Васильевич,

директор Центра физико-технических проблем энергетики Севера КНЦ РАН, д.т.н.

Россия, 184209, Мурманская область, г. Апатиты, мкр. Академгородок, д. 21А

Эл. почта: efimov@ien.kolasc.net.ru

Гумерова Натэлла Идрисовна,

доцент кафедры «Электроэнергетика, техника высоких напряжений» Санкт-

Петербургского государственного политехнического университета, к.т.н., ст.н.с.

Россия, Санкт-Петербург, Политехническая ул., д.29

Тел. 8-911-257 3809,

Эл. почта: nigumerova@mail.ru