Теплофизика высоких температур, 2020, T. 58, № 5, стр. 770-781

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

Э. Е. Сон 1*, В. С. Бондарь 2, Ю. М. Темис 3, Х. Х. Азметов 3

1 Объединенный институт высоких температур РАН
Москва, Россия

2 Московский политехнический университет
Москва, Россия

3 Центральный институт авиационного моторостроения имени П.И. Баранова
Москва, Россия

* E-mail: son.eduard@gmail.com

Поступила в редакцию 10.03.2020
После доработки 10.03.2020
Принята к публикации 18.06.2020

Полный текст (PDF)

Аннотация

Представлены методы и результаты моделирования взрывов и кинетики напряженно-деформированного состояния высоковольтного маслонаполненного электрооборудования на примере разрушения высоковольтных маслонаполненных трансформаторов при возникновении короткого замыкания. Описаны теория, методы и разработанные модули в пакетах Flow Vision, GDT (Gas Dynamic Tool) и сопряжении пакетов ANSYS и LS DYNA. Описаны физико-математические модели, лежащие в основе разработанных подходов, вычислительные схемы и интерфейсы пакетов. Проведено сравнение с экспериментальными результатами.

ВВЕДЕНИЕ

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

Высоковольные маслонапоненные трансформаторы, как и другое ВМЭО, представляют из себя объекты, надежность которых существенно ниже, чем, например, атомных станций. Если для атомных станций коэффициент надежности составляет величину порядка 10–5 и аварии случаются крайне редко, то для ВМЭО такие аварии происходят практически ежегодно как в России, так и за рубежом [1].

Взаимодействие течений газа и/или жидкости со структурами в потоке (в зарубежной литературе [26] для такого типа течений принято наименование Fluid Structure Interaction, или FSI) представляет собой взаимодействие некоторой подвижной или деформируемой структуры с внутренним или окружающим потоком жидкости или газа. Взаимодействие флюида и структуры может быть стабильным, колебательным, устойчивым или неустойчивым, приводящим в некоторых случаях к катастрофам, как ураганы и тропические циклоны приводят к разрушениям. Этот раздел механики возник еще в классических задачах Л. Эйлера, М.А. Лаврентьева, А.Ю. Ишлинского, задачах флаттера в авиации, рассмотренных М.В. Келдышем. Примерами такого типа течений являются распространение волны давления по несжимаемой жидкости в гибкой трубке (например, течение крови по гибким кровеносным сосудам), задачи аэроупругости, испарение с поверхности в пограничном слое, пленки в многофазных жидкостях, движение нефти и газа по трубам и пористым средам и многие другие явления.

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

Вследствие указанных и других трудностей даже при постоянном развитии методов прямого численного моделирования в настоящее время в основном используются пакеты численного моделирования, в одних из которых решаются газодинамические и гидродинамические задачи (ANSYS, Fluent, GDT, Comsol, Flow Vision и др.), а в других задачи деформируемых твердых тел (ABACUS, FIDESIS и др.). Пока лишь в некоторых из них (LS DYNA) решаются задачи не только деформации, но и разрушения конструкций.

Задачи о расчетах высоковольтных маслонаполненных трансформаторов в стационарных режимах и авторские пакеты программ представлены в работах [711].

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

ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ КИНЕТИКИ НАПРЯЖЕННО-ДЕФОРМИРОВАННОГО СОСТОЯНИЯ, РАЗРУШЕНИЯ И ВЗРЫВА ТРАНСФОРМАТОРОВ

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

Для анализа процесса деформирования и разрушения трансформатора была разработана методика математического моделирования напряженно-деформированного состояния (НДС), реакции конструкции на короткое замыкание, разрушение и взрыв трансформатор. На рис. 1 показано экспериментально измеренное изменение давления, действующего на стенку корпуса трансформатора, после короткого замыкания. Если короткое замыкание происходит около точки внутри трансформатора (с координатой l/lст = = 0.13), то первоначально отмечается рост давления в этой точке, затем происходит перераспределение давления по всей стенке. Такое динамическое воздействие требует создания адекватной математической модели, позволяющей учесть особенности поведения конструкционного материала и разрушения конструкции при высокоскоростном нагружении. Динамическая модель должна описывать волновые процессы в конструкции с учетом особенностей ее упругопластического деформирования. При этом необходимо рассмотреть возможные механизмы разрушения: упругопластическое предельное состояние либо квазихрупкий откол от свободной поверхности [12, 13].

Рис. 1.

Распределение давления на стенке трансформатора мощности вдоль горизонтальной плоскости на уровне т. А (см. рис. 14) для различных моментов времени: 1 – 0.1 мс, 2 – 0.2, 3 – 0.3, 4 – 0.4, 5 – 0.5, 6 – 0.75, 7 – 1, 8 – 1, 9 – 1.25, 10 – 1.5, 11 – 1.75, 12 – 2, 13 – 2.5.

Одной из наиболее эффективных систем автоматизированного инженерного анализа и нестационарных процессов деформирования является многоцелевой пакет LS-DYNA, предназначенный для решения трехмерных динамических нелинейных задач механики деформируемого твердого тела, механики жидкости и газа, теплопереноса, а также связанных задач [14]. Пакет LS-DYNA применяется для решения задач соударения тел, взрыва, обработки металлов давлением и ряда других проблем. В пакете реализованы эффективные методы решения перечисленных выше задач, в том числе явный и неявный метод конечных элементов, процедуры автоматической перестройки и сглаживания конечно-элементной сетки при вырождении элементов, произвольные лагранжево-эйлеровые сетки, высокоэффективные алгоритмы решения контактных задач, широкий набор моделей материалов, возможности пользовательского программирования, а также процедуры лагранжево-эйлерового связывания и расчета многокомпонентных течений сжимаемых сред на подвижных эйлеровых сетках. Программный код LS-DYNA оптимизирован под основные платформы и операционные системы, векторизован, распараллелен для систем с общей и распределенной памятью. Основной отличительной особенностью пакета, пока не реализованных в других пакетах, является возможность моделирования больших деформаций, включая разрушения объектов. На основании этого принято решение применить программу LS-DYNA для анализа отклика корпуса трансформатора на динамическое воздействие. К таким задачам относится моделирование процесса разрушения трансформатора при возникновении в его конструкции короткого замыкания. Моделирование проводилось для двух типов трансформаторов: трансформатора мощности и трансформатора тока. Несмотря на то, что конструктивно корпуса этих трансформаторов отличаются друг от друга, методика моделирования процессов изменения НДС и разрушения для них общая. Этапы расчета представлены на рис. 2.

Рис. 2.

Блок-схема расчета кинетики напряженно-деформированного состояния.

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

В связи с тем, что программа LS-DYNA по сути является только расчетной системой, в которой модули подготовки данных и анализа результатов расчета (препроцессор и постпроцессор) развиты слабо, то подготовка исходных данных, в частности, построение модели и задание граничных условий и анализ полученных результатов целесообразно проводить в системах с более развитыми модулями препроцессора и постпроцессора. Примером такой системы является программный комплекс ANSYS, в который программа LS-DYNA входит как дополнительный расчетный модуль [14].

Конечно-элементные модели корпусов трансформаторов (рис. 3) создавались в препроцессоре комплекса ANSYS. Давление на стенки корпуса рассчитывалось в комплексе моделирования поведения масла в условиях высоковольтного разряда, который использует метод конечных объемов. Для преобразования этого давления в формат конечно-элементной модели комплекса ANSYS создан специальный конвертор.

Рис. 3.

Конечно-элементные модели трансформаторов: (а) трансформатор мощности, (б) трансформатор тока.

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

Конечно-элементная модель трансформатора мощности состоит из Np = 183 455 узлов и Ne = = 170 874 оболочечных конечных элементов. Число уравнений составляет Neq = 1 100 730. Для трансформатора тока: Np = 81 351, Ne = 72 792 и Neq = 488 106.

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

${{\left( {\frac{{{{F}_{n}}}}{{\left[ {{{F}_{n}}} \right]}}} \right)}^{f}} + {{\left( {\frac{{{{F}_{\tau }}}}{{\left[ {{{F}_{\tau }}} \right]}}} \right)}^{s}} \geqslant 1,$

где Fn и Fτ – рассчитанные приходящиеся на связь нормальное растягивающее и сдвиговое усилия соответственно в расчете; [Fn] и [Fτ] – максимально допустимые нормальное и сдвиговое усилия связи при разрушении соответственно; f и s – показатели степени для нормального и сдвигового усилия. При задании нулевых значений [Fn] или [Fτ] для соответствующего направления связь переводится в неразрушающую. Величины [Fn] и [Fτ] определяются для каждого соединения (болтового или сварного) по известным моделям предельного состояния [15]. Условие разрушения связи по пластической деформации определяется из выражения ${{\varepsilon }_{p}} \geqslant [{{\varepsilon }_{p}}]$, где ${{\varepsilon }_{p}}$ – пластическая деформация в связи, определяемая по связанным элементам; [${{\varepsilon }_{p}}$] – максимально допустимая пластическая деформация в связи [15]. При задании нулевого значения максимальной пластической деформации $[{{\varepsilon }_{p}}]$ разрушение по пластической деформации не рассматривается.

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

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

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

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

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

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

На рис. 4 показаны модели трансформаторов с указанием связей.

Рис. 4.

Модели трансформатора с показанными связями: (а) трансформатор мощности, (б) трансформатор тока.

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

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

Модель материала. Если короткое замыкание в трансформаторе приводит к динамическому воздействию на его конструкцию, вызывающему деформирование ее элементов за пределами упругости, для описания характеристик материала необходимо применять кривые деформирования, зависящие от скорости деформации. Влияние скорости деформации на кривую деформирования конструкционных материалов исследовались различными авторами. Детальные обзоры можно найти в работах [5, 6, 12].

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

Рис. 5.

Кривые деформирования стали А36 при различных скоростях деформаций: 1 – 1.95 × 10–4, 2 – 1.06, 3 – 3.6 × 102 1/с.

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

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

Наиболее распространенной является модель Купера–Симонда. Однако для моделирования влияния скорости деформации на кривую деформирования можно применять более современную зависимость Джонсона–Кука [15, 16]:

${{\sigma }_{Y}} = [A + B{{(\varepsilon _{{{\text{eff}}}}^{p})}^{N}}]\left( {1 + C{\kern 1pt} {\text{ln}}\dot {\varepsilon }} \right)[1 - {{({{T}_{H}})}^{M}}],$
где $\varepsilon _{{{\text{eff}}}}^{p}$ – интенсивность пластических деформаций; A, B, C, N и M – константы материала;
$\dot {\varepsilon } = \frac{{\dot {\varepsilon }_{{{\text{eff}}}}^{p}}}{{{{{\dot {\varepsilon }}}_{0}}}},$
${{\dot {\varepsilon }}_{0}}$ – скорость деформации, при которой определялись A, B, N; ${{T}_{H}} = \frac{{T - {{T}_{R}}}}{{{{T}_{M}} - {{T}_{R}}}}$ – приведенная температура; TM – температура материала; TR – температура, при которой определялись A, B, N.

Данная модель требует определения пяти параметров материала A, B, C, N и M.

Сравнение ее с другими моделями для стали A36 показано на рис. 6 [13]. Нетрудно видеть, что модель Джонсона–Кука при ${{\dot {\varepsilon }}_{0}}$ = 1.0 с–1 в интересующим диапазоне скоростей деформирования приводит к результатам, незначительно отличающимся от более простой модели Купера–Симонда.

Рис. 6.

Сравнение результатов расчетов по моделям Купера–Симонда (1) и Джонма–Кука (2) с экспериментальными данными (3) для стали А36 при: (а) ${{\dot {\varepsilon }}_{0}}$ = 1.0 с–1, (б) ${{\dot {\varepsilon }}_{0}}$ = 10–4.

Поэтому в настоящей работе использована модель Купера–Симонда, имеющая вид

(1)
${{\sigma }_{{{\text{eff}}}}} = \left[ {1 + {{{\left( {\frac{{\dot {\varepsilon }}}{C}} \right)}}^{{\frac{1}{P}}}}} \right]\left( {{{\sigma }_{0}} + \beta {{E}_{P}}\varepsilon _{P}^{{{\text{eff}}}}} \right),$
где ${{\sigma }_{0}}$ – начальный предел текучести, $\dot {\varepsilon }$ – скорость деформации; C и P – параметры скорости деформации модели; $\varepsilon _{P}^{{{\text{eff}}}}$ – интенсивность пластической деформации; ${{E}_{P}}$ – модуль пластического упрочнения, определяемый по выражению
${{E}_{P}} = \frac{{{{E}_{{\tan }}}E}}{{E - {{E}_{{\tan }}}}},$
где E – модуль упругости, Etan – касательный модуль.

Свойства материала приведены в таблице.

Таблица 1.  

Свойства материала

Плотность ρ, кг/м3 7830
Модуль упругости E, ГПа 200
Коэффициент Пуассона 0.3
Предел текучести ${{\sigma }_{0}}$, МПа 340
Касательный модуль Etan, МПа 763
Параметр C 40
Параметр P 5
Максимальная пластическая деформация 0.25

Кривые деформирования, рассчитанные по зависимости (1) для различных скоростей деформации (с шагом 10 с–1), показаны на рис. 7.

Рис. 7.

Кривые деформирования материала при различных значениях скорости деформации: 1 – 10–5 1/с, 2 – 10–4, 3 – 10–3, 4 – 10–2, 5 – 10–1, 6 – 1, 7 – 10, 8 – 102, 9 – 103, 10 – 104.

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

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

(2)
$\left[ {M\left( u \right)} \right]\left\{ {\ddot {u}} \right\} + \left[ {C\left( u \right)} \right]\left\{ {\dot {u}} \right\} + \left\{ {{{F}^{{{\text{int}}}}}\left( u \right)} \right\} = \left\{ {{{F}^{{{\text{ext}}}}}\left( u \right)~} \right\},$
где $\left[ {M\left( u \right)} \right]$ и $\left[ {C\left( u \right)} \right]$ – зависящие от текущей конфигурации матрица масс и матрица демпфирования конструкции соответственно; $\left\{ {{{F}^{{\operatorname{int} }}}\left( u \right)} \right\}$ – зависящий от истории процесса деформирования вектор внутренних сил; $\left\{ {{{F}^{{{\text{ext}}}}}\left( u \right)} \right\}$ – вектор внешних сил; $\left\{ {\ddot {u}} \right\}$ – вектор ускорений точек конструкции; {$\dot {u}$} – вектор скоростей точек конструкции; $\left\{ u \right\}$ – вектор перемещений точек конструкции.

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

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

Зависимости матрицы $\left[ C \right]$ от вектора скоростей {$\dot {u}$} и вектора $\left\{ {{{F}^{{\operatorname{int} }}}\left( u \right)} \right\}$ от вектора перемещений $\left\{ u \right\}$ определяются принятой моделью неупругого поведения материала. Наиболее часто применяются модели, основанные на уравнениях состояния, которые связывают приращения напряжений с приращениями деформаций – модели типа “течения”.

Наиболее просто применить модели пластического течения, если использовать явную схему интегрирования, реализованную в модуле LS-DYNA.

Вектор ускорений в момент времени t определяется как

$\left\{ {{{{\ddot {u}}}_{t}}} \right\} = {{\left[ M \right]}^{{ - 1}}}\left( {\left\{ {F_{t}^{{{\text{ext}}}}} \right\} - \left\{ {F_{t}^{{{\text{int}}}}} \right\} - \left[ {C\left( u \right)} \right]\left\{ {\dot {u}} \right\}} \right),$
где $\left\{ {F_{t}^{{{\text{ext}}}}} \right\}$ – вектор приложенных внешних и объемных сил; $\left\{ {F_{t}^{{{\text{int}}}}} \right\}$ – вектор внутренних сил, который определяется как
${{F}^{{{\text{int}}}}} = \sum\limits_e {\left( {\int\limits_{{{\Omega }_{e}}} {{{{\left[ B \right]}}^{T}}\left\{ {{{\sigma }_{n}}} \right\}d\Omega } + {{F}_{e}}^{{hg}} + {{F}_{e}}^{{{\text{cont}}}}} \right)} .$
Здесь $\left[ B \right]$ – матрица, связывающая перемещения с деформациями, $\left\{ {{{\sigma }_{n}}} \right\}$ – вектор напряжений в элементе, ${{F}_{e}}^{{hg}}$ – сила сопротивления хрупкому разрушению элемента, ${{F}_{e}}^{{{\text{cont}}}}$ – контактные усилия в элементе.

Затем определяются векторы скоростей и перемещений как

$\left\{ {{{{\dot {u}}}_{{t + \Delta t}}}} \right\} = \left\{ {{{{\dot {u}}}_{{t - \Delta t}}}} \right\} + \left\{ {{{{\ddot {u}}}_{t}}} \right\}\Delta {{t}_{t}},$
$\left\{ {{{u}_{{t + {{\Delta t} \mathord{\left/ {\vphantom {{\Delta t} 2}} \right. \kern-0em} 2}}}}} \right\} = \left\{ {{{u}_{t}}} \right\} + \left\{ {{{{\dot {u}}}_{{t + {{\Delta t} \mathord{\left/ {\vphantom {{\Delta t} 2}} \right. \kern-0em} 2}}}}} \right\}\Delta {{t}_{{t + {{\Delta t} \mathord{\left/ {\vphantom {{\Delta t} 2}} \right. \kern-0em} 2}}}},$
где
$\Delta {{t}_{{t + {{\Delta t} \mathord{\left/ {\vphantom {{\Delta t} 2}} \right. \kern-0em} 2}}}} = 0.5\left( {\Delta {{t}_{t}} + \Delta {{t}_{{t + \Delta t}}}} \right),$
$\Delta {{t}_{{t - {{\Delta t} \mathord{\left/ {\vphantom {{\Delta t} 2}} \right. \kern-0em} 2}}}} = 0.5\left( {\Delta {{t}_{t}} - \Delta {{t}_{{t + \Delta t}}}} \right).$
Координаты узлов в новом положении для момента времени $t + \Delta t$ изменяются прибавлением приращения перемещений к начальной геометрии $\left\{ {{{x}_{0}}} \right\}$:
$\left\{ {{{x}_{{t + \Delta t}}}} \right\} = \left\{ {{{x}_{0}}} \right\} + \left\{ {{{u}_{{t + \Delta t}}}} \right\}.$
Напряженное состояние в каждом элементе определяются по зависимости
${{\left\{ {{{\sigma }_{e}}} \right\}}_{{t + \Delta t}}} = {{\left\{ {{{\sigma }_{e}}} \right\}}_{t}} + {{\left\{ {\Delta {{\sigma }_{e}}} \right\}}_{{\Delta t}}}.$
Приращения напряжений ${{\left\{ {\Delta {{\sigma }_{e}}} \right\}}_{{\Delta t}}}$ на каждом шаге определяются по зависимости
${{\left\{ {\Delta {{\sigma }_{e}}} \right\}}_{{\Delta t}}} = \left[ {{{D}_{\tau }}} \right]\left( {{{{\left\{ {\Delta \varepsilon } \right\}}}_{{\Delta t}}} - {{{\left\{ {\Delta {{\varepsilon }_{\theta }}} \right\}}}_{{\Delta t}}}} \right),$
где $\left[ {{{D}_{\tau }}} \right]$ – матрица связи приращений напряжений и приращений деформаций, определяемая выбранной моделью теории течения, ${{\left\{ {\Delta \varepsilon } \right\}}_{{\Delta t}}}$ = $ = {{\left[ B \right]}_{e}}{{\left\{ {\Delta {{u}_{e}}} \right\}}_{{\Delta t}}}$ – приращение деформаций в элементе, ${{\left\{ {\Delta {{\varepsilon }_{\theta }}} \right\}}_{{\Delta t}}}$ – приращение деформаций, вызванное изменением температуры в элементе, ${{\left\{ {\Delta {{u}_{e}}} \right\}}_{{\Delta t}}}$ – вектор приращений перемещений узлов элемента.

При решении нелинейных задач:

1. Требуется создание расчетной схемы, обеспечивающей обратимость матрицы сосредоточенных масс.

2. Уравнения становятся несвязанными и могут решаться напрямую (явно).

3. Не требуется построения и обращения матрицы жесткости. Все нелинейности, включая контакт, входят в вектор внутренних сил.

4. Главные вычислительные затраты приходятся на определение внутренних сил.

5. Не требуется контроля сходимости, так как уравнения не связаны.

6. Требуются очень маленькие шаги по времени для обеспечения устойчивости счета.

Явная схема решения устойчива, если только шаг по времени меньше критического шага $\Delta {{t}^{{{\text{crit}}}}}$:

$\Delta t \leqslant \Delta {{t}^{{{\text{crit}}}}} = \frac{2}{{{{\omega }_{{\max }}}}},$
где ${{\omega }_{{\max }}}$ – наибольшая собственная частота.

Размер критичного временного шага определяется по критерию Куранта–Фридрихса–Леви следующим образом.

Собственная частота определяется как

${{\omega }_{{\max }}} = 2\frac{c}{l},$
где c – скорость распространения волны и определяется из выражения
$c = \sqrt {\frac{E}{\rho }} ,$
l – характерный размер элемента.

Подстановка ${{\omega }_{{\max }}}$ в уравнение для размера критического временнóго шага дает

$\Delta t = \frac{l}{c}.$
Таким образом, $\Delta t$ – время необходимое для прохождения волны через элемент длиной l, т.е. размер критического шага по времени для явной схемы интегрирования зависит от размеров элементов и свойств материала (скорости звука).

Система LS-DYNA проверяет все элементы и затем вычисляет необходимый шаг по времени. Для увеличения устойчивости расчета используется коэффициент 0.9 для уменьшения шага по времени

$\Delta t = 0.9\frac{l}{c}.$

Для четырехугольных оболочечных конечных элементов характерная длина определяется как

$l = \frac{S}{{\max \left( {{{l}_{1}},{{l}_{2}},{{l}_{3}},{{l}_{4}}} \right)}},$
где S – площадь элемента, l1, l2, l3, l4 – характерные размеры.

Для треугольных оболочечных конечных элементов

$l = \frac{{2S}}{{\max \left( {{{l}_{1}},{{l}_{2}},{{l}_{3}}} \right)}}.$

Для всех типов оболочечных конечных элементов

$c = \sqrt {\frac{E}{{\rho (1 - {{\mu }^{2}})}}} .$

Из-за очень малого шага по времени явная схема пригодна только для очень быстрых переходных процессов, к которым относится рассматриваемая проблема.

РЕЗУЛЬТАТЫ РАСЧЕТОВ УДАРНОГО ВОЗДЕЙСТВИЯ

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

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

На рис. 8–10 показаны распределения суммарных перемещений деталей, интенсивности эквивалентных напряжений в различные моменты времени для описанных вариантов расчета. Эквивалентные напряжения соответствуют вторым инвариантам тензоров напряжений.

Рис. 8.

Распределение суммарных перемещений (мм) в момент времени t = 0.025 с для различных вариантов расчета: (а) вариант 1, (б) вариант 2, (в) вариант 3.

Рис. 9.

Распределение суммарных перемещений (мм) в момент времени t = 0.05 с для различных вариантов расчета: (а)–(в) – см. рис. 8.

Рис. 10.

Распределение интенсивности эквивалентных напряжений (Па) в момент времени t = 0.025 с для различных вариантов расчета: (а)–(в) – см. рис. 8.

Представленные на рис. 8 и 9 результаты расчета показывают, что учет связей в соединениях существенно влияет на поведение модели. Так, в первом варианте корпус не разрушается, так как максимальные эквивалентные напряжения (рис. 10) не превосходят предела текучести, а соответствующие им максимальные деформации не превосходят величины максимальных пластических деформаций. Характер деформирования корпуса в моменты времени t = 0.025 с и t = 0.05 с показывают, что стенка трансформатора колеблется. При этом максимальные пластические деформации возрастают.

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

УДАРНОЕ ВОЗДЕЙСТВИЕ ПРИ КОРОТКОМ ЗАМЫКАНИИ

Моделирование заданного воздействия при коротком замыкании производилось для трансформаторов мощности и тока (рис. 3).

При моделировании трансформатора мощности была использована модель, соответствующая третьему варианту из предыдущего параграфа. Свойства материала определялись по зависимости (1).

Основным отличием являлось то, что при коротком замыкании скорость приложения давления к стенкам корпуса трансформатора ниже, чем при моделировании взрыва.

На рис. 11 показано изменение давления в характерных точках A, B и C корпуса трансформатора мощности, где красным ромбом обозначен момент начала разрушения.

Рис. 11.

Изменение давлений и напряжений в точке A корпуса трансформатора мощности в зависимости от времени.

Давление определено расчетным путем по газодинамическим расчетам. По полученным результатам можно сделать вывод о том, что максимальное давление ∼55 МПа в точке A возникает через 0.00025 с после короткого замыкания. Давление действует на стенки корпуса с определенной задержкой в точках B и C. При этом максимальные величины давлений в этих точках в три раза меньше, чем в точке A. Это приводит к тому, что разрушение развивается вблизи точки A.

На рис. 11 показаны графики изменения давления p и эквивалентного напряжения в точке A в зависимости от времени от момента короткого замыкания до момента выполнения условия разрушения.

Момент разрушения определялся при достижении эквивалентной пластической деформации предельного значения εplmax = 0.25. Резкий рост напряжений в точке A связан с особенностями высокоскоростного деформирования. Скорость деформирования в точке A показана на рис. 12, где приведены для сравнения скорости деформирования в точках B и C.

Рис. 12.

Динамические кривые деформирования для точек A, B, и C корпуса трансформатора мощности.

Отличие в скоростях деформирования и условиях приложенных давлений приводит к различию в уровнях деформаций и характере изменения напряжений в точках A, B и C (рис. 12). Нетрудно видеть, что предельные деформации возникают в точке A, вблизи которой и разрушается корпус.

На рис. 13, 14 показаны изменения давления и интенсивности эквивалентных напряжений и интенсивности пластических деформаций в корпусе трансформатора мощности в различные моменты времени.

Рис. 13.

Изменение давления в характерных точках конструкции трансформатора тока в зависимости от времени; на вставке – положения точек A, B, C, D, E и F.

Рис. 14.

Изменение эквивалентных напряжений в характерных точках корпуса трансформатора тока в зависимости от времени.

Разрыв корпуса происходит в зоне взрыва за счет отрыва панели вертикальной стенки у дна корпуса.

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

На рис. 15 показано распределение суммарных перемещений корпуса трансформатора тока в различные моменты времени. Различные условия приложенных давлений приводят к различию в скоростях деформирования и, соответственно, уровнях деформаций и характере изменения эквивалентных напряжений в характерных точках корпуса трансформатора тока.

Рис. 15.

Распределение суммарных перемещений корпуса трансформатора тока в различные моменты времени: (а) 0.0005 с, (б) 0.001, (в) 0.0015.

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

На рис. 15 также показаны распределения суммарных перемещений, интенсивности эквивалентных напряжений и интенсивности пластических деформаций в корпусе трансформатора тока в различные моменты времени соответственно.

Как видно из представленных рисунков, происходит разрыв корпуса масляного бака.

Характер изменения скоростей деформирования, эквивалентных напряжений и пластических деформаций в точке разрушения (точка D на рис. 13) корпуса трансформатора тока аналогичен характеру изменений в точке A корпуса трансформатора мощности.

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

ЗАКЛЮЧЕНИЕ

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

Работа выполнена при поддержке гранта РНФ № 19-49-02031.

Список литературы

  1. Сон Э.Е., Фортов В.Е., Азметов Х.Х., Бондарь В.С., Бабаева Н.Ю., Гаджиев М.Х., Зицерман В.Ю., Киверин А.Д., Коробейников С.М., Куликов Ю.М, Минцев В.Б., Панов В.А., Смирнов Б.М., Темис Ю.М., Тюфтяев А.С. Электрофизика и взрывобезопасность высоковольтного маслонаполненного электрооборудования. М.: Изд-во РАН, 2019. 417 с.

  2. Fluid-structure Interaction: Modelling, Simulation, Optimization. Eds. Bungartz H.; Schäfer M. Springer, 2006.

  3. Son E.E., Tarasova E.N., Zubkov P.T. Fully Developed Two-phase Liquid-liquid Flow in Finned Duct. MIT Conference Fluid-Structure Interaction, 2005.

  4. Sigrist J.F. Fluid-Structure Interaction: An Introduction to Finite Element Coupling. Wiley, 2015.

  5. Heil M. An Efficient Solver for the Fully Coupled Solution of Large-displacement Fluid-structure Interaction Problems // Computer Methods in Applied Mechanics and Engineering. 2004. V. 193. № 1–2. https://doi.org/10.1016/j.cma.2003.09.006

  6. Bathe K.-J., Zhang H. Finite Element Developments for General Fluid Flows with Structural Interactions // Int. J. Numerical Methods in Engineering. 2004. V. 60. № 1. P. 213.

  7. Kulikov Y.M., Son E.E. Stability of Thermoviscous Fluid Flow Under High Temperature Gradients // High Temp. 2017. V. 55. № 1. P. 131.

  8. Аксенов А.А., Жлуктов С.В., Кудимов Н.Ф., Сон Э.Е., Таран М.Д., Третьякова О.Н., Шишаева А.С. О моделировании сложного теплообмена в силовых трансформаторах большой мощности // Известия РАН. Серия Энергетика. 2013. № 2. С. 131.

  9. Фортов В.Е., Сон Э.Е., Зибаров А.В., Бондарь В.С., Темис Ю.М., Дарьян Л.А., Долуденко А.Н., Горюшин Ю.А., Дементьев Ю.А., Иванов М.Ф. Моделирование взрыва высоковольтных маслонаполненных трансформаторов // Известия РАН. Серия Энергетика. 2011. № 5. С. 64.

  10. Дарьян Л.А., Козлов А.В., Полищук В.П., Поварешкин М.Н., Сон Э.Е., Фортов В.Е., Шурупов А.В. Бездуговой метод испытания высоковольтного маслонаполненного оборудования на взрывобезопасность // Известия РАН. Серия Энергетика. 2011. № 5. С. 74.

  11. Бердников Р.Н., Фортов В.Е., Сон Э.Е., Горюшин Ю.А., Дементьев Ю.А., Зибаров А.В. Программный комплекс моделирования взрыва высоковольтного маслонаполненного электрооборудования в закрытых камерах трансформаторов ПС (ExploTraP). 2 011 617 913 А.с. 2 012 611 186.

  12. Трощенко В.Т., Лебедев А.А., Стрижало В.А. и др. Механическое поведение материалов при различных видах нагружения. Киев: Логос, 2000. 571 с.

  13. Прочность материалов и конструкций / Под ред. Трощенко В.Т. (отв. ред.) и др. Киев: Академпериодика, 2005. 1088 с.

  14. ANSYS LS-DYNA User’s Guide. SAS IP, Inc., 2009.

  15. Бондарь В.С. Неупругость. Варианты теории. Физматлит, 2004. 186 с.

  16. Schwer D. Optional Strain-rate Forms for Johnson–Cook Constitutive Model and the Role of Parameter epsilon_01 // The 6th European LS-DYNA Conference, 2007.

Дополнительные материалы отсутствуют.