Численное решение нелинейной сингулярно возмущенной краевой задачи с пограничными и внутренними переходными слоями Текст научной статьи по специальности «Математика»

ФИЗИКО-МАТЕМАТИЧЕСКИЕ НАУКИ

Численное решение нелинейной сингулярно возмущенной краевой задачи с пограничными и внутренними переходными слоями Петухова Н. Ю.

Петухова Наталья Юрьевна /Petukhova Natalya Yurievna — кандидат физико-математических наук, доцент,

Центр математического образования, Московский политехнический университет, г. Москва

Аннотация: рассматривается сингулярно возмущенная нелинейная краевая задача для дифференциального уравнения. Уравнение исследовано в случае, когда в решении имеется внутренний переходный слой. Предложен численный метод, основанный на вычислении функций, составляющих асимптотическое разложение решения. Показана сходимость и устойчивость метода. Ключевые слова: сингулярно возмущенная задача, разложение по малому параметру, контрастная структура, внутренний переходный слой.

DOI: 10.20861/2312-8267-2016-30-001

1. Постановка задачи.

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

s2y"(x) = F(x, y(x)), 0 < x < 1

y(0) = A, y(1) = B. (1)

S > 0: малый параметр. Такая задача может моделировать распространение звука в среде с малой вязкостью или распространение тепла в тонком стержне.

Задача рассматривается в случае, когда F(x,y) нелинейно зависит от y и в ее решении y(x, s) имеется внутренний переходный слой - контрастная структура. Асимптотическое разложение решения подобной задачи по малому параметру получено в [1; 2]. Цель этой работы - получить приближенное решение (1) при малом е, построив численное приближение нескольких первых функций из асимптотического ряда.

Для наличия контрастной структуры типа ступеньки уравнение F(x,y)=0 должно иметь 3 корня; при этом на фазовой траектории эквивалентной (1) системы

Sz'Jx) = F(x, z ), 0 < x < 1

(2)

Sz2(x) = zi (Z2 = y(x) zi = y,(x))

должны быть два седла, соединенных сепаратрисами. Асимптотическое разложение задачи (1) в [2] построено при следующих условиях:

1) F (x, z2 ) бесконечно дифференцируема в некоторой области G.

2) Уравнение F(x, z2 ) = 0 имеет 3 корня z = Щ (x), i = 1, 2, 3; Щ (x) G G, 0 < x < 1

при этом Щ (x) < (p2 (x) < Щ (x).

3) Fz (x, щ (x)) > 0, i = 1,3; Fz (x,( (x)) < 0, 0 < x < 1. (3)

(2( x)

4) Существует решение уравнения I(x) = 0: x0 G (0; 1), где I(x) = J F(X, z)dz.

Щ( x)

(2(x0)

5) I'(x0) = J Fx (x0, z)dz * 0.

Щ(. xa~)

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

переходит из окрестности ф1 (х) в окрестность ф3 (х) . Это решение имеет следующую структуру:

у(хе) = !Г"~}(Х'£) + ^(х^ХО < X < х* ХХ'£) = |у„(+)(х,е) + ^е), х* < х < 1 (4) Для остаточного члена (х, е) справедлива оценка:

(х, е)| < Мпеп+1,0 < х < 1, Мп не зависит от е. В (4): х - внутренняя точка отрезка [0; 1], около которой происходит переход решения у(х, е) с функции ((х) на ((х); X заранее неизвестно и вычисляется в виде ряда: х* = х0 + еху + е2х2 +..., где х0 из (3); (5)

п п п

У"(х,е) = ин(х) + ПУ(т,е) + ОУн(т,е) = ^ЧЧх) + ^екПУкТ) + ХеО^Т). (6)

к=0 к=0 к=0

В (6): и( ^ (х) - регулярные функции асимптотического ряда, и0 (х) = ((х); ПУк ("Т0 ) -функции пограничного слоя в окрестности точки х = 0, Т0 = х / е, е > 0; ОУ(—Чт) -функции пограничного слоя в левой окрестности точки х , Т = (х — х ) / е.

Аналогично: У(+) (х, е) = и(+) (х) + ЯУ(т , е) + QY(+) (Т, е), и(х) - регулярные функции асимптотического ряда при х £ [х ; 1], и0(+)(х) =( (х); ЯУ (т ,е) - функции

пограничного слоя в окрестности точки х = 1, Т = (х — 1) / е; ОГ^+) (Т, е) - функции

*

пограничного слоя в правой окрестности точки х .

Теперь для получения приближенного решения задачи (1) с точностью Д^ = Мые^+1, N > 0, надо вычислить с такой точностью функцию

У^" (х,е) = -|

N N N

xеки[-)к(х) + xекПУк"Т) + xекОУ(-)к(т)., 0 < х < х*

к=0 к=0 к=0 (7) N N N .

ХекиГ (х) + ХекЯУ""Т) + ХекОУ(+)"(т), х* < х < 1

При этом для достижения данной точности Дд, каждую из функций приближенного решения (7) надо вычислять с точностью

бк =ДN /ек, к = 0,1,...,N. (8)

Построенная таким образом функция У" (х, е) будет равномерным по е на всем отрезке [0; 1] численным приближением решения задачи (1). При больших значениях N (N>3) реализация такого численного метода станет достаточно громоздкой; если же 0 < N < 3 , его применение оправдано.

2. Численное решение задач для регулярных и пограничных функций.

Далее рассмотрим подробно вычисление с требуемой точностью левого приближения решения у(х,е) - функции У^ —)) (х, е) ; функция " (х,е) вычисляется аналогично.

2.1. Вычисление регулярных функций.

Уравнения для регулярных функций И( ) (X) в (6) получаем, приравнивая коэффициенты при £к в разложении

да да . . к. (-У'^л т?/.. V к. (-) л „ .. „ ..» тг„„ „_п _______ СУ-V 1/(-)

) "(X) = Р(X,)(X)), 0 < X < X*. При к=0 имеем Р(X,х)) = 0,

к=0 к=0

согласно (3) выбираем и( ) (X) = (1 (X). Далее, раскладывая функцию Р(X, (( (X) + А)

да

(А = £кик (X)) в ряд по степеням е, получим уравнения для следующих регулярных функций:

к=1

к = 1: 0 = Р (X, ( (X)) • и{ ) (X). ^ и\ ) (X) = 0. (9) к = 2: ((X) = Р/ (X, (1 (X)) • и2-) (X) + >2Р/(X, ((X)) • (и" (X))2. о Р/ (X, (1 (X)) • и" (X) = ((X)

к = 3 : и> "(X) = Р/ (X,(> (X)) • и3-) (X) + Р;(x, (> (X)) • и>-) (X) • и" (X) + >6Р/3)(X, ((X)) • (и>(-) (X))3.

^ и (-)(X) = 0.

к=4: и2" (X) = Ру' (X, ( (X)) • и4-) (X) + >2 Р;( X, ( (X)) • (2 • и>-) (X) • и3-) (X) + (и2-) (X))2) +

>6 Ру(3) (X, (> (X)) • (3 • и>2 (X) • и2 (X) + и>3 (X)). о

и2 " (X) = Ру" (X, (X)) • и4-) (X) + >2 РД X, (X)) • (и2-) (X))2.

Уравнения (9) для функций и(-) (X) не являются дифференциальными; входящая в них функция

П , ч (-) / N. (-) П , ч

( (X) известна точно; при п>1 входящую в уравнение для Ы^ (X) функцию И^-г (X) можно

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

Так будет получена с заданной точностью регулярная часть асимптотического разложения решения

(4). Функции и(-) (X) вычисляем на отрезке 0 < X < X,), где X известно из (3).

2.2. Вычисление пограничных функций. А. Пограничные функции на левом конце отрезка.

Пусть X = X / £> 0, П^, ("Т0) = П^ (X). Уравнения для пограничных функций П^ (X) получаем, приравнивая коэффициенты при £к в уравнении

X екП/ (X) = Р(ей X ек (ик") (е() + Пк (X)) - РX ек (и<") (е()).

к =0 к=0 к=0

При £° имеем задачу для П0 (X) :

УО т ГП^ + И 0 < Г < оп

(10)

Ц^) = Р(0,т>(0) + П0^)), 0 < X < да

П0(0) = у(0)-т(0), П0(х) ^ 0.

Уравнения для следующих пограничных функций Пк (X), к > 1 будут линейными общего вида:

[Пк " (X) = Ру' (0, (1 (0) + П0 (X)) •Пк (X) + Ок (X), 0 < X < да Пк (0) = -«^(0), Пк (X) ^ 0.

^ г^+да

Здесь:

С,(/)=/• ср[(0)• (Р/(М)-Р'(М)) + /(Р'(М)-р'(М)); М = (0,й(0) + П0(/)),М = (0,я(0));

со=(у2(2- щ (0)+и2 (0)) • (р; (м) - г; (м»+у2и-щ т2 (р/ (м> - р/ (м» + (12) х /2 ср" (м) - р" (м))+е ■ (о)^; <м) - <му)+/ • ^ <м) ■ п, «+(м) • (о

(Уравнения написаны с учетом: Щ (х) = 0 ).

Уравнение для функции П0 (г) является нелинейным; его численное решение рассмотрим отдельно. Из условия (3) следует, что Р (0, ( (0) + П0 (г)) > к > 0 по крайней мере при г > ¿0 . Тогда, на основании метода дифференциальных неравенств [4], для функции П0 (г) справедлива оценка: |П0 (г)| < Се Ы . Для получения численного решения (10) с точностью 80 = Дм из (8), задачу для П0 (г) достаточно решать на отрезке 0 < г < Т, где Т определяется из условия:

е кТ° = 8,

^ = 1п(1 / 80 ) / к . При численном решении краевой задачи (10) сведем ее к задаче Коши для 2 (г) = (П0 (г), П0 (г)) . Для 2 (г) имеем следующую систему уравнений:

й2х / Жг = 22

, (13)

\й22 / Ж = Р(0,( (0) + 2), 0 < г < Т, 2(0) = (П0(0), П0 (0))

Легко найти первый интеграл (13):

2 • /Жг = Р(0,((0) + 2)• /Жг — >222(г) = |р(0,((0) + 2)• 2'(г)Жг —

ш 0

>2(222(^) — 222(0)) = |Р(0,(1(0) + 2>)• 2>(г)Жг = | Р(0,(>(0) + 21)Ж21 .

0 г (0)

0

Тогда (П0 (0))2 = —2 | Р(0, ( (0) + 2)Ж2 . (14)

П>(0)

Значение последнего интеграла можно без труда вычислить с требуемой точностью. Если для численного решения (13) применять какой-либо стандартный устойчивый метод порядка р, р > 1: многошаговый или метод типа Рунге-Кутта, то для его погрешности справедлива оценка [6]:

\х(гт) — II < 112(0) — 20" | • ещ>(гтЬ) + (р+1° (г„ )|| • (ехр^Ь) — 1) / ь,0 < гт < Т. (15)

Здесь 2" - численное решение в точке г = т" ; Ь - норма матрицы

8/ / 82 (/ = (22, Р(0, ( (0) + 2 )) . Норма матрицы А, согласованная с нормой

п

1211 = шах(|21) , вычисляется по формуле: Ь = шах(^^ |) . Для (13):

' ] =1

д/ / 82 =

Г0 Р'г (0,(1(0) + 2>)Л V1 0 у

и L = шах(1; шах

2£[ у (0),((0)]

F' (0, z) ) ограничена при

всех 1

Погрешность численного решения в (15) \2{гт ) — 2) || = Ж + Ж.

Ж =| 2(0) — 20" | ехр^Ь) <| 2(0) — 20" | ехр(Т Ь) < 12(0) — 2к0\\ • 8—Ь/к. Тогда <8^ если р(0) — 2" < 8>+Ь/к, т.е. значение (0) в (14) надо вычислять с точностью 8>

Ь/к ^0 .

й2 •||^+1)(^)||(ехр(^)-1) Для функции

г(р+1) ^ справедлива та же оценка, что и для П0 (/) ; (/) ~ в Ы, р > 1. Тогда: при

условии к>Ь: для выполнения требования < О0, шаг сетки выбирается только из условия необходимой точности: Нр < д{); если к<Ь и с12 ~ Нр • ехр((Х — , для выполнения

а2 < О0, шаг сетки надо подчинить условию: Ъ < О0 = О0 .

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

Краевые задачи для следующих пограничных функций Пк (X), к > 1 (11) являются линейными;

при этом для функций Ск (X) в (12) легко получить оценку \Ск (X) < С ехр(-Ы / 2) , такие же оценки на основании метода дифференциальных неравенств [4; 9] будут справедливы и для Пк (X), к > 1. Численное решение подобного типа задач подробно рассмотрено в [5]. Путем построения неравномерной сетки задачу (11) надо свести к решению краевой задачи на отрезке; далее для вычисления П ^ (X) можно использовать стандартную разностную схему второго порядка с шагом

сетки Нк , выбираемым только из условия требуемой точности (8): Ъ < 8к . Так же в [5] обоснована устойчивость такого численного алгоритма, позволяющая использовать значения ранее вычисленных пограничных и регулярных функций ПЪ (X), иЪ (X), 0 < / < к -1 при формировании правой части

и начальных условий для функции П ^ (X) .

Так будет построена с требуемой точностью пограничная составляющая численного решения (7) у( )Ъ (X, е) на левом конце отрезка при х=0.

Б. Пограничные функции слева от внутренней точки перехода. Здесь Т = (X - X*) / £, X* из (5), тогда X = X* + Т£ = X0 + ^, где

= £X^ + £ X2 + ... + £ X + ... + еТ . Уравнения для пограничных функций

^(т) = о- (т) в (6), получаем, приравнивая коэффициента: при £ в равенстве:

да да да

Хека" (Т) = Р (X + 4, Хек иЧX + 4) + О (Т)) - Р (Xo + 4, Хек (и«( Xo + 4)).

к =0 к =0 к=0

да

Здесь: р (^+а,, х е (и(-) (xo + 4)+а (т))=Р (X + ^1,(1 (X)+а, (т)+=р (м)+

к=0

+р;(м)а2 + рдмц + х(а•а/аx+а•а/Эу)2Р(м)+...; м = (Xo,ml(Xo)+&(т)),

а = u(_)'(x м +12 и(-) "(X м2 +...+е(и^ )+и^'^ м + ■■■)+£2(u(,-)(x0)+u(_)'(x м + ■■■)+ ...+£Й(Т)+£2е2(т)+... .

да

Аналогично: Р(X0 + ^, X £ (Мк ) (\ + ^ )) =

к=0

= Р(М) + Р'(М)«52 +ХМ -5/аг + ^2 •Э/ф)2Р(М) + ...;М = (х0,^1(х0)), ¿2 =

и,

(м + х и( м2 +...+е(и1< )+и ''(X м +...)+£(и( '(X)+и м +...)+...

Получим задачу для О (т) :

о; (т) = Р(х0, ( (х0 ) + о0 (т)), — ш < т < 0

о (0) = ( (х0 ) — ( (х0 ), О (—Ш) = Задачи для следующих пограничных функций снова будут линейными общего вида

О (т) = К (х0, ( (х0) + О0 (т)) • Ок (т) + Ик (т), — ш < т < 0

Ок (0) = Ч, Ок (т) — 0.

Значения определяются, приравнивая коэффициента: при еК в равенстве

ш

У )(х*) = Хек (ик )(х0 + ех1 + е2 х2 +...) + О (0)) = (2 (х*) = (2(х0 + ех1 + е2 х2 +...).

к=0

При к=1 имеем следующую задачу для О (т) :

'й'(т) = Р/ (х0,(( х0) + 00(т)) • О>(т) + И> (т), — ш<т< 0

(18)

О> (0) = х> (( (х0) — ((х0)), О> (т) — 0.

г —^—ш

И1(т) = (>'(х0)(х1 +т)(РДх0,((х0) + О0(т)) — Ру'(х0,((х0)) +

(х' + т)(р;(х0 , ((х0)+а (т))—Р;( х0, ((х0)).

Задача (18) содержит неизвестный параметр х . Чтобы найти его, надо приравнять коэффициенты при ек в равенстве У)(х*,е) = У(+\х*,е) — по (6):

ш ш ш ш

]Текик-)(х*) +ХекО«(0) = ]Хекик+)(х*) +ХекО«(0). (19)

к=0 к=0

учетом их(~) (х) = и(+) (х) = 0 ) имеем:

((х0) + О'(—)(0) = (з'(х0) + О'(+)(0) .

Подставляя в это равенство общий вид решения задачи (18) для о^ )(т) и аналогичной задачи для (т) , можно получить уравнение для х . В [2] для х1 получена следующая формула:

к=0 1

к =0

В (19) при е* (с

х =■

(з( х0) /

(1( х0)

¥'х(х0,2) | (2 | Р(х0, р)ф)—1/2 Л

(2 ( х0 ) (1 ( х0 )

Ж.2

(з(х0 )

(20)

| Рх(Xо, 2)Ж2

(' ( х0 )

Знаменатель этой дроби не равен нулю в силу (3) и для реализации нашего численного метода х можно определить из (20) с требуемой в (8) точностью 8 . При к>1 задача для О (т) содержит коэффициенты х.,1 < I < к для вычисления х надо приравнять коэффициенты при ек в (19).

Задача (16) для О (т) полностью аналогична задаче (10) для П0 (г) . При |т| >Т имеем оценки: Р'у (х0, ( (х0 ) + О (т)) > к > 0, |О0 (т)| < С ехр(—к | т|). Получить численное решение (16) с требуемой точностью 80 можно вышеописанным методом вычисления П0 (г) . Также при к=1: после определения х задача (18) аналогична задаче (11) для Пх (г) ; для функции

(г) справедлива оценка: (г)| < С ехр(—к Ц /2) , такая же оценка выполняется и для О (г) . Далее получить численное решение (18) можно по схеме из [5]. Этим же методом после определения X можно вычислить с требуемой точностью функции О (г) , к>1. Так будет

определена с требуемой точностью Ад, по формуле (7) функция У^ ^ (X, £) - численное

*

приближение решения задачи (1) слева от точки перехода 0 < X < X . Функцию Уг+)А (X, £), X < X < 1 можно получить аналогичным образом, взяв за начальное приближение й0 (х) = (ръ (х) . Обобщение результатов этого раздела приведем в следующей теореме.

Теорема 1. Функция У^ (X, £), полученная по формулам (7) - (8) является на всем отрезке 0 < X < 1 численным решением задачи (1), (3), удовлетворяющим условию: шах|у^,£) — У^ (X, £)| <МЫ£Н+1, N = 0; 1; 2. Константа Мы не зависит от е. 3. Вопросы практической реализации.

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

п0 (^), О) (г) при решении системы вида (13):

П0 ~ Т0 / к, где к = 8^1к , р > 1; Т0 = 1п(1 / С>0). Метод вычисления пограничных функций Пк ({), Ок (г), к > 1, описанный в [5], требует число действий Пк = 1 / Нк, Нк = ^8к .

Тогда общее число действий П = 0(П0 + Пк ) .

Традиционные численные методы решения сингулярно возмущенной краевой задачи на отрезке требуют при своей реализации условия Р' (X, у) > 0 , которое в нашем случае (условия (3)) не выполнено. При решении нелинейной задачи методы обычно имеют первый порядок сходимости О^), реже: 0(И\п(\ / ИУ)2 [3; 7; 8]. Далее, в случае решения нелинейной задачи эти методы приводят к системе нелинейных алгебраических уравнений большой размерности т ~ 1 / /7 = 1/ 8. Решать такую систему надо итерационным методом, где число итераций Т) ~ 1 / £ : велико при малом е.

Таким образом, здесь для нелинейной задачи общее число действий П = (1 / 8) • (1 / £) . В нашем случае реализация численного метода оказалась проще.

В заключении отметим еще одно полезное свойство вышеописанного численного алгоритма. Функции Ык (X), П (0, Ок (г) , составляющие численное решение У^ (X, £) в (7), не зависят от е.

Это позволяет, один раз вычислив У^ (X, £0), использовать эти же функции

Ык (X), П (), О (г) для получения численных решений серии задач (1), (3) при различных

значениях малого параметра £ Ф £0, легко проводить сравнительный анализ в моделях распространения звука или тепла.

Литература

1. Васильева А. Б., Бутузов В. Ф. Асимптотические методы в теории сингулярных возмущений. М.: Высшая Школа, 1990. 208 с.

2. Васильева А. Б., Бутузов В. Ф., Нефедов Н. Н. Контрастные структуры в сингулярно возмущенных задачах // Журнал вычисл. математики и матем. физики, 2001. Том 41. № 7. С. 799-851.

3. Дулан Э., Миллер Дж., Шилдерс У. Равномерные численные методы решения задач с пограничным слоем. М.: Мир, 1983. 200 с.

4. Нефедов Н. Н. Метод дифференциальных неравенств для некоторых сингулярно возмущенных задач в частных производных // Дифференц. уравнения, 1995. Том 31. № 4. С. 719-722.

5. Петухова Н. Ю. Численное решение краевой задачи для эллиптического уравнения с малым параметром на основе асимптотических представлений // Актуальные проблемы гуманитарных и естественных наук, 2016. № 6 (89). Часть 1. С. 15 - 21.

6. Хайер Э., Нерсетт С., Ваннер Г. Решение обыкновенных дифференциальных уравнений. Нежесткие задачи. М.: Мир, 1990. 512 с.

7. Шишкин Г. И. Сеточные аппроксимации сингулярно возмущенных эллиптических и параболических уравнений. Екатеринбург. УрО РАН, 1992. 215 с.

8. Шишкин Г. И. Обусловленность и устойчивость разностных схем на равномерных сетках для сингулярно возмущенного параболического уравнения конвекции диффузии // Журнал вычисл. математики и матем. физики, 2013. Том 53. № 4. С. 575-599.

9. Pao C. V. Nonlinear parabolic and elliptic equations. NY: Plenum Press, 1992. 807 p.

Энергия отделения протонов Sp и S2p в изотонических цепочках

с N=48, 49, 50, 51, 52 Кульков А. В.1, Казаков А. В.2

'Кульков Алексей Владимирович /Kulkov Alexej Vladimirovich — магистрант;

2Казаков Алексей Владимирович /Kazakov Alexej Vladimirovich — магистрант, кафедра прикладной математики, физико-математический факультет, Смоленский государственный университет, г. Смоленск

Аннотация: с использованием экспериментальных данных для дефекта масс ядер вычислены энергии отделения протонов Sp и S2p для семейства изотонических цепочек около магической с N=50. Для Sp явно проявляется чёт-нечёт осцилляция, для S2p таких осцилляций не наблюдается. Для семейства изотонических цепочек около магической с N=50 найдены значения спаривательной щели Ap. Ключевые слова: энергия отделения протонов, чет-нечет осцилляция, ширина спаривательной щели.

Энергии отделения нуклонов являются важными характеристиками внутренней структуры атомных ядер, сил взаимодействия между нуклонами в ядре. В данной работе вычисляются энергии отделения протонов Sp, S2p и ширина спаривательной щели Др для изотонических цепочек с N=48, 49, 50, 51, 52.

Энергии отделения протонов определяются известным образом [1, 3] через энергии связи соответствующих ядер B(Z,N):

SP = B(Z, N) - B(Z -1, N ) (1) S2P = B(Z, N) - B(Z - 2, N) (2)

Энергии связи ядер B(Z,N) можно найти, используя экспериментальные данные [4] для масс соответствующих нуклидов (M-A)ZN:

M(N,'Z) = Z • + N • m - B(N, Z), (3)

где M(N,Z) - масса нуклида с числом нейтронов N и числом протонов Z, MH - масса нуклида водорода, mn - масса нейтрона. Все величины выражаются в эВ. Данные также можно переписать через соответствующие дефекты масс.

Ширина спаривательной щели вычисляется [5] по формуле

Л = 1 [B(Z, N) - 2B(Z -1, N) + B(Z - 2, N)]. (4)

Численные значения Sp (в кЭв) приведены в таблице 1 «Энергии отделения протонов Sp».