РЕШЕНИЕ ОБРАТНОЙ ЗАДАЧИ НАХОЖДЕНИЯ ПОВЕРХНОСТНОЙ ПЛОТНОСТИ ЗАРЯДА ПО ИЗВЕСТНОМУ ЭЛЕКТРИЧЕСКОМУ ПОТЕНЦИАЛУ СРЕДСТВАМИ COMSOL MULTIPHYSICS Текст научной статьи по специальности «Физика»

Павел Игоревич Пичугин Pavel I. Pichugin

магистрант кафедры

техники и электрофизики высоких напряжений, Национальный исследовательский университет «Московский энергетический институт», Москва, Россия

УДК 621.317 DOI: 10.17122/1999-5458-2021-17-3-4-39-48

РЕШЕНИЕ ОБРАТНОЙ ЗАДАЧИ НАХОЖДЕНИЯ ПОВЕРХНОСТНОЙ ПЛОТНОСТИ ЗАРЯДА ПО ИЗВЕСТНОМУ ЭЛЕКТРИЧЕСКОМУ ПОТЕНЦИАЛУ СРЕДСТВАМИ COMSOL MULTIPHYSICS

Актуальность

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

Цель исследования

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

Методы исследования

Задача решается с применением численного моделирования в программном обеспечении на базе метода конечных элементов Comsol Multiphysics. С помощью созданной 3D модели, содержащей воспроизведённую конструкцию ёмкостного зонда и исследуемый заряженный материал, поясняется принцип формирования матрицы весовых коэффициентов, через которую впоследствии устанавливается взаимосвязь между дискретным набором измеренных потенциалов и неизвестным распределением поверхностной плотности заряда.

Результаты

Разработан эффективный метод пересчёта известного распределения электрического потенциала от заряженной поверхности в её поверхностную плотность заряда. Показан принцип построения и задания граничных условий в необходимой для пересчёта 3D модели. На примере одного заданного распределения электрического потенциала была показана высокая эффективность метода получения распределения

- 39

Электротехнические и информационные комплексы и системы. № 3-4, т. 17, 2021

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

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

SOLUTION OF AN INVERSE PROBLEM OF FINDING SURFACE CHARGE DENSITY WITH A KNOWN ELECTRICAL POTENTIAL USING COMSOL

MULTIPHYSICS

Relevance

Polymer materials capable of electric charge accumulation on their surface and in volume have wide application in a great number of applied problems of electrophysics and high-voltage engineering due to their excellent dielectric property. This property can be both an advantage (in electret manufacture), and a disadvantage (when using polymers as an isolating medium). Fundamentally important here is the need to characterize and control the degree of charge using the value of surface charge density which, however, can be determined indirectly through the electric potential only. When measuring the potential using a capacitive probe, difficulties arise in solving an inverse problem of converting the potential into surface charge density, since the relationship between these values largely depends on the geometry of the charged body and the capacitive probe. An effective and universal conversion method will greatly simplify the measurement of the surface charge density that will be helpful in practical tasks of electrical technology and operation of high-voltage equipment with elements made of polymer materials.

Aim of research

Aim of the research is to develop a technique for restoring surface charge density distribution according to a given measured distribution of electric potential created by the surface charge; give recommendations on how to achieve an optimal balance between accuracy and labor intensity of numerical calculations.

Methods

The problem is solved using numerical modeling in finite element method based software Comsol Multiphysics. With the help of the created 3D model containing the reproduced design of the capacitive probe and the charged examined material, the principle of the weight coefficients matrix formation is explained, through which the relationship between a discrete set of measured potentials and the unknown distribution of the surface charge density is subsequently established.

Results

An effective method has been developed to recalculate the known distribution of the electric potential from a charged surface into its surface charge density. The principle of constructing and setting the boundary conditions in the 3D model required for recalculation is shown. On the example of one given distribution of electric potential, the high efficiency of the method for obtaining the distribution of the charge density was shown. It has been established that the accuracy of the method increases with a finer partition of the charged surface into elementary charged areas. Recommendations are given on the use of methods for solving the arising systems of linear algebraic equations. It was also shown how to simplify the technique used in the presence of an axisymmetric distribution of the charge density.

Keywords: inverse problem, electric potential, surface charge density, capacitive probe, floating potential, matrix of weighting coefficients, generalized minimal residual method, finite element method

Постановка проблемы

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

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

Данная задача часто возникает в следующих областях:

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

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

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

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

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

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

Теоретические положения метода

Пусть имеется заряженная поверхность с поверхностной плотностью заряда а(х, y). Поверхностный заряд создаёт в пространстве электрическое поле, которое характеризуется распределением потенциала ф(х, y, z).

Одним из распространённых способов измерения данного потенциала является использование ёмкостных зондов (Capacitive Probe) [5]. Находясь в электрическом поле на некотором расстоянии от заряженной поверхности, согласно закону электростатической индукции, на поверхности зонда наводится заряд противоположного знака Q, а сам зонд приобретает некоторый потенциал (р относительно заранее заданного нулевого потенциала (земли). Зная ёмкость зонда C, можно из соотношения Q = Сф отыскать потенциал зонда. Если заряд равномерно распределён по поверхности (а = const), то в большинстве случаев по измеренному потенциалу зонда ф можно сразу получить искомую поверхностную плотность заряда.

Если же заряд распределён неравномерно по поверхности, то чувствительность данного метода уже недостаточна для точного измерения потенциала. В таких случаях зонд соединяется с вибрирующим или вращающимся вольтметром. Вибрации или вращение элементов вольтметра создают изменение его ёмкости C по заранее известному закону, а нахождение вольтметра вместе с зондом под одним и тем же потенциалом ф создаёт в его цепи протекание ёмкостного тока IC,

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

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

(1)

где й — расстояние между обкладками;

е — диэлектрическая проницаемость среды между обкладками;

е0 — электрическая постоянная, е0 = 8,85-10-12 Ф/м.

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

Поскольку зонд и заряженная поверхность представляют собой аналогичным образом некоторую ёмкость, то для него справедливо соотношение (1) в следующем виде:

(р = Ао, (2)

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

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

Решить данную проблему позволяет дискретизация задачи [6], заключающаяся в разбиении заряженной поверхности на N элементов, каждый из которых заряжен с постоянной поверхностной плотно-

стью заряда ау (/ = 1,..., К). Пусть потенциал измеряется аналогично в N точках пространства. Тогда каждый потенциал (/ = 1, ..., К) по принципу суперпозиции является суммой потенциалов, создаваемых каждым заряженным элементом поверхности с поверхностной плотностью заряда ау:

n

ъ = (3)

7=1

Расписывая выражение (3) для всех N точек измерения, можно получить следующее матричное уравнение:

M

Am

Ai ¿22 M

¿N2

A A О A

a\n ¿in M

A-rn

~ <Р\

<Р2

M M

<JN <Pn

(4)

Из выражения (3) видно, что функция А заменяется матрицей коэффициентов А у, определяющих вклад заряженного элемента в потенциал /-ой точки пространства ф.

Таким образом, обратная задача сводится к нахождению матрицы весовых коэффициентов Ау. Её нахождение возможно с помощью программного обеспечения для решения инженерных и научных задач на базе МКЭ Comsol Multiphysics. Покажем методику отыскания весовых коэффициентов Ау.

Создание расчётной модели в Comsol Multiphysics

Создадим в программе Comsol Multiphysics 3D модель, состоящую из заряженной квадратной плёнки со стороной а = 40 мм, диэлектрической проницаемостью е = 2,3 и толщиной h = 100 мкм, а также произвольно начерченного ёмкостного зонда, служащего для измерения созданного поверхностным зарядом потенциала электрического поля на высоте Н = 10 мм от поверхности плёнки (рисунок 1). При этом исключается 2-координата, измеренный потенциал

Рисунок 1. 3D-модель для нахождения весовых коэффициентов в программе Comsol Multiphysics

Figure 1. 3D model for finding weight coefficients in Comsol Multiphysics program

ф(х, y) зависит теперь только от двух координат.

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

Поскольку МКЭ требует для расчётов задания искусственных границ, система «зонд — плёнка» была помещена в цилиндр больших размеров с нулевыми

Таблица 1. Граничные условия 3D-модели Table 1. Boundary conditions of a 3D model

граничными условиями 2-го рода на нём. Приведём в таблице 1 заданные граничные условия.

Разобьём поверхность плёнки на N квадратных элементов, N = 25 (рисунок 2).

Принцип нумерации заряженных элементов и точек измерения потенциала (точки на рисунке) показан на рисунке 2.

П V-» и

Зарядим в модели первый квадратный элемент поверхностной плотностью заряда а1 = 1 нКл/см2, заряды же всех остальных квадратов примем равным

Граничное условие Место Значение

Ground 1 нижняя сторона плёнки плёнка заземлена (ф = 0 В)

Ground 2 весь объём экрана зонда экран заземлён (ф = 0 В)

Floating Potential весь объём стержня зонда плавающий потенциал на стержне зонда (электростатическая индукция)

Surface Charge Density выбранный квадратный заряженный элемент поверхностная плотность заряда квадрата (а = 1 нКл/см2)

20

15

10

га

? О

S

о о

-10

-15

-20

<Ръ Ü5 + + On

* + +

+ + +

Ol * + +

+<pl Ol +<Р6 Об + + +

-20 -15 -10 -5 0 5 10

Координата х, мм

15 20

Рисунок 2. Разбиение заряженной поверхности на элементы Figure 2. Partition of a charged surface into elements

нулю а] = 0 нКл/см2 (] = 2, ..., Ы). Измерим после полевого расчёта в программе Comsol Multiphysics потенциалы во всех точках г (г = 1, ..., Ы), которые обусловлен^! только зарядом первого элемента. Потенциалы ф(х, у) будем измерять в точках, координаты х, у которых соответствуют центрам квадратных элементов (рисунок 2).

Тогда выполняя матричное умножение в (4), можно получить следующий ряд выражений:

¿и ¿ю А Ал

41 ¿21 м

¿ni

12 ^22 м

¿n2

А А О А

4n ¿21V

м

¿nn

<Р1

0 <Р2

м м

0 <Pn

Л11 ¿21 М

¿n1

12

А

¿22 М

¿n2

А А О А

4n ¿2 n

м

¿nn

0 <Р1

<Р2

м м

0 <Pn

'¿и ¿12 А ¿1n -1 >1

а2 ¿21 ¿22 А ¿2 n <Р2

м м м О м м

oN_ ¿n1 ¿n2 А ¿nn _ _<Pn

(5)

<Р1 =Ап<71> Ч>2 =А21°\> •••> <Рк = ¿та1

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

Пусть далее заряжен второй квадратный элемент а2 = 1 нКл/см2, заряды же всех остальных квадратов примем равным нулю а] = 0 нКл/см2 (] = 1, 3, ..., Ы). Проделав ту же самую процедуру измерения потенциалов точек можно получить значения коэффициентов второго столбца:

Л

■ (6)

&1 - ^12°2> Ф2 ~ А22&2> Ры ~ ¿N2^2

Заряжая далее последовательно все остальные квадратные элементы, можно получить всю матрицу весовых коэффициентов А.

Зная матрицу А и имеющееся реальное распределение потенциала ф(х, у) в данных точках, можно решить систему уравнений (4) и получить интересующее рас-

пределение a(x, y) поверхностной плотности заряда в образованных квадратных элементах:

(7)

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

Особое внимание следует уделить решению матричного уравнения (7). Значительная доля обратных задач являются вследствие дискретизации изначально непрерывного решения математически некорректно поставленными. Часто это обстоятельство приводит к неустойчивости решения системы линейных алгебраических уравнений (СЛАУ). Опыт решения рассматриваемой обратной задачи показывает, что из-за очень малого определителя матрицы весовых коэффициентов А получить адекватное решение традиционными методами решения СЛАУ не представляется возможным. В этой связи при решении обратной задачи следует прибегнуть к использованию итеративных решателей, среди которых эффективными оказываются обобщённый метод минимальных невязок (Generalized Minimal Residual Method) или его подвид метод квази-минималь-ных невязок (Quasi-Minimal Residual Method). Используя для решения данной обратной задачи программу MATLAB,

можно обратиться к соответствующим решателям Gmres или Qmr [7]. Сходимость к решению наблюдается при задании относительной погрешности вычисления 10-2:

погт {<р-Ао)

<10

-2

(8)

попп(^)

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

Пример использования метода

Продемонстрируем применение изложенной методики на конкретном примере.

Пусть измеренное распределение потенциала ф(х, у) на высоте Н = 10 мм от заряженной плёнки является симметричным относительно центра системы координат (0; 0) и описывается аналитическим выражением (х, у в мм):

= (9)

Проведём два расчёта, различающиеся числом квадратных элементов, на которые разбита заряженная поверхность: N = 25 (рисунок 2) и N = 36. Подставляя координаты центров квадратов в выражение (9), можно получить потенциалы точек ф1 (/ = 1, ..., К).

Подготовив по вышеизложенной методике матрицы весовых коэффициентов А и используя какой-либо из перечисленных итеративных решателей, можно получить распределения поверхностной плотности заряда по квадратам а] (/ = 1, ..., Щ. Приведём полученные поверхности (рисунок 3), дополнительно построим график распределения поверхностной плотности заряда вдоль одной из диагоналей плёнки (с 5 до 21 элемента

для N = 25, с 6 до 31 элемента для N = 36). Поскольку распределение потенциала (9) задаёт параболоид, то распределения плотностей заряда вдоль диагоналей квадратной плёнок должно быть также симметричным относительно центра плёнки. Для возможности сравнения двух результатов плотности заряда при более мелком разбиении интерполированы на узлах более крупной сетки элементов.

Из рисунка 3, Ь видно, что при более мелком дроблении сетки полученный результат численно мало отличается от результата при дроблении сетки на N = 25 элементов (рисунок 3, а). Относительная погрешность максимальной поверхностной плотности заряда ам=25 по сравнению с ам=36 составляет около 2,5 %. Однако из рисунка 3, Ь можно заметить, что более мелкое дробление сетки привело к повышению симметричности распределения, что и должно наблюдаться в данной задаче. Таким образом, начиная с числа элементов N = 36, можно получить удовлетворительное по точности решение.

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

Расчёты показывают, что довольно хорошая точность может быть достигнута при числе колец N = 20.

Номер элемента I

a) N = 25; b) N = 36; c) вдоль диагонали заряженной поверхности a) N = 25; b) N = 36; c) along the diagonal of the charged surface Рисунок 3. Результаты вычислений поверхностной плотности заряда Figure 3. Results of calculating surface charge density

-20 -15 -10 -5

10 15 20

X, MM

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

Figure 4. Partition of a charged surface into ring elements in case of axisymmetric

distribution of the charge density

Выводы

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

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

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

1. Кондратюк А.Д. Электретные свойства диэлектриков // Ярославский педагогический вестник. 1998. № 3. С. 109-112.

2. Бойцов В.Г., Рычков Д.А. Полимерные электреты в инновационных технологиях // Известия РГПУ. Естественные и точные науки: Научный журнал. 2002. № 2(4). С. 118-132.

3. Smili K., Herous L. Corona Charging and Charge Decay on Polyethylene Terephthalate Films (PET) // Electrotehnica, Electronica, Automatica (EEA). 2019. Vol. 67. No. 3. P. 81-90. ISSN 1582-5175.

4. Кушелев В.П., Орлов Г.Г., Сорокин Ю.Г. Охрана труда в нефтеперерабатывающей и нефтехимической промышленности: учебник для вузов. М.: Химия, 1983. 472 с.

5. Noras M.A. Non-Contact Surface Charge/Voltage Measurements // Capacitive Probe — Principle of Operation. 2002. No. 3001. P. 1-8.

6. Takuma T., Yashima M., Kawamoto T. Principle of Surface Charge Measurement for Thick Insulating Specimens // IEEE Transactions on Dielectrics and Electrical Insulation. 1998. Vol. 5. No. 4. P. 497-504. doi: 10.1109/94.708264.

7. MATLAB: Iterative Methods for Linear Systems. URL: https://www. mathworks.com/help/matlab/math/iterative-methods-for-linear-systems.html (дата обращения 01.12.2021).

References

1. Kondratyuk A.D. Elektretnye svoi-stva dielektrikov [Electret Properties of

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

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

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

Dielectrics]. Yaroslavskii pedagogicheskii vestnik — Yaroslavl Pedagogical Bulletin, 1998, No. 3, pp. 109-112. [in Russian].

2. Boitsov V.G., Rychkov D A. Poli-mernye elektrety v innovatsionnykh tekhno-logiyakh [Polymer Electrets in Innovative Technologies/ Izvestiya RGPU. Estest-vennye i tochnye nauki: Nauchnyi zhurnal — Izvestiya RSPU. Natural and Exact Sciences: Scientific Journal, 2002, No. 2(4), pp. 118-132. [in Russian].

3. Smili K., Herous L. Corona Charging and Charge Decay on Polyethylene Tere-phthalate Films (PET). Electrotehnica, Electronica, Automatica (EEA), 2019, Vol. 67, No. 3, pp. 81-90. ISSN 1582-5175.

4. Kushelev V.P., Orlov G.G., Soro-kin Yu.G. Okhrana truda v neftepereraba-tyvayushchei i neftekhimicheskoi promysh-lennosti: uchebnik dlya vuzov [Labor Protection in Oil Refining and Petrochemical Industry: Textbook for Universities]. Moscow, Khimiya Publ., 1983. 472 p. [in Russian].

5. Noras M.A. Non-Contact Surface Charge/Voltage Measurements. Capacitive Probe — Principle of Operation, 2002, No. 3001, pp. 1-8.

6. Takuma T., Yashima M., Kawamoto T. Principle of Surface Charge Measurement for Thick Insulating Specimens. IEEE Transactions on Dielectrics and Electrical Insulation, 1998, Vol. 5, No. 4, pp. 497-504. doi: 10.1109/94.708264.

7. MATLAB: Iterative Methods for Linear Systems. URL: https://www.math-works.com/help/matlab/math/iterative-methods-for-linear-systems.html (accessed 01.12.2021).