CC BY
4
Вестник КРАУНЦ. Физ.-мат. науки. 2024. Т. 47. №2. C.35-57. ISSN 2079-6641
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ " https://doi.org/10.26117/2079-6641-2024-47-2-35-57 Научная статья
Полный текст на русском языке УДК 519.642.2, 517.972.7
О задаче оптимизации для определения вида функциональной зависимости переменного порядка дробной производной типа Герасимова-Капуто
Д. А. Твёрдый*, Р. И. Паровик
Институт космофизических исследований и распространения радиоволн ДВО РАН, 684034, е. Паратунка, ул.Мирная, д. 7, Россия
Аннотация. При решении задач математического моделирования часто приходится обращаться к теории интегрально-дифференциального исчисления. С ее помощью можно описывать динамические процессы самой разной природы. Использование аппарата дробных производных позволяет уточнить некоторые из этих моделей за счет учета в уравнениях эффекта памяти. Данный эффект выражается в зависимости текущего состояния динамической системы от предыдущих состояний, то есть нелокальности. Интенсивность этого эффекта будет определяться значением показателя степени дробной производной. Классически это некое значение а является нецелым и постоянным. Однако существуют обобщения дробных производных на случай переменной во времени нелокальности а^) и других функциональных зависимостей. Подобные дробно-дифференциальные модели все чаще находят свое применение в теории и практике физико-математических, а также технических наук. Однако, учитывая понимание природы моделируемого процесса, подбор различных параметров таких моделей приходится осуществлять эмпирически. Например, модельные параметры уточняются путем перебора значений и сопоставления временных рядов: результатов моделирования и экспериментальных данных, представляющих процесс. Это продолжается до тех пор, пока результаты моделирования не начнут качественно аппроксимировать данные. Такой подход трудоемок, что неизбежно приводит нас к идеям о решении обратных задач. Цель данной работы - показать, что с помощью методов безусловной оптимизации возможно решение обратных задач для определения вида функциональной зависимости а(1;). Прямая задача определяется как задача Коши для дробного уравнения, где производная понимается в смысле Герасимова-Капуто с переменным показателем степени дробной производной аЩ. Прямая задача решается численно с помощью нелокальной неявной конечно-разностной схемы. Обратная задача определяется как задача дискретной минимизации функции а^) на основе экспериментальных данных. В качестве метода для решения был выбран итерационный метод Левенберга-Марквардта. На тестовых примерах было показано, что метод Левенберга-Марквардта действительно может быть использован для безусловной оптимизации с целью определения вида функции и её оптимальных значений в конкретных моделях.
Ключевые слова: Обратные задачи, безусловная оптимизация, Ньютоновские методы минимизации функции, алгоритм Левенберга-Марквардта, дробные производные, Герасимов-Капуто, эффект памяти, нелокальность, неявные конечно-разностные схемы.
Получение: 29.04.2024; Исправление: 30.05.2024; Принятие: 09.06.2024; Публикация онлайн: 25.08.2024
Для цитирования. Твёрдый Д. А., Паровик Р. И. О задаче оптимизации для определения вида функциональной зависимости переменного порядка дробной производной типа Герасимова-Капуто // Вестник КРАУНЦ. Физ.-мат. науки. 2024. Т. 47. № 2. C. 35-57. EDN: PVTXPV. https://doi.org/10.26117/2079-6641-2024-47-2-35-57. Финансирование. Исследование выполнено за счет гранта Российского научного фонда № 23-71-01050, https://rscf.ru/project/23-71-01050/
Конкурирующие интересы. Конфликтов интересов в отношении авторства и публикации нет. Авторский вклад и ответственность. Авторы участвовали в написании статьи и полностью несут ответственность за предоставление окончательной версии статьи в печать.
* Корреспонденция: А E-mail: [email protected]
Контент публикуется на условиях Creative Commons Attribution 4.0 International License © Твёрдый Д. А., Паровик Р. И., 2024
© ИКИР ДВО РАН, 2024 (оригинал-макет, дизайн, составление)
Vestnik KRAUNC. Fiz.-Mat. nauki. 2024. vol. 47. no. 2. P. 35-57. ISSN 2079-6641
MATHEMATICAL MODELING
" https://doi.org/10.26117/2079-6641-2024-47-2-35-57 Research Article Full text in Russian MSC 26A33, 49N45
The optimization problem for determining the functional dependence of the variable order of the fractional derivative of
the Gerasimov-Caputo type
D. A. Tverdyi*, R.I. Parovik
Institute of Cosmophysical Research and Radio Wave Propagation FEB RAS, 684034, Paratunka village, Mirnaya str., 7, Russia
Abstract. When solving mathematical modeling problems, it is often necessary to turn to the theory of integral and differential calculus. This theory can be used to describe dynamic processes of various types. The use of fractional derivatives allows us to refine some models by taking into account the memory effect, which is expressed in the equations depending on the current state of the system from previous states. This effect is called non-locality and its intensity is determined by the value of the exponent in the fractional derivative. Classically, this value a a noninteger constant, but there are also generalizations for time-varying nonlocality and other functional dependencies. Fractional differential models are finding increasing application in the physical, mathematical, and technical sciences. However, given the nature of the modeled process, the selection of various parameters for such models must be carried out empirically. Model parameters are refined by iterating through values and comparing simulation results with experimental data representing the process. This process continues until the results begin to qualitatively approximate the data, which is a time-consuming process that inevitably leads to ideas about solving inverse problems. The purpose of this work is to demonstrate that it is possible to use methods of unconditional optimization to solve inverse problems and determine the type of functional dependence a(t). The direct problem is formulated as a Cauchy problem for a fractional differential equation, where the derivative is interpreted in the sense of Gerasimov-Caputo with a variable exponent a(t) for the fractional derivative. The direct problem is solved numerically using a nonlocal, implicit finite difference scheme. The inverse problem is defined as the problem of discrete minimization of the function a(t) based on experimental data. To solve this problem, we have chosen the Levenberg-Marquardt iterative method. Through test examples, we have shown that this method can be used for unconstrained optimization to determine the shape of the function a(t) and its optimal values in various models.
Key words: Inverse problems, non-conditional optimization, Newton methods of function minimization, Levenberg-Marquardt algorithm, fractional derivatives, Gerasimov-Caputo, memory effect, non-locality, implicit finite-difference schemes.
Received: 29.04.2024; Revised: 30.05.2024; Accepted: 09.06.2024; First online: 25.08.2024
For citation. Tverdyi D.A., Parovik R.I. The optimization problem for determining the functional dependence of the variable order of the fractional derivative of the Gerasimov-Caputo type. Vestnik KRAUNC. Fiz.-mat. nauki. 2024, 47: 2, 35-57. EDN: PVTXPV. https://doi.org/10.26117/2079-6641-2024-47-2-35-57.
Funding. The research was funded by a grant from the Russian Science Foundation, project number 23-71-01050, which can be found at https://rscf.ru/project/23-71-01050/
Competing interests. There are no conflicts of interest regarding authorship and publication.
Contribution and Responsibility. All authors contributed to this article. Authors are solely responsible for providing the final version of the article in print. The final version of the manuscript was approved by all authors.
* Correspondence: A E-mail: [email protected] ^
The content is published under the terms of the Creative Commons Attribution 4-0 International License © Tverdyi D. A., Parovik R. I., 2024
© Institute of Cosmophysical Research and Radio Wave Propagation, 2024 (original layout, design, compilation)
Введение
Интегрально-дифференциальное исчисление активно развивается уже на протяжении более 300 лет, и за это время зарекомендовало себя надёжным и мощным инструментом решения теоретических и прикладных задач. На протяжении всей истории этого направления математической науки исследователей увлекал вопрос о возможности обобщения понятий интеграла и производной от целочисленного к вещественному порядку [1]. Со времён Г. В. Лейбница и Г. Ф. Лопиталя (1965) были разработаны самые разные определения понятия оператора дробного дифференцирования: Римана-Лиувилля, Вейля, Грюнвальда-Летникова и многие другие. Также всё чаще находят приложения дробных производных и интегралов в тех областях, где целочисленного определения было недостаточно. В квантовой физике [2], вязкоупругости [3, 4] и механике твердого тела [5], аномальных диффузионных процессов [6, 7], фрактальных осцилляторах [8] и т. д.
В нашем исследовании дробные производные (FD) интересны тем, что они полезны для иллюстрации памяти и наследственных свойств многих динамических процессов [9]. Этот эффект также называют эредитарностью, и одним из первых ученых, которые занялись выведением этой теории, стал итальянский математик и физик В. Вольтерра [10]. Эредитарность и ее связь с дробными производными можно показать следующим образом. Рассмотрим следующее интегро-дифференциальное уравнение (IDE):
K(t - a)ü(a)da = F(u(t),t),
(1)
где, согласно принципам эредитарности В. Вольтера [10], память учитывает всю возможную предысторию. Так как на практике мы не обладаем данными о всей истории процесса, то следует рассматривать интеграл на подинтервале (0,t) где, t G [0, T] - текущее время моделирования, T > 0 - общее время моделирования. Главным моментом является то, что выбор вида ядра K(t — а) как степенной функции K(t — а) (р—1—— позволяет нам переходить к операторам дробного исчисления [11,12].
Переходя к формализму FD при работе с такими IDE, в нашем исследовании мы будем рассматривать интеграл в левой части (1) в терминах оператора дробного дифференцирования Герасимова-Капуто [13,14] постоянного 0 < а < 1 порядка:
:(-
1
г (1 - a) J
U(g)
(t - а)'
:da - 9£tu(ff),
(2)
где, Г(■) - гамма-функция Эйлера; й = ^й.
Эффект памяти выражается в зависимости текущего состояния динамической системы от предыдущих состояний, то есть нелокальности. Интенсивность этого эффекта будет определяться значением показателя степени РБ. Классически это некое а значение является нецелым и постоянным [15]. В науке и технике свойства памяти широко распространены в сложных системах [16]. Однако
t
t
о
существуют обобщения дробных производных на случай переменной во времени нелокальности о("Ь) и других функциональных зависимостей [17]. Дробные производные переменного порядка (УОРБ) находит свое применение в теории и практике физико-математических наук. Например, диффузионные уравнения [18,19] или модель динамического вязкоупругого осциллятора [20].
В работе [21] подробно исследуется обобщение РБ вида (2) до УОРБ при 0 < о(1") < 1.
1 Г1 й(ст) а(1) , л /0ч
^ = Э01 )й( а), (3)
Г(1 - a(t))J
о (1 - а)а(1)
где, а("Ь) € С1 (0,1) - непрерывно-дифференцируемая (гладкая) функция.
В таком случае, уравнение (1) в области П = |("Ь) : 0 < " < Т} можно представить как дробное уравнение (РБЕ) вида:
ЭО^Уа) = ?(й("),"), (4)
где, й("Ь) € С2(0,Т) - дважды непрерывно-дифференцируемая функция решения, ?(й("Ь),"Ь) - функционал, определяющий вид РБЕ.
Цикл работ автора данного исследования, включающий публикации [22-24], посвящён исследованию вопросов, связанных с математическим моделированием объемной активности газа радона (КУЛ) в накопительных камерах с газоразрядными датчиками. Накопительные камеры располагаются в пунктах наблюдения, которые находятся над зонами динамического влияния разломов в земной коре. Эти зоны имеют повышенную проницаемость, что способствует стоку подпочвенных газов в атмосферу [25]. Предполагается, что параметр о("Ь) описывает фрактальность геосреды и связан с такими её характеристиками как пористость, проницаемость и трещиноватость, которые прямо влияют на интенсивность процесса. В результате был разработан ряд математических моделей на основе нелинейных РБЕ вида (4). В [22] показано, что нелинейные модели с РБ вида (2) постоянного порядка о могут хорошо описывать накопительные режимы КУЛ. В работе [23] показано, что обобщение до УОРБ вида (3) на случай о("Ь) позволяет также описывать некоторые аномальные режимы КУЛ. Непрерывный мониторинг вариаций КУЛ с целью выявления аномалий в ее значениях является одним из эффективных методов исследования напряженно-деформированного состояния геосреды. КУЛ считается информативным и оперативным предвестником сейсмических событий, что обуславливает актуальность таких исследований.
При решении многих фундаментальных и прикладных задач, связанных с исследованием определенных свойств динамических систем, часто возникают ситуации, когда параметры модели точно не определены. Однако если известны экспериментальные данные или информация о диапазонах значений параметров модели, то можно осуществлять подбор различных параметров таких моделей эмпирически, учитывая имеющееся понимание природы моделируемого процесса. В работах [22-24] модельные параметры уточняются путем перебора значений и сопоставления временных рядов: результатов моделирования и
экспериментальных данных, представляющих процесс. Этот процесс продолжается до тех пор, пока результаты моделирования не начинают качественно аппроксимировать данные, то есть до достижения высоких значений корреляции и И2 детерминации [26]. Такой подход является трудоемким, что неизбежно приводит нас к идеям различных способов автоматизации подбора оптимальных параметров.
Например, в исследовании [27] для дробных уравнений применяют численные методы получения интервальных оценок решений по известным интервальным значениям параметров для задач аномальной диффузии.
Для решения проблемы поиска оптимальных параметров можно решать обратную задачу. Такие задачи часто являются некорректно поставленными [?], т. е. решения может не существовать, решение может быть не единственно или решение неустойчиво [29]. Однако обратная задача — это распространенный тип задач во многих научных областях [30], где необходимо определить значения параметров модели на основе наблюдаемых данных [31]. Несмотря на то что этот раздел математики стал активно развиваться только в Х1Х-ХХ веке, его по праву можно считать одним из наиболее важных. Необходимость такого подхода часто возникает при работе с геологическими данными [32], в геофизике и сейсмологии [33], в задачах компьютерной томографии [34,35] и многих других. Особенно остро этот вопрос стоит для геофизики, так как невозможно провести прямые измерения глубоко под землей [36].
В данной работе исследуется обратная задача идентификации функции о(Ь) в РБЕ вида (4). Статья имеет следующую структуру: сначала определяется прямая задача и её разностный аналог для неоднородного дробного уравнения с дробной производной Герасимова-Капуто 0 < о(Ь) < 1 и описывается численный метод её решения; далее на основе прямой задачи определяется соответствующая обратная задача и её разностный аналог по определению вида функции о(Ь); после этого для решения поставленной обратной задачи описывается итерационный метод безусловной оптимизации Левенберга-Марквардта ньютоновского типа; далее по пунктам описывается алгоритм Левенберга-Марквардта; далее приводятся тестовые примеры и результаты работы описанного алгоритма, где на графиках полученные результаты сравниваются с эталонным решением прямой задачи и экспериментальными данными; в заключении делаются выводы о проделанной работе и определяются возможные направления дальнейших исследований.
Постановка прямой задачи
Рассмотрим случай когда функционал ?(й(Ь),Ь) как неоднородное РБЕ, для которого в области П = задачу Коши вида:
30,"1)й(а) = 1, й(0) = й0
где,
• Э0"")й(а) - РБ переменного порядка 0 < о(Ь) < 1 вида (3);
1 . Тогда (4) определяется {(Ь) : 0 < Ь < Т} рассмотрим
• u(t) - функция решения;
• a(t) - показатель степени дробной производной (3),
• u(t) G C2(0,T) = U - класс дважды непрерывно-дифференцируемых функций;
• a(t) G C(0,1) = A - класс непрерывно-дифференцируемых функций;
• u0 - некоторая константа, начальное условие задачи Коши.
Определение 1. Прямая задача представляет собой задачу Коши (5), состоящую в поиске u(t) G U при известной a(t) G A.
Для решения прямой задачи (5) воспользуемся ранее разработанной и теоретически обоснованной нелокальной неявной конечно-разностной схемой (IFDS), которая была апробирована на ряде тестовых задач [21]. Пусть задана равномерная сеточная область О., в которой определены классы сеточных функций U и A с шагом дискретизации т:
т = T/N, Q = {(ti = гт) : 0 < i < N}, U G ГГ, A G Q,
(6)
u(t) = u(ti) = ui, a(t) = a(ti) = аг, 0 < аг < 1,
Неоднородное FDE (5) может быть аппроксимировано в области (6) по схеме IFDS следующим образом:
i-1
A^ wi (ui_j - ui-j-i) - 1 = 0, u = C,
j=l (7)
Ai = rfT ai л, wi = (j + 1)1-ai - j1-ai, 1 < i<N,
I (2 - ai) j
где C - известная константа.
Определение 2. Разностная прямая задача (7) состоит в том, чтобы в области (6) найти сеточную u(ti) G U при известных значениях сеточной функции a(ti) G A.
IFDS схема (7) может решаться методом Ньютона или модифицированным методом Ньютона. В работе [21] были сформулированы и доказаны ряд теорем о сходимости и устойчивости IFDS для FDE более общего вида (4), а также аппроксимации дробной производной (3). Подробнее о численном методе можно узнать из работ [21,37].
Замечание 1. Разностная прямая задача (7) на основе схемы IFDS и решаемая
методом Ньютона, согласно [21], безусловно устойчива и сходится с порядком
Q /2-mfx(a(t))
Постановка обратной задачи
Пусть значения функции а("Ь) С А (или её сеточного аналога а("Ц) € А) неизвестны. Однако о решении разностной прямой задачи Коши (7) в области П
известна дополнительная информация (экспериментальные данные):
и(ъ) = еио, (8)
Пусть а("") € А - функция некоего известного класса, и её вид однозначно определяется неким набором параметров. Тогда, согласно книге [38] Тихонова А. Н., Самарского А. А., решение обратной задачи сводится к поиску значений этих параметров. Следуя этому подходу, определим множество векторов X - как пространство решений обратной задачи (5). Вектор X = [Хо,..., ХК-|] - набор неизвестных параметров характеризующих вид функции а("Ь), где X € X, X с КК, К - число компонент. Тогда:
Определение 3. Обратная задача для (5) определяется как восстановление (определение) функции а("Ь) = а( X) € А:
Э0а(* 'и(а) = 1, и(0) = ео, (9)
Определение 4. Следовательно, разностная обратная задача для (7) определяется как восстановление сеточной функции а ("О = а( X) = а(Х0,..., ХК-1) € А, по известной (8) дополнительной информации:
г—1
л^ ^ (иг—) — иг—5—1) — 1 = 0, ио = ео, иг = ег,
77^ (10)
Лг = ^-7-^, Ц = (з + 1)1—(а^))— ]1—(а^)), 1 < г < N,
Г (2 — (а(Х))) 5
Пусть е = [е0,...,еN—1 ] - вектор экспериментальных данных согласно (8), а € и, ) = X) = [ш0, 1] - вектор модельных данных, т. е. решение
разностной прямой задачи (7), полученное относительно приближения а ^Х^ при
некоторых заданных значениях X.
Пусть — = п( X) = [п0, .",Лм-1] - вектор невязки размерности N > К, такой что КК —> т. е. вещественное подмножество К, зависящее от решения X) размерности N и, соответственно, зависящее от а( X) от К параметров X:
— = — — X), (11)
Тогда, в терминах теории безусловной оптимизации, согласно [39,40], решение разностной обратной задачи (10) сводится к минимизации ¥ функционала невязки:
ш ¥ (X) = М(л(XЛ , п : КК-> М : К,
/ -> х 1 1 N—1 1 N—1 (12)
м (Л(Х)) = ^ —Пт — = 2 ^ п2 =1X (ег — ^)2,
г=0 г=0
Замечание 2. Так как непрерывная обратная задача дискретизируется, часто случается, что нарушается именно условие устойчивости полученного
численным методом, определяя разностную обратную задачу (10) как некорректно поставленную [41]. Однако в рассмотренном случае для решения прямой задачи (5) в дискретной постановке (7) используется безусловно устойчивый алгоритм ^Б [21].
Метод минимизации целевой функции
Для решения такой задачи минимизации (12) обратимся к методам математической оптимизации [40, 42]. Существует класс методов называемый итерационными методами спуска или градиентного спуска [43,44], но у них есть неприятная особенность. На практике значения тП могут описывать «овражную» структуру, т. е. иметь локальные минимумы, помимо искомого глобального минимума. И бывает так, что в зависимости от выбора начального приближения метод может сойтись именно к локальному минимуму функционала ¥( X), а не к глобальному.
Как альтернатива итерационным методам спуска, существует класс ньютоновских (квазиньютоновских) методов, одних из самых эффективных на практике [39]. Можно воспользоваться методом безусловной оптимизации Ньютона, но он требует, чтобы функционал X) € С2(С С Ек) т. е. был дважды непрерывно-дифференцируем в открытом выпуклом множестве С. Данное условие требуется для расчёта Гессиана (матрицы Гессе) вида:
Н = У2¥(Х)= У2¥(Х= ддх^дХ^), I > 0, к < к - 1,
где, У¥(Х)
- градиент ^ в X. Подробнее в [39,40].
Метод безусловной оптимизации Ньютона реализуется в виде итерационной процедуры [40]:
АХ = (-Н-1) х )'
где, приближение X - восстанавливаемые значения, вычисленные на текущей итерации процедуры; АХ = [АХ0... АХк-1 ] - оптимальные приращения X для следующей итерации.
Замечание 3. Важной характеристикой численных методов оптимизации является скорость сходимости, характеризующая эффективность используемого метода для решения обратной задачи. Поэтому желательно использовать метод с самой высокой скоростью сходимости.
Поэтому воспользуемся модификацией метода Ньютона, называемой итерационным методом безусловной оптимизации Левенберга-Марквардта [45-47], который позволяет избавиться от условия на существование вторых производных. Полагая что как минимум выполняется ¥( X) € С1 (С С Ек), и X) = |т х тП, а матрица Гессе размерности к х к рассчитывается следующим образом:
Н = |т х | + уЕ,
тогда метод Левенберга-Марквардта реализуется в виде следующей итерационной процедуры для системы линейных алгебраических уравнений:
АХ = (- (Г х | + уЕ)-1) х (Г х т?) , (13)
где,
• Е - единичная матрица размерности К х К;
• т| = п = п ^X(n)J - вектор невязки;
• Т = Т ^X ^ = Т ^ X - матрица Якоби размерности N х К, элементы которой имеют вид:
зл(п)
|1к = , 1 = 0..N - 1, к = 0..К - 1,
тг,к эхkn), , ,
причём производные в Тг,к аппроксимируются конечной разностью, вида:
(п)6 (п)
т = п( ) - п() (14)
т1,к = ^к , (14)
где,
• - X = [^0... 6X^1 - заданные малые приращения X;
(п) (п)
П(п) = [п0 ).. 11 - вычисляется по формуле (11) на п итерации метода, на
/(п) у(п) т ^о ,..., XK—1
основе решения разностной прямой задачи относительно X ( п) = ^0^, ...^К-^ -
приближения;
• аналогично получаем - (п)6 = [Л0п)г.. . ЛгП-'1] для решения прямой задачи относительно X (п)6 = X (п) + X X - приближения с заданными приращениями.
Замечание 4. Как правило, так называемый параметр регуляризации у -некоторое число, принимающее разное значение при разных задачах оптимизации, т.е. при разной дополнительной информации 0 .
Параметр регуляризации у имеет ключевое значение для метода Левенберга-Марквардта и напрямую связан с шагом и направлением сходимости метода. Поэтому стоит отметить важные свойства, накладываемые на у, согласно [40]:
1. ¥ = Е>о = {у|у € К, у > 0}
2. если у € ¥ и матрица Гессе Н положительно определена, то АX - направление спуска (градиента) для «хорошего» (оптимального) шага метода;
3. если у —> 0 очень мало, то АX - оптимальный шаг метода в направлении градиента. Это также означает, что скорость сходимости (13) метода близка к квадратичной;
4. если у € ¥,у » 0, то из (13) следует что АX ^ -^ (Тт х тП) - шаг метода в направлении антиградиента.
Алгоритм итерационного метода Левенберга-Марквардта
1) На п = 0 итерации, передаём методу:
1.а) X0 - начальное приближение для компонент X; 1.Ь) _ X - заданные малые приращения; 1.с) е - дополнительную информацию;
1^) у(0) = V • тах Х^ - стартовое значение параметра регуляризации, где: V - заданная константа, —Х - диагональ матрицы В(0) = I х I (X ( ^;
1.е) с - заданная константа для пересчёта параметра регуляризации.
2) На п = 0 итерации, вычисляем з0 - начальное значение функционала невязки:
2.а) Решаем разностную прямую задачу (7) относительно - начального приближения, получая о—0) - вектор модельных данных;
2.Ь) Вычисляем согласно (11) компоненты вектора невязки: п = е — о(0);
_V _V 1 N—1 ^
2.с) Вычисляем согласно (12) значение: з0 = ¥ (Xх) = М = ^ ^
2 г=0
3) На итерации п > 0, вычисляем значения матрицы Якоби I = I (X'
относительно текущего приближения X(n':
3.а) Решаем разностную прямую задачу (7) относительно X(n)г = X(n' + _ X приближения с заданными приращениями, получая о(п)г;
3.Ь) Вычисляем по формуле (11) компоненты вектора невязки: П = е — о(п)г.
3.с) Получаем компоненты вектора невязки: цп = е — о(п). Так как о(п) будет уже получена на п — 1 итерации, то:
(при п = 0): _Х - известна с шага (2.Ь) этого алгоритма;
(при п > 1): п п = е — о(п) = е — о(п 1)А - известна с шага (4.с);
3.d) Вычисляем J (X(n)) по формуле (14).
4) На итерации п > 0, вычисляем в1 - текущее значение функционала невязки: 4.а) Вычисляем ДX - оптимальные приращения X, решая систему (13);
4.Ь) Решаем прямую задачу (7) относительно X(n)А = X(n' + AX - приближения с вычисленными оптимальными приращениями, получая о(п)А;
4.с) Вычисляем по формуле (11) компоненты вектора: п п = е — о(п)А;
4.а) Вычисляем согласно (12) значение: в1 = ¥ (X (п) = М (п)А^.
5) Проверка условий выхода из цикла итерационной процедуры:
5.а) Метод сходится к оптимальному решению - если е < I, где I - заданная точность решения разностной обратной задачи, е - среднеквадратичная ошибка модельных данных:
N—1
1 'N—' г
i=0
< = n^K14
5.Ь) Метод сходится к локальному минимуму - если не даёт существенного
изменения значений ДX в ходе итераций. Иначе говоря А ^АXJ —> 0, где А ^А^ - оценка скорости «приращения приращений» определяемая:
д И=£ К—^r i)
\ k=0 /
0.
6) Проверка условий продолжения цикла. Если условия (5.а) и (5.Ь) не выполнены,
то возможны два сценария:
6.а) Если з0 > 51 - считая что алгоритм вычислил оптимальный шаг, то:
§0 = $1;
X(n) = X(n)А;
г = *;
п = п + 1;
Переходим к (3)-му шагу представленного алгоритма.
6.Ь) Иначе, если з0 < в1 - алгоритм вычислил не оптимальный шаг, то:
У = су;
Переходим к (4)-му шагу представленного алгоритма.
Представленный алгоритм для решения обратной задачи (10) был реализован в виде набора подпрограмм на языке С в ходе работ по гранту Российского научного фонда № 23-71-01050 по теме «Разработка программного комплекса для моделирования и анализа объемной активности радона как предвестника сильных землетрясений Камчатки».
Замечание 5. Особенностью описанных методов ньютоновского типа минимизации функционала (12) для решения разностной обратной задачи является то, что необходимо множество раз решать разностную прямую задачу (7), при различных значениях X. Это требует использования самых эффективных реализаций прямой задачи (7) и методов её решения. Для этого были использованы ранее реализованные автором подпрограммы [37] на языке С, с использованием ОрепМР [48] и/или СиЭЛ [49], для ускорения расчётов за счёт распараллеливания.
Тестовые примеры и результаты численных экспериментов
Покажем на примерах, что описанный метод решения разностной обратной задачи (10), при определённых заданных X (о), 6X, c,v параметрах связанных с направлением и шагом метода Левенберга-Марквардта, действительно позволяет восстановить вид функциональной зависимости a(t). Структура численного эксперимента такова:
• Положим для всех примеров далее: T = 30, N = 300;
• Решим «эталонную» разностную прямую задачу (7) при всех заданных и известных параметрах модели, в том числе и восстанавливаемой функции a(t) и определяющих её параметрах X = [X0,..., XK-1]. Синяя (Blue) модельная кривая на (рис.1-4 b);
• На основе «эталонного» решения сгенерируем 0 - псевдослучайные экспериментальные данные (8). Чёрная (Black) кривая на (рис.1-5 b);
• Решим разностную обратную задачу (10) описанным алгоритмом метода Левенберга-Марквардта и восстановив значения a(t) = а( X) = a(X0, ...,XK-1) € A по данным 0 . Красная (Red) и/или зелёная (Green) модельная кривая на (рис.1-5 b);
Замечание 6. На (рис.1-4) приведены только наиболее удачные результаты. Это связано с тем, что в ходе численных экспериментов по подбору
X (0)
, 6X, c, v получен
огромный объем данных для анализа работы алгоритма с дробными уравнениями, многие из которых приводили к тому что метод разваливался.
Для генерации 0 будем использовать данные «эталонного» решения, к которым, поэлементно, внесём псевдослучайное возмущение в диапазоне [-0.25,0.25], сгенерированное на компьютере:
0(tO = u(ti)+ p(s), i = 0..N - 1, s = ±0.25,
где, p(s) - функция распределения дискретной случайной величины s. Данный способ генерации экспериментальных данных реализован на языке программирования С на основе rand() - функции для генерации псевдослучайных чисел, стандартной библиотеки <stdlib.h>.
Расчёты для решения прямых и обратных задач в тестовых примерах проводились на персональной ЭВМ приобретенной в для работ по гранту Российского научного фонда № 23-71-01050 по теме «Разработка программного комплекса для моделирования и анализа объемной активности радона как предвестника сильных землетрясений Камчатки». Персональная ЭВМ имеет следующие характеристики системы: CPU - AMD Ryzen 9 7950X 16 х 4.5 GHz (32 Threads), cache L2 16 Mb & L3 64 Mb; RAM - 96 Gb; GPU - GeForce RTX 4090, 24 Gb, 2235 MHz, ALU 16384.
Рис. 1. Пример 1. a(t0 = X0 ^-T^J - линейно-возрастающая. a) Графики ct(t0:
(Blue) - эталонный X, (Red) - восстановленный X; b) Модельные кривые: (Blue) - эталонная, (Red) - рез. оптимизации при X(o) = [0.01], 6X = [0.005] [Figure 1. Example 1. ct(t0 = X0 is a linear-increasing. a) Graphs ct(t0: (Blue) - reference X, (Red) - restored X; b) Model curves: (Blue) - reference, (Red) -optimization result at X-^ = [0.01], deltaX = [0.005]]
Рис. 2. Пример 2. a(tO = X0
(N - i)T
T
- линейно-убывающая. a) Графики a(ti):
(Blue) - эталонный X, (Red) - восстановленный X; b) Модельные кривые: (Blue) - эталонная, (Red) - рез. оптимизации при X (0 = [0.01], 6X = [0.005]
[Figure 2. Example 2. a(ti) = X0
(N-i)T
linear-decreasing. a) Graphs a(ti): (Blue)
- reference X, (Red) - restored X; b) Model curves: (Blue) - reference, (Red) -optimization result at X^ = [0.01], deltaX = [0.005]]
tine cak: (Red) -14.0070 sec,, (Blue) - 0.006 sec. aJ algorithm DP: (Red) - parallel(OperiMP), CPU = 24 I (Blue) - parallel(OperiMP), CPU = 24
(Red) - /count = 62 /n =25 / step = 4/1 = 0.1 / c = 2 / v = 50 / [Xc = 0.498, Xn - 0.100, X= = 0, X3 = 0] (Blue) - / court = 0 / n = 0 / step = 0 / I = 0 / c = 0 / v = 0 / [Xc = 0.5, Xi = 0.105, X: = 0, X3 = 0]
0.50
2
Рис. 3. Пример 3. a(ti) = X0cos(X1iT) - периодическая функция. a) Графики a(ti): (Blue) - эталонный X, (Red) - восстановленный X; b) Модельные кривые: (Blue) - эталонная, (Red) - результат оптимизации при X (0) = [0.1,0.02], 6X = [0.05,0.01] [Figure 3. Example 3. a(ti) = X0cos (X1iT) - periodic function. a) Graphs a(ti): (Blue) - reference X, (Red) - restored X; b) Model curves: (Blue) - reference, (Red) - result optimization at X0 = [0.1,0.02], 6X = [0.05,0.01]]
2.99 2.15 1.31 0.48
rANfWw^^^
Ь)
3.00
6.00
1.00
12.00
15.00 t
18.00
21.00
24.00
27.00
30.00
2
Рис. 4. Пример 4. ct(t0 = Xosin (XiIt) - периодическая. a) Графики a(t0: (Blue) - эталонный X, (Red, Green) - восстановленные X; b) Модельные кривые:
(Blue) - эталонная, (Red, Green) - результаты оптимизации при X(o) = [0.005,0.001], 6X = [0.005,0.001] [Figure 4. Example 4. ct(t0 = X0sin (X1It) - periodic. a) Graphs a(t0: (Blue) -reference X, (Red, Green) - restored X; b) Model curves: (Blue) - reference, (Red,
Green) - results optimization at X-^ = [0.005,0.001], deltaX = [0.005,0.001]]
Заключение
Из результатов можно сделать вывод о возможности решения обратной задачи для нелинейного дробного уравнения с оператором дифференцирования типа Герасимова-Капуто как переменного, так и постоянного порядка. Показано, что с помощью методов математической безусловной оптимизации, в частности итерационного метода Левенберга-Марквардта, имея математическую модель и экспериментальные данные процесса, можно восстанавливать близкие к оптимальным значения a(t) и а показателя степени дробной производной. В том числе и функции a(t), зависящие от двух неизвестных параметров.
Однако, чем более сложную структуру имеют экспериментальные данные и вид функции, определяющий a(t), тем сложнее подобрать параметры c,v метода Левенберга-Марквардта и стартовые значения , 6X, при которых бы выполнялись условия, накладываемые на параметр регуляризации у.
Продолжение работы состоит в том, чтобы перейти от решения обратных задач для нелинейного дробного уравнения к обратным задачам для квадратично нелинейного дробного аналога уравнения Риккати. В свою очередь, это позволит решать различные обратные задачи по определению тех или иных параметров модельного уравнения объемной активности газа радона [23] на основе регистрируемых экспериментальных данных.
Рис. 5. Пример 5. a(t0 = ( X,
(N - i)T T
cos(XiIt) - периодическая функция
с линейно-убывающей амплитудой. а) Графики a(ti): (Blue) - эталонный X, (Red) - восстановленный X; b) Модельные кривые: (Blue) - эталонная, (Red) - результат оптимизации при X^ = [0.22,0.0404], 6X = [0.15,0.03]
[Figure 5. Example 5. ct(t0 = (X0 ( ^ i)T j j c 0S(x1It)2 is a periodic function
T
with a linear-decreasing amplitude. a) Graphs a(t0: (Blue) - reference X, (Red) -reconstructed X; b) Model curves: (Blue) - reference, (Red) - optimization result at
X^ = [0.22,0.0404], 6X = [0.15,0.03] ]
Аббревиатуры
FD Fractional Derivatives
VOFD Variable-Order Fractional Derivatives
IDE Integro-Differential Equations
RVA Radon Volumetric Activity
FDE Fractional Derivative Equation
IFDS Implicite Finite-Difference Sheme
Список литературы
1. Kilbas A.A., Srivastava H.M., Trujillo J.J. Theory and Applications of Fractional Differential Equations. Amsterdam: Elsevier, 2006. 540 pp. ISBN 9780444518323.
2. Iomin A. Fractional-time quantum dynamics, Physical Review E, 2009. vol. 80, no. 2, pp. 1-4 DOI: 10.1103/PhysRevE.80.022103.
3. Bagley R. L., Torvik P. J. A theoretical basis for the application of fractional calculus to viscoelasticity, Journal of rheology, 1983. vol.27, no. 3, pp. 201-210 DOI: 10.1122/1.549724.
4. Coimbra C. F. M. Mechanics with variable-order differential operators, Annalen der Physik, 2003. vol. 12, no. 11-12, pp. 692-703 DOI: 10.1002/andp.200310032.
5. Rossikhin Y. A., Shitikova M. V. Application of fractional calculus for dynamic problems of solid mechanics: novel trends and recent results, Applied Mechanics Reviews, 2010. vol.63, no. 1:010801, pp. 1-52 DOI: 10.1115/1.4000563.
6. Metzler R., Klafter J. The restaurant at the end of the random walk: recent developments in the description of anomalous transport by fractional dynamics, Journal of Physics A Mathematical and General, 2004. vol.37, no. 31, pp. 161-208 DOI: 10.1088/0305-4470/37/31/R01.
7. Moroz L.I., Maslovskaya A. G. Numerical Simulation of an Anomalous Diffusion Process Based on a Scheme of a Higher Order of Accuracy, Mathematical Models and Computer Simulations, 2021. vol. 13, no. 3, pp. 492-501 DOI: 10.1134/S207004822103011X.
8. Parovik R. I. Mathematical modeling of linear fractional oscillators, Mathematics, 2020. vol. 8, no. 11:1879, pp. 1-26 DOI: 10.3390/math8111879.
9. Sun H.G., Chen W., Wei H., Chen Y. Q.A comparative study of constant-order and variable-order fractional models in characterizing memory property of systems, The European Physical JournalSpecial Topics, 2011. vol.193, no. 1, pp. 185-192 DOI: 10.1140/epjst/e2011-01390-6.
10. Volterra V. Sur les équations intégro-différentielles et leurs applications, Acta Mathematica, 1912. vol. 35, no. 1, pp. 295-356 DOI: 10.1007/BF02418820.
11. Нахушев А. М. Дробное исчисление и его применение. Москва: Физматлит, 2003. 272 с. ISBN 5-9221-0440-3.
12. Рехвиашвили С.Ш., Псху А. В. Дробный осциллятор с экспоненциально-степенной функцией памяти, Письма в ЖТФ, 2022. Т. 48, №7, С. 33-35 DOI: 10.21883/PJTF.2022.07.52290.19137.
13. Gerasimov A. N. Generalization of linear deformation laws and their application to internal friction problems, Applied Mathematics and Mechanics, 1948. vol. 12, pp. 529-539.
14. Caputo M. Linear models of dissipation whose Q is almost frequency independent - II, Geophysical Journal International, 1967. vol.13, no. 5, pp. 529-539 DOI: 10.1111/j.1365-246X.1967.tb02303.x.
15. Uchaikin V. V. Fractional Derivatives for Physicists and Engineers. Vol. I. Background and Theory. Berlin: Springer, 2013.373 DOI: 10.1007/978-3-642-33911-0 pp.
16. Westerlund S.Dead matter has memory!, Physica Scripta, 1991. vol.43, no. 2, pp. 174-179 DOI: 10.1088/0031-8949/43/2/011.
17. Patnaik S., Hollkamp J. P., Semperlotti F. Applications of variable-order fractional operators: a review, Proceedings of the Royal Society A, 2020. vol.476, no. 2234, pp. 20190498 DOI: 10.1098/rspa.2019.0498.
18. Lin R., Liu F., Anh V., Turner I. W. Stability and convergence of a new explicit finite-difference approximation for the variable-order nonlinear fractional diffusion equation, Applied Mathematics and Computation, 2009. vol.212, no. 2, pp. 435-445 DOI: 10.1016/j.amc.2009.02.047.
19. Fang Z.W., Sun H. W., Wang H.A fast method for variable-order Caputo fractional derivative with applications to time-fractional diffusion equations, Computers & Mathematics with Applications, 2020. vol.80, no. 5, pp. 1443-1458 DOI: 10.1016/j.camwa.2020.07.009.
20. Sahoo S., Saha Ray S., Das S., Bera R. K.The formation of dynamic variable-order fractional differential equation, International Journal of Modern Physics C, 2016. vol. 27, no. 07, pp. 1650074 DOI: 10.1142/S0129183116500741.
21. Tverdyi D. A., Parovik R. I. Investigation of Finite-Difference Schemes for the Numerical Solution of a Fractional Nonlinear Equation, Fractal and Fractional, 2022. vol.6, no. 1:23, pp. 1-27 DOI: 10.3390/fractalfract6010023.
22. Tverdyi D. A., Parovik R. I., Makarov E. O., Firstov P. P. Research of the process of radon accumulation in the accumulating chamber taking into account the nonlinearity of its entrance, E3S Web Conference, 2020. vol.196, no. 02027, pp. 1-6 DOI: 10.1051/e3sconf/202019602027.
23. Tverdyi D. A., Makarov E. O., Parovik R. I. Hereditary Mathematical Model of the Dynamics of Radon Accumulation in the Accumulation Chamber, Mathematics, 2023. vol.11, no. 4:850, pp. 1-20 DOI: 10.3390/math11040850.
24. Tverdyi D.A., Makarov E. O., Parovik R. I. Research of Stress-Strain State of Geo-Environment by Emanation Methods on the Example of alpha(t)-Model of Radon Transport, Bulletin KRASEC. Physical and Mathematical Sciences, 2023. Т.44, №3, С. 86-104 DOI: 10.26117/2079-6641-2023-443-86-104.
25. Рудаков В. П. Эманационный мониторинг геосред и процессов. Москва: Научный мир, 2009.175 с. ISBN 978-5-91522-102-3.
26. Cox D. R. Hinkley D.V. Theoretical Statistics, 1st edition. London: Chapman & Hall/CRC, 1979.528 pp. ISBN 9780412161605.
27. Ревизников Д. А., Морозов А. Ю. Алгоритмы численного решения дробно-дифференциальных уравнений с интервальными параметрами, Сибирский журнал индустриальной математики, 2023. Т. 26, №4, С. 93-108 DOI: 10.33048/SIBJIM.2023.26.407.
28. Hadamard J. S. Sur les problèmes aux derivees partielles et leur significa tion physique, Princeton University Bulletin, 1902. vol. 13, no. 4, pp. 49-52.
29. Morozov V. A. Methods for Solving Incorrectly Posed Problems. Springer: New York, 1984. 257 pp. ISBN 978-0-387-96059-3 DOI: 10.1007/978-1-4612-5280-1.
30. Mueller J.L., Siltanen S. Linear and Nonlinear Inverse Problems with Practical Applications. Philadelphia, USA: Society for Industrial and Applied Mathematics, 2012.351 pp. ISBN 978-1-61197233-7 DOI: 10.1137/1.9781611972344.
31. Tarantola A. Inverse problem theory : methods for data fitting and model parameter estimation. Amsterdam and New York: Elsevier Science Pub. Co., 1987. 613 pp. ISBN 0444427651.
32. Tahmasebi P., Javadpour F., Sahimi M. Stochastic shale permeability matching: Three-dimensional characterization and modeling, International Journal of Coal Geology, 2016. vol. 165, no. 1, pp. 231242 DOI: 10.1016/j.coal.2016.08.024.
33. Lailly P. The seismic inverse problem as a sequence of before stack migrations, Conference on Inverse Scattering, Theory and application, 1983, pp. 206-220.
34. Mohamad-Djafari A. Inverse Problems in Vision and 3D Tomography. New-York: ISTE-WILEY, 2010. 480 pp. ISBN 9781848211728 DOI: 10.1002/9781118603864.
35. Hayotov A.R., Jeon S., Shadimetov K. M. Application of optimal quadrature formulas for reconstruction of CT images, Journal of Computational and Applied Mathematics, 2021. vol. 388, pp. 113313 DOI: 10.1016/j.cam.2020.113313.
36. Gubbins D. Book reviews. Inverse Problem Theory. Methods for Data Fitting and Model Parameter Estimation Albert Tarantola. Elsevier, Amsterdam and New York, 1987, Geophysical Journal International, 1988. vol.94, no. 1, pp. 167-168 DOI: 10.1111/j.1365-246X.1988.tb03436.x.
37. Tverdyi D.A., Parovik R. I. Hybrid GPU-CPU Efficient Implementation of a Parallel Numerical Algorithm for Solving the Cauchy Problem for a Nonlinear Differential Riccati Equation of Fractional Variable Order, Mathematics, 2023. vol.11, no. 15:3358, pp. 1-21 DOI: 10.3390/math11153358.
38. Тихонов А.Н., Самарский А. А. Уравнения математической физики. Москва: Наука, 1977. 736 с.
39. Dennis J.E., Robert Jr., Schnabel B. Numerical methods for unconstrained optimization and nonlinear equations. Philadelphia: SIAM, 1983.378 pp.
40. Иващенко Д. С. Численные методы решения прямых и обратных задач для уравнения диффузии дробного порядка по времени, дисс. ... канд. физ.-матем. наук. Томск, 2008.187 с.
41. Тихонов А. Н. О решении некорректно поставленных задач и методе регуляризации, Докл. АН СССР, 1963. Т. 151, №3, С. 501-504.
42. Кабанихин С. И., Искаков К. Т. Оптимизационные методы решения коэффициентных обратных задач. Новосибирск: Новосибирский государственный университет, 2001.315 с. ISBN 5-94356-022-Х.
43. Калиткин Н.Н. Численные методы. 2-е изд.. Санкт-Петербург: БХВ, 2011.592 с. ISBN 9785-9775-0500-0.
44. Arridge S. R., Schweiger, M. A General Framework for Iterative Reconstruction Algorithms in Optical Tomography, Using a Finite Element Method, Computational Radiology and Imaging: Therapy and Diagnostics, 1902. vol. 13, no. 4, pp. 40-70 DOI: 110.1007/978-1-4612-1550-9_4..
45. Levenberg K. A method for the solution of certain non-linear problems in least squares, Quarterly of applied mathematics, 1944. vol.2, no. 2, pp. 164-168 DOI: 10.1090/qam/10666.
46. Marquardt D. W. An algorithm for least-squares estimation of nonlinear parameters, Journal of the society for Industrial and Applied Mathematics, 1963. vol.11, no. 2, pp. 431-441 DOI: 10.1137/0111030.
47. More J. J. The Levenberg-Marquardt algorithm: Implementation and theory, In: Watson, G.A. (eds) Numerical Analysis. Lecture Notes in Mathematics, 1978. vol.630, pp. 105-116 DOI: 10.1007/BFb0067700.
48. Борзунов С. В., Кургалин С. Д., Флегель А. В. Практикум по параллельному программированию: учебное пособие. Санкт-Петербург: БХВ, 2017.236 с. ISBN 978-59909805-0-1.
49. Sanders J., Kandrot E. CUD A by Example: An Introduction to General-Purpose GPU Programming. London: Addison-Wesley Professional, 2010.311 pp. ISBN 978-0-13-138768-3.
Информация об авторах
Твёрдый Дмитрий АлександровичА физико-математических наук,
лаборатории элетромагнитных космофизических исследований
радиоволн ДВО РАН, с. © СЖСГО 0000-0001-6983-5258.
- кандидат научный сотрудник
излучений, Институт и распространения
Паратунка, Россия,
Паровик Роман Иванович А - доктор физико-математических наук, доцент, ведущий научный сотрудник лаборатории моделирования физических процессов института космофизических исследований и распространения радиоволн ДВО РАН, с. Паратунка, Россия, © СЖСГО 0000-0002-1576-1860.
References
[1 [2
[3 [4 [5
[6
[7
[8] [9]
10' 11 12
Kilbas A.A., Srivastava H.M., Trujillo J.J. Theory and Applications of Fractional Differential Equations. Amsterdam, Elsevier, 2006, 540 pp., isbn: 9780444518323. Iomin A. Fractional-time quantum dynamics, Physical Review E, 2009, vol. 80, no. 2, pp. 1-4. DOI: 10.1103/PhysRevE.80.022103.
Bagley R. L., Torvik P.J. A theoretical basis for the application of fractional calculus to viscoelasticity, Journal of rheology, 1983, vol. 27, no. 3, pp. 201-210. DOI: 10.1122/1.549724. Coimbra C. F. M. Mechanics with variable-order differential operators, Annalen der Physik, 2003, vol. 12, no. 11-12, pp. 692-703. DOI: 10.1002/andp.200310032. Rossikhin Y. A., Shitikova M.V. Application of fractional calculus for dynamic problems of solid mechanics: novel trends and recent results, Applied Mechanics Reviews, 2010, vol. 63, no. 1:010801, pp. 1-52. DOI: 10.1115/1.4000563.
Metzler R., Klafter J. The restaurant at the end of the random walk: recent developments in the description of anomalous transport by fractional dynamics, Journal of Physics A Mathematical and General, 2004, vol. 37, no. 31, pp. 161-208. DOI: 10.1088/0305-4470/37/31/R01.
Moroz L. I., Maslovskaya A.G. Numerical Simulation of an Anomalous Diffusion Process Based on a Scheme of a Higher Order of Accuracy, Mathematical Models and Computer Simulations, 2021, vol. 13, no. 3, pp. 492-501. DOI: 10.1134/S207004822103011X. Parovik R. I. Mathematical modeling of linear fractional oscillators, Mathematics, 2020, vol. 8, no. 11:1879, pp. 1-26. DOI: 10.3390/math8111879.
Sun H.G., Chen W., Wei H., Chen Y. Q. A comparative study of constant-order and variable-order fractional models in characterizing memory property of systems, The European Physical Journal-Special Topics, 2011, vol. 193, no. 1, pp. 185-192. DOI: 10.1140/epjst/e2011-01390-6.
Volterra V. Sur les equations integro-differentielles et leurs applications, Acta Mathematica, 1912, vol. 35, no. 1, pp. 295-356. DOI: 10.1007/BF02418820.
Nahushev A. M. Drobnoe ischislenie i ego primenenie [Fractional calculus and its application]. Moscow: Fizmatlit, 2003, 272 pp., isbn: 5-9221-0440-3 (In Russian) Rekhviashvili S.S., Pskhu A.V. Drobnyj oscillyator s eksponencial'no-stepennoj funkciej pamyati [Fractional oscillator with exponential-power memory function], Pis'ma v ZHTF [Letters to ZhTF], 2022, vol. 48, no. 7, pp. 33-35. DOI: 10.21883/PJTF.2022.07.52290.19137,(In Russian)
Gerasimov A. N. Generalization of linear deformation laws and their application to internal friction problems, Applied Mathematics and Mechanics, 1948, vol. 12, pp. 529-539. Caputo M. Linear models of dissipation whose Q is almost frequency independent - II, Geophysical Journal International, 1967, vol. 13, no. 5, pp. 529-539. DOI: 10.1111/j.1365-246X.1967.tb02303.x.
Uchaikin V. V. Fractional Derivatives for Physicists and Engineers. Vol. I. Background and Theory. Berlin, Springer, 2013, 373 pp. DOI: 10.1007/978-3-642-33911-0. Westerlund S. Dead matter has memory!, Physica Scripta, 1991, vol. 43, no. 2, pp. 174-179. DOI: 10.1088/0031-8949/43/2/011.
Patnaik S., Hollkamp J. P., Semperlotti F. Applications of variable-order fractional operators: a review, Proceedings of the Royal Society A, 2020, vol. 476, no. 2234, pp. 20190498. DOI: 10.1098/rspa.2019.0498.
Lin R., Liu F., Anh V., Turner I.W. Stability and convergence of a new explicit finite-difference approximation for the variable-order nonlinear fractional diffusion equation,
Applied Mathematics and Computation, 2009, vol. 212, no. 2, pp. 435-445. DOI: 10.1016/j.amc.2009.02.047.
[19] Fang Z. W., Sun H. W., Wang H. A fast method for variable-order Caputo fractional derivative with applications to time-fractional diffusion equations, Computers & Mathematics with Applications, 2020, vol. 80, no. 5, pp. 1443-1458. DOI: 10.1016/j.camwa.2020.07.009.
[20] Sahoo S., Saha Ray S., Das S., Bera R. K. The formation of dynamic variable-order fractional differential equation, International Journal of Modern Physics C, 2016, vol. 27, no. 07, pp. 1650074. DOI: 10.1142/S0129183116500741.
[21] Tverdyi D. A., Parovik R. I. Investigation of Finite-Difference Schemes for the Numerical Solution of a Fractional Nonlinear Equation, Fractal and Fractional, 2022, vol. 6, no. 1:23, pp. 1-27. DOI: 10.3390/fractalfract6010023.
[22] Tverdyi D.A., Parovik R.I., Makarov E.O., Firstov P.P. Research of the process of radon accumulation in the accumulating chamber taking into account the nonlinearity of its entrance, E3S Web Conference, 2020, vol. 196, no. 02027, pp. 1-6. DOI: 10.1051/e3sconf/2020196020278.
[23] Tverdyi D.A., Makarov E.O., Parovik R.I. Hereditary Mathematical Model of the Dynamics of Radon Accumulation in the Accumulation Chamber, Mathematics, 2023, vol. 11, no. 4:850, pp. 1-20. DOI: 10.3390/math11040850.
[24] Tverdyi D.A., Makarov E. O., Parovik R.I. Research of Stress-Strain State of Geo-Environment by Emanation Methods on the Example of alpha(t)-Model of Radon Transport, Bulletin KRASEC. Physical and Mathematical Sciences, 2023, vol. 44, no. 3, pp. 86-104. DOI: 10.26117/2079-6641-2023-44-3-86-104.
[25] Rudakov V. P. Emanacionnyj monitoring geosred i processov [Emanational monitoring of geoenvironments and processes]. Moscow: Science World, 2009, 175 pp., isbn: 978-5-91522102-3 (In Russian)
[26] Cox D. R. Hinkley D. V. Theoretical Statistics, 1st edition. London, Chapman & Hall/CRC, 1979, 528 pp., isbn: 9780412161605.
[27] Reviznikov D.L., Morozov A.Y. Algorithms for the numerical solution of fractional differential equations with interval parameters, Siberian journal of industrial mathematics, 2023, vol. 26, no. 4, pp. 93-108. DOI: 10.33048/SIBJIM.2023.26.407.
[28] Hadamard J. S. Sur les problemes aux derivees partielles et leur significa tion physique, Princeton University Bulletin, 1902, vol. 13, no. 4, pp. 49-52.
[29] Morozov V. A. Methods for Solving Incorrectly Posed Problems, New York: Springer, 1984, 257 pp., isbn: 978-0-387-96059-3, DOI: 10.1007/978-1-4612-5280-1.
[30] Mueller J. L., Siltanen S. Linear and Nonlinear Inverse Problems with Practical Applications. Philadelphia, USA, Society for Industrial and Applied Mathematics, 2012, 351 pp., isbn: 978-1-61197-233-7. DOI: 10.1137/1.9781611972344.
[31] Tarantola A. Inverse problem theory : methods for data fitting and model parameter estimation, Amsterdam and New York: Elsevier Science Pub. Co., 1987, 613 pp., isbn: 0444427651.
[32] Tahmasebi P., Javadpour F., Sahimi M. Stochastic shale permeability matching: Three-dimensional characterization and modeling, International Journal of Coal Geology, 2016, pp. 231-242, vol. 165, no. 1. DOI: 10.1016/j.coal.2016.08.024.
[33] Lailly P. The seismic inverse problem as a sequence of before stack migrations, Conference on Inverse Scattering, Theory and application, 1983, pp. 206-220.
[34] Mohamad-Djafari A. Inverse Problems in Vision and 3D Tomography. New-York, ISTE-WILEY, 2010, 480 pp., isbn: 9781848211728. DOI: 10.1002/9781118603864.
[35] Hayotov A. R., Jeon S., Shadimetov K. M. Application of optimal quadrature formulas for reconstruction of CT images, Journal of Computational and Applied Mathematics, 2021, vol. 388, pp. 113313. DOI: 10.1016/j.cam.2020.113313.
[36] Gubbins D. Book reviews. Inverse Problem Theory. Methods for Data Fitting and Model Parameter Estimation Albert Tarantola. Elsevier, Amsterdam and New York, 1987, Geophysical Journal International, 1988, vol. 94, no. 1, pp. 167-168. DOI: 10.1111/j.1365-246X.1988.tb03436.x.
[37] Tverdyi D.A., Parovik R.I. Hybrid GPU-CPU Efficient Implementation of a Parallel Numerical Algorithm for Solving the Cauchy Problem for a Nonlinear Differential Riccati Equation of Fractional Variable Order, Mathematics, 2023, vol. 11, no. 15:3358, pp. 1-21. DOI: 10.3390/math11153358.
[38] Tihonov A.N., Samarskij A.A. Uravneniya matematicheskoj fiziki [Mathematical physics equations]. Moskva: Nauka, 1977, 736 pp. (In Russian)
[39] Dennis J.E., Robert Jr., Schnabel B. Numerical methods for unconstrained optimization and nonlinear equations. Philadelphia, SIAM, 1983, 378 pp.
[40] Ivashchenko D. S. Chislennye metody resheniya pryamyh i obratnyh zadach dlya uravneniya diffuzii drobnogo poryadka po vremeni [Numerical methods of solution of direct and inverse problems for the fractional order time-dependent diffusion equation]. Diss. ... PhD (Phys.-Math.). Tomsk. 2008. 187 p. (In Russian)
[41] Tihonov A. N. O reshenii nekorrektno postavlenny'x zadach i metode regulyarizacii [On the solution of ill-posed problems and the method of regularization], Dokl. Akad. Nauk SSSR [Report AS USSR], 1963, vol. 151, no. 3, pp. 501-504 (In Russian)
[42] Kabanihin S.I., Iskakov K.T. Optimizacionnye metody resheniya koefficientnyh obratnyh zadach [Optimization methods for solving coefficient inverse problems]. Novosibirsk: Novosibirskij gosudarstvennyj universitet, 2001, 315 pp., isbn: 5-94356-022-X (In Russian)
[43] Kalitkin N. N. Chislennye metody. 2-e izd. [Numerical methods. 2nd ed.]. Saint Petersburg: BVH, 2011, 592 pp., isbn: 978-5-9775-0500-0 (In Russian)
[44] Arridge S. R., Schweiger, M. A General Framework for Iterative Reconstruction Algorithms in Optical Tomography, Using a Finite Element Method, Computational Radiology and Imaging: Therapy and Diagnosticsn, 1999, vol. 110, pp. 40-70. DOI: 110.1007/978-1-4612-1550-9_4.
[45] Levenberg K. A method for the solution of certain non-linear problems in least squares, Quarterly of applied mathematics, 1944, vol. 2, no. 2, pp. 164-168. DOI: 10.1090/qam/10666.
[46] Marquardt D. W. An algorithm for least-squares estimation of nonlinear parameters, Journal of the society for Industrial and Applied Mathematics, 1963, vol. 11, no. 2, pp. 431-441. DOI: 10.1137/0111030.
[47] More J. J. The Levenberg-Marquardt algorithm: Implementation and theory, In: Watson, G.A. (eds) Numerical Analysis. Lecture Notes in Mathematics, 1978, vol. 630, pp. 105-116. DOI: 10.1007/BFb0067700.
[48] Borzunov S.V., Kurgalin S.D., Flegel A.V. Praktikum po parallel'nomu programmirovaniyu: uchebnoe posobie [Workshop on Parallel Programming: A Study Guide]. Saint Petersburg: BVH, 2017, 236 pp., isbn: 978-5-9909805-0-1 (In Russian)
[49] Sanders J., Kandrot E. CUDA by Example: An Introduction to General-Purpose GPU Programming. London, Addison-Wesley Professional, 2010, 311 pp.
Information about the authors
Tverdyi Dmitrii AlexsandrovichA - PhD (Phys. & Math.), Researcher at the Electromagnetic Radiation Laboratory, Institute of Cosmophysical Research and Radio Wave Propagation, FEB RAS, Paratunka village, Russia, ©ORCID 0000-0001-6983-5258.
Parovik Roman IvanovichA - D. Sci. (Phys. & Math.), Associate Professor, Leading researcher laboratory of modeling physical processes Institute of Cosmophysical Research and Radio Wave Propagation FEB RAS, Paratunka village, Russia, ©ORCID 0000-0002-1576-1860.
