CC BY
70
МАТЕМАТИЧЕСКАЯ ФИЗИКА
УДК 519.63
ИТЕРАЦИОННОЕ РЕШЕНИЕ ОДНОЙ НЕКЛАССИЧЕСКОЙ ЗАДАЧИ ДЛЯ УРАВНЕНИЯ КОЛЕБАНИЙ СТРУНЫ
В.И. Васильев, В.В. Попов, М.С. Еремеева, А.М. Кардашевский
Северо-Восточный федеральный университет им. М.К. Аммосова, Якутск, Российская Федерация e-mail: [email protected]
Рассмотрена неклассическая задача для гиперболичеекого уравнения второго порядка, в которой кроме граничных условий на концах струны заданы дополнительные условия: в начальный момент времени — скорость движения струны, в конечный момент времени — смещение струны. Для численного решения поставленной задачи предложен итерационный метод, ранее использованный авторами настоящей статьи для решения ретроспективной задачи теплопроводности. Приведены примеры расчетов для модельных задач, в том числе со случайными погрешностями во входных данных.
Ключевые слова: обратные задачи, неклассическая задача, конечно-разностный метод, метод сопряженных градиентов, случайные погрешности.
ITERATIVE SOLUTION OF A NONCLASSICAL PROBLEM FOR THE EQUATION OF STRING VIBRATIONS
V.I. Vasilyev, V.V. Popov, M.S. Eremeeva, A.M. Kardashevskiy
North-Eastern Federal University n.a. M.K. Ammosov, Yakutsk, Russian Federation e-mail: [email protected]
The article considers a nonclassical problem of the second-order hyperbolic equation, where some additional conditions are given as a supplementary to the boundary conditions at the ends of the strings: at the initial time it is the speed of the string, at the final time — displacement of the string. The iterative method is proposed for the numerical solution of the problem. This method was previously used by the authors of this article to address the retrospective thermal conductivity problem. Sample calculations of model problems including the ones with random errors in the input data are given.
Keywords: inverse problems, nonclassical problem, finite difference method, conjugate gradient method, random errors.
Введение. Неклассическая задача с заданной скоростью точек струны в начальный момент времени и смещением струны на конечный момент времени для гиперболического уравнения, как и задача Дирихле, относится к классу условно корректных задач математической физики [1, 2]. Общий подход к решению таких некорректных задач базируется на использовании градиентных итерационных методов. Основу градиентных итерационных методов составляет нахождение приближенного численного решения некорректных задач из
итерационной минимизации соответствующего функционала. В цикле работ С.И. Кабанихина и его учеников, например в работе [3], численно исследовано решение задачи Дирихле для волнового уравнения методами итераций Ландвебера и наискорейшего спуска, проведено теоретическое исследование устойчивости и сходимости приближенного решения к решению исходной дифференциальной задачи для одномерного и двумерного уравнений гиперболического типа второго порядка.
Для восстановления начального условия задачи теплопроводности в работе [4] предложен быстро сходящийся итерационный метод минимальных невязок, продемонстрирована его высокая эффективность при решении модельных задач. В работе [5] для численного решения задачи Дирихле для гиперболического уравнения второго порядка использован итерационный метод сопряженных градиентов. В настоящей работе для численного решения неклассической задачи для уравнения колебаний струны применен аналогичный подход с итерационным уточнением начального условия. Для поставленной неклассической задачи по заданному решению на конечный момент времени итерационным методом сопряженных градиентов уточнено значение решения в начальный момент времени. Результаты вычислительного эксперимента демонстрируют высокую эффективность предлагаемого итерационного метода.
Постановка неклассической задачи. В прямоугольнике (0, /) х х (0,Т) ищем функцию и(х,г) — решение одномерного гиперболического уравнения второго порядка
д2и д , ч ди\ , т
^-0<х<;, 0<г<Т (1)
Предположим, что концы колеблющейся струны х = 0, х = I закреплены, т.е. выполняются однородные граничные условия первого рода:
и(0, г) = 0, и(1, г) = 0, 0 < г < Т. (2)
В начальный момент времени задана скорость точек струны, а в конечный момент времени — смещение струны
ди
— (х, 0) = V(х), и(х, Т) = р(х), 0 < х < 1. (3)
Пусть коэффициент к(х) дифференциального уравнения (1) является достаточно гладкой, т.е. имеет непрерывные производные до третьего порядка включительно, положительной, ограниченной функцией и удовлетворяет условиям 0 < к\ < к(х) < к2 < то. Рассматриваемая неклассическая начально-краевая задача относится к классу
условно корректных задач математической физики [6]. В соответствии с приведенным ниже простейшим примером эта задача может иметь неединственное решение. При задании v(x) = 0, p(x) = 0, k(x) = 1 задача (1)-(3) имеет семейство решений, описываемых формулой u(x,t) = csin(n1nx) cos(n1nt), l = nt, T = n3 — (2n4 — 1) /(2n1), c G R, n1,nt,n3,n4 G N.
Разностная задача. Обозначим whT множество внутренних узлов пространственно-временной сетки whT = wh х wT, где wh = {x¿ = ih, i = 1, 2,... ,m — 1}; wT = {tj = jr, j = 1, 2,... ,n — 1}; wh = {x¿ = ih, i = 0,1,..., m}. На множестве сеточных функций y G H таких, что y(x) = 0, x G wh, определим сеточный оператор A соотношением
АУ = —(a(x)yx)x =
1( ( , h\y(x +h) — y(x) í Mx) — y(x — h)\
= —— a(x + h)-;--a(x)-;- , x G wh
h h h
приняв, например, a(x) = k(x — 0, 5h).
В конечномерном сеточном гильбертовом пространстве H скалярное произведение и норму зададим соотношениями (y, v) = ^^ yvh,
11у11 = (у,у)1/2. В гильбертовом пространстве H оператор A является самосопряженным, положительно определенным (A = A* > 0) при 0 < k1 < k(x) < kt < ж и имеет место двухстороннее неравенство 8ki < ||A|| < 4k2/h2.
Сначала от исходной неклассической задачи (1)-(3), применив дискретизацию по пространственной переменной, перейдем к задаче определения решения дифференциально-операторного уравнения
djty + Ay = 0, x G wh, 0 <t<T (4)
при заданных дополнительных условиях
dy;(0) = v, y(T) = <p, x G Wh. (5)
При использовании симметричной трехслойной разностной схемы с весовым множителем а дискретный аналог задачи (4), (5) имеет вид
y — 2y + У + A(ay + (1 — 2а)у + ау) = 0, t G wT; (6)
т
2
у(т) = у(0) + ти, у (и) = X Е (7)
где 0 < а < 1, применена безиндексная система обозначений, введенная А.А. Самарским [7]:
у = у3+1, у = у3, у = у3-1, X Е Шн, 3 = 1, 2,... ,п - 1.
Для численной реализации системы линейных алгебраических уравнений (6), (7) воспользуемся итерационным методом сопряженных градиентов, позволяющим уточнить смещение струны в начальный момент времени. На каждой итерации решается прямая задача для уравнения колебаний струны с помощью стандартной симметричной трехслойной разностной схемы, обладающей вторым порядком аппроксимации по к и т [7]. Итак, пусть вместо обратной задачи рассматривается прямая задача для этого же уравнения, когда второе условие в (5) заменено начальным условием у(х, 0) = <^(х), х € шь,. Дискретный аналог прямой задачи имеет вид
у - 2у + у + А(ау + (1 - 2а)у + ау) = 0, г € шт; (8)
т
2
y(0) = р, У (т) = y(0) + TV, x G Wh. (9)
* >7 " (10)
Разностная схема (8), (9) устойчива при выполнении условия [7]
1 1
4 - т
и справедлива следующая априорная оценка:
у(г) = ||у(г + т)||,<||У(г)1< ... < 11^(т)||„, г € , (11)
где
У
норма,
У
2 1 1
= 4 НУ + yllA+т2 bllLa/^ R = 4E+
Тем самым норма решения задачи Коши со временем убывает. Отметим, что для явной разностной схемы (а = 0) условие устойчивости (10) принимает вид т
Итерационный метод. Для нахождения решения сеточной задачи (6), (7) используем наиболее быстро сходящийся итерационный метод вариационного типа — метод сопряженных градиентов [8], основанный на последовательном уточнении искомого начального условия V(х) с дальнейшим решением на каждой итерации дискретного аналога прямой задачи (8), (9). Придадим этой задаче соответствующую операторную формулировку. Последовательно исключая промежуточные значения сеточной функции у (г), г € шт при заданных значениях и V на конечный момент времени, получаем у(гп) = Л^ + Bv, где А, В — операторные полиномы п-й и (п — 1)-й степеней от положительно определенного и самосопряженного оператора А. С учетом этого приближенному решению обратной задачи естественно сопоставить решение следующего операторного уравнения: А^ = <^(х) — Bv = В силу самосопряженности оператора А самосопряженным является и операторный полином Л.
Для решения операторного уравнения используем трехслойный итерационный метод сопряженных градиентов, записанный в каноническом виде [8]
^/к+1 = + (1 — 0+1)^-1 — а^вк+Л, к = 0,1,..., (12)
где ;к — невязка, ;к = — Для вычисления итерационных параметров вк+1, ак+1 используются следующие рекуррентные формулы
вк+1=(А;г;)), к=0,1,...; (13)
Л в+1 (;к ,;к) М-1 7 10 1 Л/1Л
а^+1 = 1--—--- — , к = 1, 2,..., «1 = 1. (14)
\ в/ (;к-1,;к-1) а/
В этом итерационном методе вычисление вектора уП осуществляется численной реализацией прямой задачи
— 2ук + + Л(а^ + (1 — 2а)ук + ау/) = 0, г € Шт, (15) т2
с начальными условиями
У0 = ^, У1 = ^ + тV, х € (16)
Следовательно, невязка определяется по формуле ;к = уП — <£, х € Ш^. Аналогично выполняется вычисление вектора = А;к :
¿/к — +
т2
+ A(ay + (1 - 2a)zfe + ) = 0, (x,t) G , (17)
с начальными условиями
20 = , = , х € (18)
Таким образом, в случае неклассической задачи для уравнения колебаний струны итерационный метод сопряженных градиентов реализуется в следующем порядке.
1. Полагаем к = 0 и задаем начальное приближение искомой функции ^, х €
2. Запускаем счетчик итераций к = к + 1, последовательно решая прямые задачи (15), (16) и (17), (18), определяем невязку ;к = у^ — <£, х € ш^, и вектор = А;к, х € ш^.
3. По рекуррентным формулам (13), (14) находим значения итерационных параметров вк+1, ак+1.
4. По формуле (12) рассчитываем очередное приближение искомого начального условия ^(х)к+1, х € Ш^.
5. Процесс продолжается до тех пор, пока не выполнится критерий выхода из итерационного цикла ||;к || < е, иначе продолжаем процесс, возвращаясь к пункту 2.
Примеры расчетов. Вычислительную эффективность предложенного метода удобно проиллюстрировать на примере численного решения простейшей неклассической задачи для уравнения колебания струны. Будем искать приближенное решение задачи для уравнения
д2и д2и ^
— — — = 0, 0 < х < 1, 0 <г<Т, дЬ2 дх2
с однородными граничными условиями и(0,Ь) = 0, и(1,Ь) = 0, 0 < Ь < Т, и дополнительными условиями в начальный и конечный моменты времени
ди
— (х, 0) = 0, и(х, Т) = ф(х), 0 < х < 1.
В рамках квазиреального вычислительного эксперимента ограничимся примером численного решения обратной задачи, которая соответствует решению классической прямой задачи для уравнения колебаний с теми же граничными условиями и двумя условиями в начальный момент времени:
и(х, 0) = р(х), ди(х 0) =0, 0 < х < 1.
Из решения прямой задачи определяем функцию р(х), которая и присутствует в постановке обратной задачи: и(х,Т) = ц>(х), 0 < х < 1. Приведем результаты решения обратной задачи в условиях, когда искомая функция ц>(х) имеет разный вид для некоторых значений Т. В первом случае зададим искомую функцию в виде финитной функции (р(х) = е-50(х-1/2) , 0 < х < 1. В расчетах принято т = 100, е = 0,001.
В качестве начального приближения искомого начального условия возьмем функцию, отдаленно напоминающую искомую функцию = э1п(пх), 0 < х < 1. График (кривая 1) искомого начального условия у = ц>(х), 0 < х < I, найденного предложенным итерационным методом, приведен на рис. 1, а, график заданного дополнительного условия в конечный момент времени (кривая 2) при п = 80, Т = 0,8 — на рис. 1, б. Число итераций 11, время расчета на персональном компьютере с процессором 17 0,7 с. Аналогичные графики представлены на рис. 1, в при п = 180, Т = 1, 8, в этом случае число итераций осталось неизменным (11), а время расчета составило 1,55 с.
Результаты расчета, когда искомая функция задана в виде более сложной функции ф) = —е-100(ж-г/4)2 + е-20(х-1/2)2 - е-100(ж-Зг/4)2,
0 < х < 1, приведены на рис.1, б и г. В расчетах также принято т = 100, е = 0,001. В этом случае итерации сходятся примерно в 2,5 раза медленнее, при п = 80, Т = 0,8 (см. рис. 1, б) и п = 180, Т = 1, 8
Рис. 1. Решения неклассической задачи в условиях различного задания функции (р(х)
(рис. 1, г) соответственно, число итераций 26, 27 и время расчета 1,66 и 3,83 с.
Зависимость к(— ^ е) на последовательности сгущающихся сеток приведена на рис. 2, а. При сгущении сетки наблюдается увеличение числа итераций для всех значений е. С уменьшением значения е число итераций возрастает.
Зависимость нормы невязки от числа итераций на последовательности сгущающихся сеток приведена на рис. 2, б. Ордината показывает критерий выхода из итерационного цикла, а абсцисса — число итераций. При сгущении сетки наблюдается увеличение числа итераций. С уменьшением значения е число итераций существенно возрастает.
Рис.2. Зависимости к(—1§е) (а) и е(к) (б) на последовательности сгущающихся сеток для щ = 50, п2 = 40 (1), п± = 100, п2 = 80 (2), щ = 200, п2 = 160 (5)
Следует отметить, что кривые показывают немонотонную зависимость нормы невязки от числа итераций.
Приведем характерные результаты решения неклассической задачи в условиях, когда функция ^(х) задана с погрешностью. В экспериментах функция ^(х), х Е и, возмущалась следующим образом: (х) = ^+8а(х), х Е и, где а(х) — случайные величины, равномерно распределенные на отрезке [-1, 1]. Итерационный процесс сопряженных градиентов обрывался по достижению нормы невязки значения 8, т.е. при ||гщ) || < 8. Полученное приближенное решение (кривая 1) для уровня погрешности во входных данных, определяемых величиной 8 = 0,05 (20%), и исходная функция ^(х) (кривая 2) приведены на рис. 3, а. Найденная сеточная функция ^ (искомое начальное условие) и точное начальное условие (кривая 2) представлены на рис. 3, в. Согласно этим зависимостям, искомое начальное смещение струны найдено с сохранением начального уровня возмущений.
Для выделения более гладкого решения, принадлежащего вместо пространства Ь2 (и) сеточному пространству W\ (и), зададим сглаживающий оператор Б^ = —+ х Е и, в виде сеточных функций, обращающихся в нуль в граничных узлах.
Возмущенный вариант входящей функции ^ (х) (кривая 1) и функция ^(х) (кривая 2), которая определена в результате сглаживания функции (х), приведены на рис. 3, б. Искомое точное решение (кривая 2) и решение полученное с помощью предлагаемого итерационного процесса (кривая 1) представлены на рис. 3, г. Здесь даны результаты расчета на пространственной сетке М = 100, Н = 0,01 с помощью
в г
Рис. 3. Характерные результаты решения неклассической задачи в условиях, когда функция (ж) задана с погрешностью
явной схемы для волнового уравнения а = 0 с временным шагом, удовлетворяющим условию устойчивости т = 0,01 < h, T = 0,8, N = 80. Поскольку каждый раз генератор случайных чисел выдает новый набор случайных чисел итерационный процесс сходится за 1030 итераций. В расчетах принято е = 0,05. Очевидно, рассчитывать на точное восстановление начальной функции нельзя в силу пониженной гладкости входных данных. Приближенное решение исходной задачи с хорошей точностью удается получить при задании более гладкого дополнительного условия в конечный момент времени, т.е. при сглаживании входных данных. В двумерном случае число итераций значительно возрастает. Тогда для существенного ускорения скорости сходимости итерационного процесса возникает необходимость использования предобусловливателя.
Заключение. Приведенные результаты вычислительного эксперимента подтверждают работоспособность предлагаемого итерационного метода решения неклассической задачи для гиперболического уравнения, а также его достаточно высокую эффективность.
Авторы благодарны профессору П.Н. Вабищевичу за конструктивные замечания и советы.
Работа выполнена при финансовой поддержке РФФИ и РФФИ Дальний Восток (гранты № 13-01-00719а, № 12-01-98514).
ЛИТЕРАТУРА
1. Samarskii A.A., Vabishchevich .PN.Numerical Methods for Solving Inverse Problems of Mathematical Physics. De Gruyter, 2007.
2. Kabanikhin S.I. Inverse and ill-posed problems: theory and applications. De Gruyter, 2012.
3. An optimization method in the Dirichlet problems for the wave equation / S.I. Kabanikhin, M.A. Bektemesov, D.B. Nurseitov, O.I. Krivorotko, A.N. Alimova // Journal of inverse and ill-posed problems. 2012. No. 2 (20). P. 193-211.
4. Самарский А.А., Вабищевич П.Н., Васильев В.И.Итерационное решение ретроспективной обратной задачи теплопроводности // Матем. моделирование. 1997. Т. 9. № 5. C. 119-127.
5. Вабищевич П.Н., Васильев В.И.Итерационное решение задачи Дирихле для гиперболического уравнения // Труды Х Межд. конф. "Сеточные методы для краевых задач и приложения - 2014". Казань: Изд-во КФУ, 2014. С. 162-166.
6. Тихонов А.Н., Самарский А.А. Уравнения математической физики. М.: Наука, 1972.
7. Самарский А.А. Теория разностных схем. М.: Наука, 1989.
8. Самарский А.А., Николаев Е.С. Методы решения сеточных уравнений. М.: Наука, 1978.
REFERENCES
[1] Samarskii A.A., Vabishchevich P.N. Numerical Methods for Solving Inverse Problems of Mathematical Physics. De Gruyter, 2007.
[2] Kabanikhin S.I. Inverse and ill-posed problems: theory and applications. De Gruyter, 2012.
[3] Kabanikhin S.I., Bektemesov M.A., Nurseitov D.B., Krivorotko O.I., Alimova A.N. An optimization method in the Dirichlet problems for the wave equation. J. Inverse and ill-posed Problems, 2012, no. 2 (20), pp. 193-211.
[4] Samarskii A.A., Vabishchevich P.N., Vasil'ev V.I. Iterative Solution of a Retrospective Inverse Problem of Heat Conduction. Matem. Modelirovanie [Mathematical Models and Computer Simulations], 1997, vol. 9, no. 5. pp. 119127 (in Russ.).
[5] Vabishchevich P.N., Vasil'ev V.I. Iterative Solution of the Dirichlet Problem for a Hyperbolic Equation. Proc. X Int. Conf. "Mesh Methods for Boundary-Value Problems and Applications-2014". Kazan, KFU Publ., 2014, pp. 162-166 (in Russ.).
[6] Tikhonov A.N., Samarskii A.A. Equations of mathematical physics. Courier Dover Publications, 1990.
[7] Samarskii A.A. The theory of difference schemes. Marcel Dekker, 2001.
[8] Samarskii A.A., Nikolaev E.C. Numerical Methods for Grid Equations. Springer, 1989.
Статья поступила в редакцию 24.12.2014 Васильев Василий Иванович — д-р физ.-мат. наук, профессор, заведующий кафедрой "Вычислительные технологии" СВФУ им. М.К. Аммосова. Автор более 110 работ в области вычислительной математики и математического моделирования. СВФУ им. М.К. Аммосова, Российская Федерация, 677000, Якутск, ул. Белинского, д. 58.
Vasilyev V.I. — Dr. Sci. (Phys.-Math.), professor of "Computational Technology", head of "Computational Technologies" department of the North-Eastern Federal University n.a. M.K. Ammosov (NEFU). Author of more than 110 publications in the field of computational mathematics and mathematical simulation.
North-Eastern Federal University n.a. M.K. Ammosov, ul. Belinskogo 58, Yakutsk, 677000 Russian Federation.
Попов Василий Васильевич — канд. физ.-мат. наук, доцент кафедры "Вычислительные технологии" СВФУ им. М.К. Аммосова. Автор более 40 работ в области вычислительной математики и математического моделирования.
СВФУ им. М.К. Аммосова, Российская Федерация, 677000, Якутск, ул. Белинского, д. 58.
Popov V.V. — Cand. Sci. (Phys.-Math.), assoc. professor of "Computational Technology" department of the North-Eastern Federal University n.a. M.K. Ammosov (NEFU). Author of more than 40 publications in the field of computational mathematics and mathematical simulation.
North-Eastern Federal University n.a. M.K. Ammosov, ul. Belinskogo 58, Yakutsk, 677000 Russian Federation.
Еремеева Майя Семеновна — аспирантка кафедры "Вычислительные технологии" СВФУ им. М.К. Аммосова.
СВФУ им. М.К. Аммосова, Российская Федерация, 677000, Якутск, ул. Белинского, д. 58.
Yeremeeva M.S. — post-graduate of "Computational Technology" department of the NorthEastern Federal University n.a. M.K. Ammosov (NEFU).
North-Eastern Federal University n.a. M.K. Ammosov, ul. Belinskogo 58, Yakutsk, 677000 Russian Federation.
Кардашевский Анатолий Михайлович — заведующий кабинетом кафедры "Вычислительные технологии" СВФУ им. М.К. Аммосова.
СВФУ им. М.К. Аммосова, Российская Федерация, 677000, Якутск, ул. Белинского, д. 58.
Kardashevskiy A.M. — head of the cabinet "Computational Technology" department of the North-Eastern Federal University n.a. M.K. Ammosov (NEFU). North-Eastern Federal University n.a. M.K. Ammosov, ul. Belinskogo 58, Yakutsk, 677000 Russian Federation.
Просьба ссылаться на эту статью следующим образом:
Васильев В.И., Попов В.В., Еремеева М.С., Кардашевский А.М. Итерационное решение одной неклассической задачи для уравнения колебаний струны // Вестник МГТУ им. Н.Э. Баумана. Сер. Естественные науки. 2015. № 3. C. 77-87.
Please cite this article in English as:
Vasilyev V.I., Popov V.V., Eremeeva M.S., Kardashevskiy A.M. Iterative solution of a nonclassical problem for the equation of string vibrations. Vestn. Mosk. Gos. Tekh. Univ. im. N.E. Baumana, Estestv. Nauki [Herald of the Bauman Moscow State Tech. Univ., Nat. Sci.], 2015, no. 3, pp. 77-87.
