Метод сопряжённого уравнения в задаче об определении источника диффузии Текст научной статьи по специальности «Математика»

УДК 519.6

DOI: 10.14529/mmph190303

МЕТОД СОПРЯЖЁННОГО УРАВНЕНИЯ В ЗАДАЧЕ ОБ ОПРЕДЕЛЕНИИ ИСТОЧНИКА ДИФФУЗИИ

1 2 В.А. ЛитвиновВ.В. Учайкин2

1 Барнаульский юридический институт, г. Барнаул, Российская Федерация

2 Ульяновский государственный университет, г. Ульяновск, Российская Федерация E-mail: lva201011@yandex.ru

Объектом исследования работы являются дифференциальные уравнения диффузии (теплопроводности). Предметом исследования является алгоритм определения функции источника или начальных условий задачи по экспериментально измеряемым величинам. В основу исследования положено двойственное представление функционалов, соответствующих экспериментально наблюдаемым величинам в процессах массо- и теплообмена. Обратная задача сформулирована в виде интегральных уравнений первого рода, ядром которых является сопряженная функция (функция ценности), получаемая как решение сопряженного в смысле Лагранжа уравнения диффузии (теплопроводности) с функцией чувствительности детектора в правой части. При этом решение сопряженных уравнений путем замены переменных сводится к решению прямых уравнений. Для регуляризации решения уравнения Вольтерры первого рода, соответствующего задаче восстановления зависимости граничного условия от времени, предложено использовать минимизацию невязки для переопределенной системы линейных уравнений. Задача восстановления зависимости начального условия от координаты сформулирована в виде уравнения Фредгольма I рода, для решения которого применен метод регуляризации Тихонова. Приведены результаты модельных расчётов по восстановлению временной зависимости источников, заданных гладкой функцией, ступенчатой функцией и функцией с гармонической составляющей в задаче об одномерной диффузии в однородной среде. Из этих результатов видно, что при выбранных параметрах расчетов полученные предлагаемым методом решения ведут себя регулярно и обладают вполне приемлемой точностью даже несмотря на то, что значения искомой функции на заданном интервале поиска изменяются на шесть порядков. В этом авторы видят главное отличие предложенного ими метода от других подходов к решению данной задачи

Ключевые слова: обратная задача; диффузия; теплопроводность; ценность; чувствительность.

Введение

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

Как и в классе прямых задач, методы решения задач обратного типа разделяются на две группы: аналитические и численные. К первой относятся построение функционально-инвариантных решений гиперболических уравнений, аналитические представления решений и коэффициентов параболических уравнений, представление решения и коэффициента уравнения Штурма-Лиувилля с применением в обратных задачах теории рассеяния, построение гармонических и других потенциалов для вычисления решений (скорости) и коэффициентов (давления) системы уравнений газовой динамики и др. [1-3].

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

Литвинов В.А., Метод сопряжённого уравнения в задаче

Учайкин В.В. об определении источника диффузии

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

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

Подробное описание постановки различного типа обратных задач можно найти в работах [6, 7]. В настоящей работе рассматривается частный случай обратной задачи по определению функции источника, основанный на двойственном представлении показаний детектора с использованием функции ценности [15, 16] и в этом смысле являющийся продолжением идей, сформулированных в наших ранних работах [17, 18].

Описание метода

Уравнение теплопроводности (диффузии), как и другие уравнения переноса, могут быть представлены в операторной форме

ЬФ(и) = 5 (и), (1)

где и - набор фазовых координат, например, для одномерного уравнения диффузии и = х}. Наличие известного решения Ф(и) превращает определение источника 5(и) в прямую задачу о вычислении ЬФ(и), и особых проблем не возникает. Однако в реальном физическом эксперименте измеряются не значения функции Ф(и), а некоторый функционал от неё, представляющий показание прибора (детектора). Ограничившись линейным представлением, запишем его в виде скалярного произведения

J = | ёиФ(и)Б(и) = (Ф, П). (2)

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

Ь+Ф+ (и,а) = П(и,а), (3)

оператор Ь в котором удовлетворяет соотношению Лагранжа:

(Ф+ , Ь Ф) = (Ф , Ь+ Ф+).

Согласно принципу двойственности показания такого детектора могут быть выражены также через функцию ценности:

J (а) = | йиФ+ (и,а)5 (и) г (Ф+, 5). (4)

Таким образом, использование двойственного представления функционала естественным образом сводит задачу параметрического определения функции источника (в чем и состоит специфика данной работы) к решению интегрального уравнения первого рода. Будет это уравнение Вольтерры или Фредгольма, зависит от свойств оператора переноса Ь и граничных условий задачи. Заметим также, что отмена правой части уравнения (1) и использование вместо неё начального условия

Ф (0, х) = g (х) (5)

эквивалентно задаче (1) с функцией источника:

5 (г, х) = 8(Г) g (х),

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

Если функцию чувствительности детектора выбрать в виде

Б(Г,х) = 8(Г - Го)8(х - хо), (6)

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

Рассмотрим модельную задачу для одномерной диффузии в неограниченном пространстве с различными точечными источниками:

5 (Г, х) = ф(Г )8( х). (7)

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

В случае постоянного коэффициента диффузии а2 уравнение (3) имеет вид:

дФ+ 2 д2Ф

—---a -uyi --xo

dt дх2

- a-S(t - to)S( x - xo). (8)

Данное уравнение отличается от классического «прямого» уравнения только знаком перед производной по времени. Путем замены переменных оно легко сводится к классическому уравнению диффузии (теплопроводности), для которого решение известно:

Ф+ (Г, х, Го, Хо) = ■ 1 ехр Г - (Х2- Хо)\ 1. (9)

2^п(Го - Г) ^ 4а2(Го - Г))

Пусть нам известно, что источник является точечным, а найти требуется зависимость мощности источника от времени, характеризуемую множителем ). Подставляя (7) и (9) в (4), получим интегральное уравнение Вольтерры I рода:

Го 1 Г х 2 1

Í

2ay]n(t0 -1)

exp

xo

<f)(t)dt - J(Xo,to). (10)

4a2(to -1

Заметим, что для краевой задачи теплопроводности в стержне единичной длины, рассматриваемой в работах [8, 9], функция ценности будет представляться в виде гармонического ряда с весовыми сомножителями ~ exp(-a2Пt), что соответствует виду ядра уравнений, решаемых в цитируемых работах.

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

Регуляризация решения

Вопросы регуляризации решений интегральных уравнений I рода рассматривались во многих работах. Одним из способов получить устойчивое решение для уравнения Вольтерры I рода является сведение его к уравнению II рода путем дифференцирования. В рассматриваемом случае ядро интегрального оператора таково, что после первого дифференцирования уравнение по-прежнему будет первого рода, так как Ф(^, to) = o. Дальнейшее дифференцирование приводит к такому же результату.

Другим способом решения интегральных уравнений Вольтерры I рода является применение различных методов регуляризации [8, 19, 2o], основная идея которых состоит в минимизации некоторого функционала от искомого решения.

Предположим, что нам известны значения функционала в разные моменты времени при фиксированном xo: Jk = J (xo, tk), к = 1, M.

Зададим две сетки значений моментов времени t1 < t2 < ... < tN = tM, в которые необходимо определить значения искомой функции <p(t), входящей в выражение (Ю). Представим её интерполяционным полиномом Лагранжа степени n:

т - ± t(tm+l) п T^j. (11)

l-o j-o, j-l (tm+l tm + j )

В выражении (11) значение индекса m определяет номер начального узла интерполяции и влияет на вид матрицы квадратурных коэффициентов, получающихся при подстановке выражения (11) в (1o). Подставим (11) в (Ю) и выполним интегрирование для всех значенийyk, получим

N

I АктФ((т) = ^, 1 < к <М.

(12)

т=1

При N < М система (12) переопределена и для её решения можно воспользоваться методом наименьших квадратов для минимизации невязки. Такой метод, в частности, реализован в пакете МЛТЬЛБ для решения переопределенных систем.

В случае мгновенного протяженного источника Б(Г, х) = 5(Г^(х) подстановка данного выражения и функции ценности (7) в равенство (4) приведет к уравнению:

1

2а» пГ,

ехр

(х - х0) 4а2Г0

2 Л

g (х)ёх = J (Хо, ^о ).

(13)

В данном случае момент наблюдения Г0 является параметром. Для восстановления функции g(x), определяющей пространственную зависимость мгновенного источника, потребуется набор значений функционала J в нескольких точках. Будем в дальнейшем для данной задачи использовать обозначение J(у) = J(у, Г0).

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

у/ = | ёу

| К (у, х) g (х)ёх - J (у)

+

а1 g 2 (х)ёх + а1 [ g'(х)]2 ёх,

где а- малый параметр регуляризации.

При получении результатов, приведенных ниже, значение параметра а определялось простым перебором, что является приемлемым при небольшом числе расчетов и применяется рядом авторов. В общем случае пределы интегрирования в приведенном выражении должны быть бесконечными, но при построении квадратурных формул в рассматриваемом примере использовался отрезок [-7, 7], выбранный исходя из поведения функций J(y) и g(x).

Обсуждение результатов

На рис. 1 и рис. 2 приведены результаты решения уравнения (10) для различных характеров поведения искомой функции ф(Г). На рис. 1 представлены расчеты для монотонно изменяющихся функций ф(Г) = Г3/2 и ф(Г) = ехр(-3Г/2). Значения второй функции на заданном интервале поиска изменяется на шесть порядков. Из рис. 1 и 2 видно, что при выбранных параметрах расчетов решения ведут себя регулярно и точность вполне достаточная. Расчеты проводились для шага по времени А/1 = 0,25 иМ= Ш.

ю2

0 1 2 3 4 5 6 7 8 9^10

Рис. 1. Результаты решения уравнения (10) для монотонных источников. Линии - точное исходное значение, +, о - численное решение; + - ф(Г) = Г3/2 , о - ф(Г) = ехр(-3Г/2)

I

На рис. 2 приводятся результаты решения уравнения (10) для источников с гармонической составляющей: ф(?) = 1 + 8т(5?) и ф(?) = 1 + 8т(2?) /(2?). В обоих случаях решалась переопределенная система линейных уравнений с количеством уравнений в два раза большим числа неизвестных. И в этом случае решения ведут себя регулярно и для выбранной сетки с шагом А? = 0,2 достигается хорошая точность.

2.5

2

1.5

1

0.5

0

-0.5

0 2 4 6 8 Г 10

Рис. 2. Результаты решения уравнения (10) для источников с гармонической составляющей. Линии - точное исходное значение, +, о - численное решение; + - ф(?) = 1 + 8т(2г)/(2г), о - ф(Г) = 1 + 8т(5г)

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

Заключение

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

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

Благодарность

Авторы благодарны Российскому фонду фундаментальных исследований за финансовую поддержку работы (грант 18-51-53018).

Рис. 3. Результаты восстановления пространственного мгновенного источника. Сплошная линия и пунктир - точные решения. Маленькие точки и кружки - результат решения уравнения (13). Сплошная линия и кружки - g(x) ~

exp(-x2/2). Пунктирная линия - прямоугольный импульс

Литература

1. Соболев, С. Л. Функционально-инвариантные решения волнового уравнения / С. Л. Соболев // Тр. Физ .-мат. ин-та им. В.А. Стеклова, 1934. - Т. 5. - С. 259-264.

2. Колмогоров, А.Н. Об аналитических методах в теории вероятностей / А.Н. Колмогоров // Успехи мат. наук. - 1938. - Вып. 5.- С. 5-41.

3. Лаврентьев, М.М. Одномерные обратные задачи математической физики / М.М. Лаврентьев, К.Г. Резницкая, В.Г. Яхно. - Новосибирск: Наука, 1982.- 88 с.

4. Самарский, А.А. Численные методы решения обратных задач математической физики / А.А. Самарский, П.Н. Вабищевич. - М.: Изд-во ЛКИ, 2009. - 478 с.

5. Кабанихин, С.И. Обратные и некорректные задачи / С.И. Кабанихин. - Новосибирск: Сибирское научное изд-во, 2009. - 457 с.

6. Гончарский, А.В. Численные методы решения обратных задач астрофизики / А.В. Гончарский, А.М. Черепащук, А.Г. Ягола. - М.: Наука, 1978. - 335 с.

7. Калинина, Е.А. Численное исследование обратной задачи восстановления плотности источника двумерного нестационарного уравнения конвекции-диффузии / Е.А. Калинина // Дальневосточный математический журнал. - 2004. - Т. 5, № 1. - С. 89-99.

8. Япарова, Н.М. Численный метод решения некоторых обратных задач теплопроводности с неизвестными начальными условиями / Н.М. Япарова // Вестник ЮУрГУ. Серия: Компьютерные технологии, управление, радиоэлектроника. - 2015. - Т. 15, № 2. - С. 55-65.

9. Япарова, Н.М. Метод решения обратной задачи идентификации функции источника с использованием преобразований Лапласа / Н.М. Япарова // Вестник ЮУрГУ. Серия: «Вычислительная математика и информатика». - 2016. - Т. 5, № 3. - С. 20-35.

10. Кожанов, А.И. Обратные задачи восстановления правой части специального вида в параболическом уравнении / А.И. Кожанов // Математические заметки СВФУ. - 2016. - Т. 23, № 4. -

C.31-44.

11. Hasanov, A. An analysis of inverse source problems with final-time measured output data for the heat conduction equation: A semigroup approach / A. Hasanov, M. Slodicka // Applied Mathematics Letters, 2013. - Vol. 26, Iss. 2. - P. 207-214.

12. An inverse time-dependent source problem for the heat equation / A. Hazanee, M.I. Ismailov,

D. Lesnic, N.B. Kerimov // Applied Numerical Mathematics. - 2013. - Vol. 69. - P. 13-33.

13. Inverse problem of time-dependent heat sources numerical reconstruction / L. Yang, M. Dehghan, J.-N. Yu, G.-W. Luo. Mathematics and Computers in Simulation. - 2011. - Vol. 81, Iss. 8. -P.1656-1672.

14. Марчук, Г.И. Методы вычислительной математики / Г.И. Марчук. - M.: Мир, 1980. -430 с.

15. Marchuk, G.I. Adjoint equations and analysis of complex systems / G.I. Marchuk // Mathematics and Its Applications (MAIA, volume 295). - Springer, Dordrecht, 1995. - 468 p.

16. Литвинов, В.А. Вариации ценности в проблеме изучения широких атмосферных ливней / В.А. Литвинов, В.В. Учайкин // Известия вузов. Физика. - 1986. - Т. 29, № 2. - С. 128.

17. Литвинов, В.А. Метод функциональных производных в проблеме чувствительности ШАЛ / В.А. Литвинов, В.В. Учайкин // Известия вузов. Физика. - 1986. - Т. 29, № 12. - С. 96.

18. Лаврентьев, М.М. Об интегральных уравнениях первого рода / М.М. Лаврентьев // Доклады АН СССР. - 1959. - Т. 127, № 1. - С. 31-33.

19. Магницкий, Н.А. О приближенном решении некоторых интегральных уравнений Воль-терра первого рода / Н.А. Магницкий // Вестник Московского университета. Серия 15: «Вычислительная математика и кибернетика». - 1978. - № 1. - С. 91-98.

Поступила в редакцию 21 апреля 2019 г

Bulletin of the South Ural State University Series "Mathematics. Mechanics. Physics" _2019, vol. 11, no. 3, pp. 20-27

DOI: 10.14529/mmph190303

ADJOINT EQUATION METHOD FOR SOLVING THE INVERSE DIFFUSION SOURCE PROBLEM

V.A. Litvinov1, V.V. Uchaikin2

Barnaul law Institute, Barnaul, Russian Federation 2 Ulyanovsk state University, Ulyanovsk, Russian Federation E-mail: lva201011@yandex.ru

The objects of the research are differential equations of diffusion (or thermal conductivity) kind. The subject of research is the algorithm for determining the function of the source or the initial conditions of the problem as per the experimentally measured values. The approach is based on a dual representation of functionals corresponding to experimentally observed quantities in the processes of mass and heat transfer. The inverse problem is formulated in the form of integral equations of the first kind, the core of which is the adjoint function (importance function) obtained as a solution of the adjoint (in the Lagrange sense) diffusion (thermal conductivity) equation with the detector sensitivity function in the right-hand side. Meanwhile, solving the adjoint equations by changing the variables is reduced to solving direct equations. To regularize the solution of the Volterra equation of the first kind corresponding to the problem of recovering the dependence of the boundary condition on time, the residual minimization for an overdetermined system of linear equations has been proposed to be used. The problem of reconstructing the dependence of the initial condition on the coordinate is formulated as a Fredholm equation of the first kind, the solution to which has been obtained using the Tikhonov regularization method. The results of model calculations are presented for the restoration of the time dependence of the sources set by a smooth function, a step function and a function with a harmonic component in the problem of one-dimensional diffusion in a homogeneous medium. These results prove that with the selected calculation parameters, the solutions obtained by the proposed method behave regularly and have quite an acceptable accuracy, even though the values of the sought-for function within the set search interval change by six orders of magnitude. This is seen by the authors as the main difference of their method from other approaches to solving this problem.

Keywords: inverse problem; diffusion; thermal conductivity; importance-function; sensitivity.

References

1. Sobolev S. Functionally-invariant solutions of wave equation, Travaux Inst. Physico-Math. Stekloff 1934, vol. 5, pp. 259-264. (in Russ.).

2. Kolmogorov A.N. On the analytic methods of probability theory, Uspekhi Mat. Nauk, 1938, no. 5, pp. 5-41. (in Russ.).

3. Lavrent'ev M.M., Reznitskaya K.G., Yakhno V.G. Odnomernye obratnye zadachi matematicheskoy fiziki (One-dimensional inverse problems of mathematical physics). Novosibirsk: Nauka Publ., 1982, 88 p. (in Russ.).

4. Samarskiy A.A., Vabishchevich P.N. Chislennye metody resheniya obratnykh zadach matematicheskoy fiziki (Numerical methods for solving inverse problems of mathematical physics). Moscow, Izd-vo LKI publ., 2009, 478 p. (in Russ.).

5. Kabanikhin S.I. Obratnye i nekorrektnye zadachi (Inverse and ill-posed problems). Novosibirsk, Sibirskoe nauchnoe izd-vo Publ., 2009, 457 p. (in Russ.).

6. Goncharskiy A.V., Cherepashchuk A.M., Yagola A.G. Chislennye metody resheniya obratnykh zadach astrofiziki (Numerical methods for solving inverse problems of astrophysics). Moscow, Nauka Publ., 1978, 335 p. (in Russ.).

7. Kalinina A.E. Numerical study of inverse problem of identification of source density for the two-dimensional non-stationary convection-diffusion equation. Dal'nevost. Mat. Zh., 2004, Vol. 5, no. 1, pp. 89-99. (in Russ.).

8. Yaparova N.M. Numerical Method for Solving some Inverse Heat Conduction Problems with Unknown Initial Conditions. Bulletin of the South Ural State University. Ser. Computer Technologies, Automatic Control, Radio Electronics, 2015, Vol. 15, no. 2, pp. 55-65. (in Russ.). DOI: 10.14529/ctcr150206

9. Yaparova N.M. Method for Solving an Inverse Term Source Problem Based on the Laplace Transform. Bulletin of the South Ural State University. Series: Computational Mathematics and Software Engineering, 2016, Vol. 5, no. 3, pp. 20-35. (in Russ.). DOI: 10.14529/cmse160302

10. Kozhanov A.I. Inverse problems of recovering the right-hand side of a special type of parabolic equations. Mathematical notes of NEFU, 2016, Vol. 23, no. 4, pp. 31-45. (in Russ.).

11. Hasanov A., Slodicka M. An analysis of inverse source problems with final-time measured output data for the heat conduction equation: A semigroup approach. Applied Mathematics Letters, 2013, Vol. 26, Iss. 2, pp. 207-214. DOI: 10.1016/j.aml.2012.08.013

12. Hazanee A., Ismailov M.I., Lesnic D., Kerimov N.B. An inverse time-dependent source problem for the heat equation. Applied Numerical Mathematics, 2013, Vol. 69, pp. 13-33. DOI: 10.1016/j.apnum.2013.02.004

13. Yang L., Dehghan M., Yu J.-N., Luo G.-W. Inverse problem of time-dependent heat sources numerical reconstruction. Mathematics and Computers in Simulation, 2011, Vol. 81, Iss. 8, pp. 16561672. DOI: 10.1016/j.matcom.2011.01.001

14. Marchuk G.I. Metody vychislitel'noy matematiki (Computational Mathematics Methods), Moscow, Mir Publ., 1980, 430 p. (in Russ.).

15. Marchuk G.I. Adjoint equations and analysis of complex systems. Mathematics and Its Applications (MAIA, vol. 295), Springer, Dordrecht, 1995, 468 p. DOI: 10.1007/978-94-017-0621-6

16. Litvinov B.A., Uchaykin B.B. Izvestiya vuzov. Fizika, 1986, Vol. 29, no. 2, p. 128. (in Russ.).

17. Litvinov B.A., Uchaykin B.B. Izvestiya vuzov. Fizika, 1986, Vol. 29, no. 12, p. 96. (in Russ.).

18. Lavrent'ev M.M. Doklady ANSSSR, 1959, Vol. 127, no. 1, pp. 31-33. (in Russ.).

19. Magnitskiy N.A. Vestnik Moskovskogo universiteta. Seriya 15: "Vychislitel'naya matematika i kibernetika", 1978, no. 1, pp. 91-98. (in Russ.).

Received April 21, 2019