CC BY
11ВЕСТН. САМАР. ГОС. ТЕХН. УН-ТА. СЕР. ТЕХНИЧЕСКИЕ НАУКИ. 2018. № 3 (59)
УДК 621.365
ПАРАМЕТРИЧЕСКИЙ АНАЛИЗ И ОПТИМАЛЬНОЕ ПРОЕКТИРОВАНИЕ ИНДУКЦИОННОЙ СИСТЕМЫ ПО КРИТЕРИЮ МАКСИМАЛЬНОГО КОЭФФИЦИЕНТА ПОЛЕЗНОГО ДЕЙСТВИЯ
А.И. Данилушкин, В.А. Данилушкин, Д.А. Животягин
Самарский государственный технический университет Россия, 443100, г. Самара, ул. Молодогвардейская, 244
E-mail: aidanilushkin@mail.ru; vasiliydan2013@yandex.ru; Denis.Zhivotyagin@arconic.com
Аннотация. Рассмотрены вопросы оптимального проектирования трехфазного индуктора для нагрева крупногабаритных цилиндрических заготовок. Отмечены специфические особенности рассматриваемых в работе индукционных нагревателей. Нагрев осуществляется в поперечном магнитном поле. Это обстоятельство существенно сказывается на распределении источников тепла в нагреваемой заготовке. Наличие источников тепла в торцевых частях заготовки позволяет компенсировать тепловые потери и повысить равномерность температурного распределения по аксиальной координате. Решена одна частная задача оптимального проектирования по критерию максимального общего коэффициента полезного действия индуктора. Поиск оптимальных параметров проводится на основе двумерной модели. В качестве варьируемых параметров рассматриваются материал и толщина футеровки, количество пазов и зубцов статора и соотношение их размеров. Другие параметры, оказывающие существенное влияние на величину к.п.д., рассматриваются как ограничения. Анализ влияния перечисленных выше параметров на общий коэффициент полезного действия установки основан на многократном уточнении результатов в процессе итерационного проектирования.
Ключевые слова: индукционный нагрев, магнитная проницаемость, плотность тока, температура, коэффициент полезного действия, оптимальное проектирование.
В технологических комплексах обработки металла на деформирующем оборудовании нагрев крупногабаритных заготовок производится в индукционных установках промышленной частоты. Распространенные конструкции индукционных нагревателей большой мощности, функционирующих на промышленной частоте, представляют собой систему индукторов, питающихся от трехфазной сети [1, 2]. Количество индукторов зависит от производительности комплекса и не всегда кратно трем. По этой причине, а также вследствие взаимной индуктивности катушек индукторов возникает существенная несимметрия трехфазной нагрузки, несмотря на равенство собственных сопротивлений отдельных секций многосекционного нагревателя. Следствием этого является неравномерная загрузка фаз питающей сети, превышающая установленные нормы несимметрии сети [3, 4, 5]. Снизить несимметрию напряжений можно при помощи специаль-
Александр Иванович Данилушкин (д.т.н., проф.), профессор кафедры «Электроснабжение промышленных предприятий».
Василий Александрович Данилушкин (к.т.н.), доцент кафедры «Электроснабжение промышленных предприятий».
Денис Александрович Животягин, аспирант.
ных симметрирующих устройств. В то же время их применение приводит к завышенной стоимости, так как стоимость и массогабаритные показатели симметрирующих устройств сравнимы с аналогичными показателями индукционной установки. Кроме того, в условиях изменяющихся при нагреве параметров нагрузки системы «индуктор - металл» необходимо использовать регулируемые симметрирующие устройства. Это еще более повышает капитальные затраты и габариты индукционной установки [5, 6, 7, 8]. Указанное обстоятельство приводит к необходимости совершенствования и создания новых модификаций индукционных нагревательных установок с целью снижения капитальных и эксплуатационных затрат. В работе в качестве альтернативного варианта рассматривается трехфазный индукционный нагреватель для сквозного нагрева крупногабаритных цилиндрических заготовок, отличие которого заключается в конструктивном исполнении индуктора с замкнутым магнитопроводом в форме статора трехфазного асинхронного электродвигателя [9, 10]. При симметричном расположении нагреваемого цилиндра и симметричной трехфазной системе обмоток индуктора несимметрия в системе «индуктор - загрузка» отсутствует. Кроме того, наличие замкнутого магнитопровода позволяет уменьшить величину реактивной мощности и количество конденсаторов для ее компенсации.
Для расчета параметров индукционной системы необходимо решить ряд задач, связанных с разработкой математической модели и исследованием взаимосвязанных электротепловых процессов в системе, содержащей магнитопровод, трехфазную обмотку индуктора, футеровку и нагреваемый цилиндр. Применение известных моделей и инженерных методик расчета для цилиндрических индукционных нагревателей продольного магнитного поля для расчета параметров предлагаемой конструкции не представляется возможным в силу ряда принципиальных конструктивных особенностей. Геометрическая модель исследуемой индукционной системы представлена на рис. 1.
Рис. 1. Трехфазный индуктор с магнитопроводом:
1 - магнитопровод; 2 - водоохлаждаемый виток обмотки; 3 - изоляция витка; 4 - футеровка; 5 - загрузка; 6 - воздушный зазор
Систему «трехфазный индуктор - металлический цилиндр» можно с определенным приближением рассматривать как асинхронный двигатель с короткоза-мкнутым ротором в режиме короткого замыкания [9]. Однако исследуемый объект имеет существенные отличия, заключающиеся в наличии футеровки, увеличивающей зазор между статором и ротором (индуктором и нагреваемым цилиндром), что нехарактерно для асинхронных машин. Кроме того, в процессе нагрева в широких пределах изменяется удельное сопротивление нагреваемого цилиндра. В связи с указанными обстоятельствами для исследования процессов нагрева и расчета конструктивных и режимных параметров оптимальной конструкции трехфазного индукционного нагревателя с замкнутым магнитопроводом предложена методика, включающая численное моделирование электромагнитных и тепловых полей в системе «трехфазный индуктор - цилиндрическая заготовка».
Математическая модель процесса индукционного нагрева цилиндрической заготовки в поперечном магнитном поле может быть представлена системой взаимосвязанных уравнений [10, 11]:
- для электромагнитной задачи
д(Ца (и (г, х, т), т )и (г, X, т)) 1 д ( 1 ди(г, х, т) ^ г дг ^ а(т) дг
д/ дх
1 дИ (г, х, т)' а(Т) дх ,
= 0;
(1)
с краевыми условиями
дИ (0, х, т)
дг
= 0:
И (я, Ь, т) = и(т)/ (и, т);
- для тепловой задачи
с(т Ит (г,е,х,т)=IАГ г^(т )дт ге,х,т) |+
дт г дг ( дг
+1 К* ^ Н (^ ^ 1+- (т);
с краевыми условиями
(2)
= 0; Мг)^^ = (т,т(я, х, т));
дг дг
,(т = х{т = вп (т, т (г,0, т));
дх дх
т(г,х,0) = т,(г,х); ге[0,я]; хе [0,Ь].
Здесь И (г, х, т) - напряженность магнитного поля; а(т) - удельная электропроводимость; - магнитная проницаемость; г, х - радиальная и осевая координаты объекта; т - время; т(г, х, е, т) - температура; е - угловая координата;
Х(Т) - коэффициент теплопроводности; ^ (т) - удельная мощность тепловыделения; с(Т) - удельная теплоемкость (в нелинейном случае это функция температуры, аппроксимированная кубическими сплайнами); у - плотность, - тепловой поток с торцевых поверхностей заготовки.
Для решения задачи оптимального проектирования используются методы структурной и параметрической оптимизации технических объектов, которые позволяют рассматривать задачи оптимизации как многокритериальные [12, 13, 14]. На практике часто многокритериальная задача сводится к поиску одного наиболее важного частного критерия, а остальные конкурирующие критерии рассматриваются как ограничения.
В настоящей работе в качестве частного критерия рассматривается коэффициент полезного действия индукционной установки.
Задача оптимального проектирования индукционного нагревателя по критерию максимального коэффициента полезного действия формулируется следующим образом.
Для объекта, описываемого системой уравнений (1), (2) с соответствующими начальными и граничными условиями, необходимо найти такие значения величин конструктивных и режимных параметров системы индукционного нагрева, которые за заданное время обеспечивают максимальное значение общего коэффициента полезного действия при выполнении требуемого по технологии конечного распределения температуры по объему нагреваемого цилиндра
Л = Лх -Л2 (3)
в условиях существующих энергетических и технологических ограничений.
Здесь ^ - электрический КПД установки [15]:
Пх = ^¡(Д, А, Рх, р2, *1, Х2, аъ а2); (4)
- тепловой КПД:
^2 = , аф, Ч(Т ))ф; (5)
рх, р2 (Т) - соответственно удельные сопротивления материала обмотки индуктора и металла заготовки; Х1 - ширина паза; X2 - ширина зубца; к - высота паза; Д, Д2 - глубина проникновения тока соответственно в металл индуктора и металл заготовки; а , а - эффективная длина индуктора и заготовки соответственно; А, А - диаметры индуктора и заготовки; - толщина футеровки; -коэффициент теплообмена между футеровкой и индуктором; Хф (Т) - коэффициент теплопроводности материала футеровки.
Энергетические и технологические ограничения принимают вид:
- удельная поверхностная мощность р0 <400кВт/ м ;
- максимальная температура наиболее нагретой точки по сечению заготовки Ттх (г, X, тКОН )< 480° С;
- перепад температур по сечению заготовки на выходе из нагревателя Ттх -Тц)= 15°С ;
- максимальная индукция в магнитопроводе Втах = 0,8 Тл;
г <60°С Т
пов.маги ' магмах
- температура магнитопровода Тпов.маги <60°С , Тмаг^ = 90°С ;
Варьируемые параметры:
1. Количество пазов и зубцов индуктора (определяет степень неравномерности температур под пазом и зубцом и величину потока рассеяния в пазу).
2. Соотношение ширины Х^ паза и Х2 зубца (максимальное значение магнитной индукции, максимальное сечение провода индуктора).
3. Плотность 5 тока в индукторе.
4. Параметры футеровки: толщина стенки (Sф), коэффициент Хф(Т) теплопроводности материала футеровки.
Результатом проектирования являются материалы магнитопровода и футеровки, геометрические размеры индуктора и его элементов, обеспечивающие в совокупности максимальный коэффициент полезного действия. К геометрическим размерам, определяющим общий коэффициент полезного действия, относятся прежде всего материал и толщина футеровки, число пазов, ширина паза и зубца магнитопровода. Анализ влияния перечисленных выше параметров на общий коэффициент полезного действия установки требует многократного уточнения результатов в процессе проектирования.
В качестве примера для исследования рассматривается индукционная установка для нагрева цилиндрической заготовки из алюминиевого сплава марки АД31. Размеры заготовки: диаметр 0,42 м, длина 0,8 м, масса 299,1 кг. Температура максимально нагретой точки 460 0С, температура центра 430 0С . Допустимое отклонение температуры ±10 0С . Напряжение питающей сети 380/220 В; частота 50 Гц. Максимальное значение индукции в магнитопроводе Втах находится в пределах 1,19-1,33 Тл.
Физические свойства сплава: плотность 2700 кг / М ; температура плавления 660 0С; коэффициент теплопроводности Х = 228 Вт/(м-0С) ; коэффициент теплоемкости Ср = 880 Дж/(кг-0С) ; удельное сопротивление р = 2.7 • 10-8 Ом• м. Поиск проводится на основе двумерной электромагнитной задачи. Перейдя от системы уравнений Максвелла [10] к формулировке с использованием векторного магнитного потенциала, можно записать дифференциальное уравнение с соответствующими граничными условиями:
д_
дх
( 1
1 дА
^ дх
дУ 1ц х дУ
- /шаА = -5^ ; А| = 0. (6)
Здесь А - векторный магнитный потенциал; цх, ц - относительная магнитная проницаемость материала по осям X, у; ] - мнимая единица; ю - круговая частота тока; Ф - удельная электрическая проводимость; 5^ - плотность стороннего тока. На базе уравнений (1), (2) строится конечно-элементная формулировка для плоской двумерной области [16]. В качестве программной реализации используется пакет Б1сШ; [17]. Конечная система уравнений имеет вид
К ]-М={/}. (7)
Рис. 2. Алгоритм расчета параметров индукционного нагревателя
Здесь \К ] {/}, и} - матрица жесткости и векторы источников и неизвестных величин. Геометрическая модель тепловой задачи полностью повторяет модель электромагнитной задачи, источником нагрева является мощность внутренних источников тепла в цилиндре. Граничные условия, которые задаются на внешних поверхностях, позволяют учесть тепловые потери с поверхности цилиндра.
Алгоритм расчета оптимальных по критерию максимального общего коэффициента полезного действия параметров индукционного нагревателя представлен на рис. 2.
В отличие от нагрева цилиндра в продольном магнитном поле соленоидаль-ного индуктора в рассматриваемом индукторе нагрев осуществляется в поперечном магнитном поле, причем чередование активных проводников с током в пазу индуктора с зубцами магнитопровода приводит к неравномерному распределению плотности тока в нагреваемой заготовке по окружности. Картина распределения плотности тока в цилиндре по радиальной координате под проводником и зубцом магнитопровода приведена на рис. 3. Эта неравномерность отражается на температурном распределении в цилиндре, однако, как показывают исследования, разность температур цилиндра на окружностях одного радиуса наблюдается только в пределах глубины проникновения тока в цилиндр и не вносит существенных изменений в общую картину распределения температуры по сечению.
радиус заготовки, мм область под пазом ^♦—область под зубцом
Рис. 3. Распределение плотности полного тока по радиальной координате
В таблице представлены параметры индуктора при различном соотношении размеров паза и зубца магнитопровода. Как следует из анализа полученных результатов, наиболее существенное влияние на общий коэффициент полезного действия оказывают количество пазов и зубцов магнитопровода, соотношение их размеров, материал и толщина тепловой изоляции. Наилучшие энергетические показатели индуктора соответствуют позиции 4 (см. таблицу).
Заключительным этапом расчета является определение интегральных параметров системы - числа витков и сечения проводов катушки индуктора, напря-
жения на индукторе в соответствии с выбранным напряжением источника питания. Материал футеровки - шамот. Проводник катушки индуктора: медная шинка сечением 4*12 мм; КПД 49,7 %; соБф = 0,29; магнитная индукция в магнитопроводе 0,8 Тл; мощность полезная 150,8 кВт; мощность полная 251,8 кВт; мощность индуктора 868 кВА.
Результаты проведенного анализа - максимальное значение индукции маг-нитопровода, графики зависимости мощности тепловыделения в заготовке, коэффициента полезного действия, коэффициента мощности от варьируемых параметров - позволяют сделать вывод о степени влияния геометрических параметров на энергетические характеристики. В исследуемой установке значения указанных параметров определяются прежде всего характеристиками материала магнитопровода и футеровки, толщиной тепловой изоляции, количеством пазом и соотношением ширины паза и зубца магнитопровода при фиксированной толщине футеровки.
Параметры индуктора при различном соотношении размеров паза и зубца
Параметр № позиции
1 2 3 4 5
Внутренний диаметр индуктора, мм 480 480 480 480 480
Внешний диаметр магнитопровода, мм 630 630 630 630 630
Толщина футеровки, мм 20 20 20 20 20
Ширина паза, мм 60 65 70 75 80
Ширина зубца, мм 63 58 53 48 43
Мощность в загрузке, кВт 152 150 149,8 150,8 149,8
Потери мощности, кВт 128 118 110 101 136
Мощность индуктора, кВт 280 268 259,8 251,8 285,8
Полная мощность, кВА 1037 1004 951 868 1039
КПД 0,42 0,447 0,453 0,497 0,428
Cos ф 0,27 0,267 0,273 029 0,275
Выводы
1. Эффективность установок индукционного нагрева определяется при прочих равных условиях предельно допустимым значением плотности тока в обмотке индуктора. В установках высокотемпературного индукционного нагрева максимальная плотность тока при принудительном водяном охлаждении индуктора определяется расходом воды и ее температурой на входе и составляет 15^25 A/мм2. В предлагаемой конструкции индуктора система охлаждения должна дополнительно обеспечивать охлаждение магнитопровода индуктора, поэтому максимальное значение плотности тока желательно поддерживать на нижнем уровне.
2. При нагреве цилиндрической заготовки в поперечном магнитном поле трехфазного индуктора кроме неравномерности температурного распределения по радиусу заготовки появляется дополнительная неравномерность температурного распределения по окружности заготовки. Эта неравномерность зависит от 150
числа зубцов и пазов статора и соотношения их размеров. Для снижения неравномерности желательно увеличивать их число. Однако увеличение числа пазов может привести к существенному уменьшению сечения зубца магнитопровода, а следовательно, к его насыщению. В связи с этим выбор количества пазов и зубцов должен производиться с учетом дополнительного требования к размерам сечения.
БИБЛИОГРАФИЧЕСКИЙ СПИСОК
1. Немков В.С., Демидович В.Б. Теория и расчет устройств индукционного нагрева. - Л.: Энер-гоатомиздат, 1988. - 280 с.
2. Электрооборудование и автоматика электротермических установок: Справочник / Под ред. А.П. Альтгаузена. - М.: Энергия, 1978. - 304 с.
3. Иванов B.C. Режимы потребления и качество электроэнергии систем электроснабжения промышленных предприятий. - М.: Энергоатомиздат, 1987. - 336 с.
4. Наумов И.В. Несимметричная нагрузка в сети 0,38 кВ и потери мощности // Электротехника (объединенный научный журнал). - 2002.- № 2 . - С. 50-52.
5. Аввакумов В.Г. Уравновешивание электрической нагрузки в трехфазной четырехпроводной системе // Известия вузов. Энергетика. - 1978. - № 5. - С. 94-99.
6. Наумов И.В. Методика расчета показателей несимметрии токов и напряжений в сети 0,38 кВ с симметрирующим устройством // Вестник АлтГАУ им И.И. Ползунова. - Барнаул: АлтГАУ, 2001. - № 2. - С. 49-56.
7. Современные энергосберегающие технологии: учеб. пособие для вузов / Ю.И. Блинов [и др.] - СПб.: Изд-во СПбГЭТУ «ЛЭТИ», 2000. - 564 с.
8. Немков В.С., Демидович В.Б., Руднев В.И. Влияние конструкции и режимов работы индукционных нагревателей на их энергетические показатели // Электротехника. - М., 1986. - № 3. -С. 23-27.
9. Данилушкин А.И., Никитина Е.А., Бойков Е.М. Повышение качества электроэнергии в системах электроснабжения кузнечных индукционных нагревателей промышленной частоты // Известия вузов. Электромеханика. Специальный выпуск «Электроснабжение». - Новочеркасск: ЮРГТУ (НПИ), 2009. - С. 35-36. ISSN 0136-3360.
10. Базаров А.А., Данилушкин А.И., Никитина Е.А. Моделирование и расчет внутренних источников тепла в трехфазном индукторе с вращающимся магнитным полем // Вестник Самарского государственного технического университета. Сер. Технические науки. - 2009. - Вып. 2(24). -С. 120-127.
11. ВайнбергА.М. Индукционные плавильные печи. - М.: Энергия, 1967. - 415 с.
12. ЛыковА.В. Теория теплопроводности. - М.: Высшая школа, 1967. - 599 с.
13. Уайлд Д. Оптимальное проектирование. - М.: Мир, 1981. - 272 с.
14. Рапопорт Э.Я. Оптимизация процессов индукционного нагрева металла. - М.: Металлургия, 1993. - 279 с.
15. Анищенко Л.М., Лавренюк С.Ю. Математические основы проектирования высокотемпературных технологических процессов. - М.: Наука, 1986. - 78 с.
16. Бате К., Вилсон Е. Численные методы анализа и метод конечных элементов. - М.: Стройиз-дат, 1982. - 448 с.
17. ELCUT. Моделирование двумерных полей методом конечных элементов. Руководство пользователя. Версия 5.9. - СПб.: Производственный кооператив ТОР, 2009.
Статья поступила в редакцию 10 июня 2018 г.
PARAMETRIC ANALYSIS AND OPTIMAL DESIGN OF INDUCTION SYSTEM WITH RESPECT TO THE MAXIMUM EFFICIENCY CRITERIA
A.I. Danilushkin, V.A. Danilushkin, D.A. Zhivotyagin
Samara State Technical University
244, Molodogvardeyskaya st., Samara, 443100, Russian Federation
Abstract. The problems of optimal design of a three-phase inductor for heating large-sized cylindrical billets are considered. The specific features of the induction heaters considered in the work are noted. Heating isrealized in a transverse magnetic field. This circumstance significantly affects the distribution of heat sources in the heated billet. The presence of heat sources in the end parts of the billet makes it possible to compensate the thermal losses and increase the uniformity of the temperature distribution along the axial coordinate. One particular problem of optimal design was solved by the criterion of the maximum overall efficiency of the inductor. The search for optimal parameters is based on a two-dimensional model. As variable parameters, the material and thickness of the lining, the number of slots and teeth of the stator and the ratio of their dimensions are considered. Other parameters that have a significant effect on the efficiency value are considered as limitations. Analysis of the influence of the above parameters on the overall efficiency of the installation is based on multiple refinement of the results in the iterative designing process.
Keywords: induction heating, magnetic permeability, current density, temperature, coefficient of efficiency, optimal design.
REFERENCES
1. Nemkov V.S., Demidovich V.B. Teoriya i raschet ustrojstv indukcionnogo nagreva [Theory and calculation of induction heating devices]. Leningrad, Energoatomizdat, 1988. 280 pp. (In Russian).
2. Elektrooborudovanie i avtomatika elektrotermicheskih ustanovok: Spravochnik [Electrical equipment and automation of electrothermal processes: Reference book]. Pod red. A. P. Al'tgauzena. -Moscow. Energiya, 1978. - 304 pp. (In Russian).
3. Ivanov B.C. Rezhimy potrebleniya i kachestvo elektroenergii sistem elektrosnabzheniya promysh-lennyh predpriyatij [Consumption modes and quality of electric power of power supply systems of industrial enterprises]. Moscow, Energoatomizdat, 1987. 336 pp. (In Russian).
4. Naumov I.V. Unbalanced load in the 0,38 kV electric grid and power loss// Elektrotekhnika (Ob"edinennyjnauchnyjzhurnal). 2002. No 2. Pp. 50-52. (In Russian).
5. Avvakumov V.G. Balancing the electrical load in a three-phase four-wire system // Izv. Vuzov. Ener-getika. 1978. No5. Pp. 94-99. (In Russian).
6. Naumov I. V. Method for calculating the unbalance of currents and voltages in a 0.38 kV electric grid with a balancing device // Vestn. AltGAU im I.I. Polzunova. Barnaul. AltGAU, 2001. No2. Pp. 4956. (In Russian).
7. Sovremennye energosberegayushchie tekhnologii: uchebnoe posobie dlya VUZov [Modern energy-saving technologies: a textbook for high schools] / Yu.I. Blinov [i dr.] Saint Petersburg. Izd-vo SPbGETU «LETI», 2000. 564 pp. (In Russian).
8. Nemkov V.S., Demidovich V.B., Rudnev V.I. Influence of design and operation modes of induction
Aleksandr I. Danilushkin (Dr. Sci. (Techn.)), Professor. Vasili A. Danilushkin (Ph.D. (Techn.)), Associate Professor Denis A. Zhivotyagin, Postgraduate Student.
heaters on their energy performances // Elektrotekhnika. 1986. No 3. - Pp. 23-27. (In Russian).
9. Danilushkin A.I., Nikitina E.A., Bojkov E.M. Improving the quality of electricity in power supply systems for forging induction heaters of industrial frequency // Elektromekhanika. Izvestiya VUZov. Special'nyj vypusk «Elektrosnabzhenie». Novocherkassk. YURGTU (NPI). 2009. Pp. 35-36. (In Russian).
10. Bazarov A.A., Danilushkin A.I., Nikitina E.A. Modeling and calculation of internal heat sources in a three-phase inductor with a rotating magnetic field // Vestn. Samar. Gos. Tekhn. Un-ta. Ser. Tekhn. Nauki. 2009. No2(24). Pp. 120-127. (In Russian).
11. Vajnberg A.M. Indukcionnye plavil'nye pechi [Induction Melting Furnaces]. Moscow, Energiya, 1967. 415 pp. (In Russian).
12. LykovA. V. Theory of heat conductivity. Vysshaya shkola, Moskva, 1967, 599 pp. (In Russian).
13. Uajld, D. Optimal design. Moscow. Mir. 1981. 272 pp. (In Russian).
14. RapoportE.Ya. Optimization of induction heating of metal. - Moscow. Metallurgiya. 1993. 279 pp. (In Russian).
15. Anishchenko L.M., Lavrenyuk S.Yu. Mathematical backgrounds of designing high-temperature technological processes. Moscow. Nauka. 1986. 78 pp. (In Russian).
16. Bate K., Vilson E. Numerical analysis methods and finite element method. Moscow. Strojizdat, 1982. 448 pp. (In Russian).
17. ELCUT. Simulation of two-dimensional fields by the finite element method. User guide. Version 5.9. Saint Petersburg. Proizvodstvennyj kooperativ TOR. 2009.
