Труды Петрозаводского государственного университета
Серия “Математика” Выпуск 17, 2010
УДК 517.954
И. А. Чернов
СХОДИМОСТЬ СЕТОЧНО-ИНТЕРПОЛЯЦИОННЫХ АППРОКСИМАЦИЙ РЕШЕНИЯ КВАЗИЛИНЕЙНОЙ ПАРАБОЛИЧЕСКОЙ КРАЕВОЙ ЗАДАЧИ НА ОТРЕЗКЕ
В статье рассматривается квазилинейная параболическая краевая задача III рода на отрезке: коэффициенты уравнения и правые части граничных условий нелинейно зависят от времени, точки и предыстории решения. Доказана сходимость разностной схемы к обобщенному решению задачи.
Постановка задачи. Рассматриваем следующую краевую задачу для параболического уравнения в частных производных в прямоугольнике П = [0, Т] х [0, Ь]:
54е = А (г, х, с(•, + В(4, х, е(-, •))5Жс- (1)
- В(4, х, с(•, •))5жс - £(£, х, с(4, х), с(^, •)), (2)
джс(£, х) = ( — 1)ж/ьС(£, х, с(4, х), с(^, •)), х = 0, Ь, 4 € [0,Т], (3)
с(0,х) = ^(ж) € С2([0, Ь]), 0 < ^ < 1, х € [0, Ь]. (4)
Коэффициент А = А(£, х, с(^, •)) является непрерывным функционалом на метрическом пространстве 4хС = П х С(П) (здесь и далее С — банахово пространство непрерывных функций с нормой У • У) со следующими свойствами:
1) 0 < А- < А < 1 (т. е. множество значений функционала отделено от нуля и ограничено);
2) А(С,х,с(-, •))= АС, х, м(-, •)) для любых функций с(4, х) и и(4, х) из С(П) при условии с(4, х) = и(4, х) при 4 < £ (независимость от будущего);
© И. А. Чернов, 2010
3) частные производные по 4 и х существуют на 4хС и тоже являются непрерывными на £хС функционалами с ограниченным множеством значений.
Коэффициенты В и В также являются непрерывными функционалами на £хС со свойствами, аналогичными А, но имеют ряд отличий:
1) нижняя граница множества значений может быть нулевой;
2) производная функционала по 4 является функционалом, непрерывным лишь на пространстве £хС1,0 = П х С1,0 (П) (пространство С1,0 состоит из непрерывных функций с непрерывной частной производной по 4 с нормой УсУ 1,0 = ||с|| + ||д4с||);
3) |д(В(4,х,с)| < Кцд) + Й2(.в)||д(с|| на 4хС для констант К*, не зависящих от (4, х) и с;
4) дХВ < 0, дХВ < 0;
5) 3ХВ < 0 не зависит от х.
Функционал Е = Е(4,х, с, «(•, •)) определен и непрерывен на пространстве 4хДС = 4хС х Д1 и имеет частные производные по первым трем (конечномерным) аргументам, причем
1) дХЕ непрерывна на 4хДС;
2) д4Е непрерывна на 4хДС1,0 = П х Д х С 1,0(П);
3) дсЕ непрерывна на 4хДС и множество ее значений положительно, отделено от нуля и ограничено: 0 < Кз(£) < дсЕ < Й4(£) (либо Е = 0);
4) |д(Е(4, х, с)| < К1(Е) + Й2(Е)||д(с|| на 4хДС для констант К4(Е), не зависящих от (4, х) и с.
Без ограничения общности можно считать, что д4А + К*2(в) + -^(д) + К*2(Е) + К4(Е) < 0 (этого всегда можно достичь заменой переменных 4 и х, сохраняющей остальные свойства коэффициентов).
Граничные условия (3) — нелинейные III рода. В правой части получается знак «плюс» при х = 0 и «минус» — при х = Ь. Правые части С (4, х, с, м(-, •)) и их производные и дсС при х = 0, Ь являются непрерывными функционалами на пространстве £ДС = [0, Т] х ДхС(П). При этом множества значений С, и дсС ограничены; зна-
чения функции м(£, х) при 4 > С не существенны; С(4, х, с, м(-, •)) > 0 при достаточно больших с; С(4, х, с, м(-, •) < 0 при с < 0 и х = 0, Ь;
дсС(£, х, с, «(•, •) > 0.
К таким задачам сводятся многие модели взаимодействия водорода с металлами, включая модели формирования или распада гидридов [1, 2, 3]. При этом с(£, х) — концентрация водорода в образце. Ограничения на коэффициенты оправданы физическим смыслом.
Примеры функционалов указанного вида доставляют подходящие ограниченные гладкие функции с ограниченными частными производными от 4, х и с(£, £) при фиксированном £ или от интеграла от с(£, х) (или от гладкой функции от с) по х € [0, Ь] (физический смысл
— количество водорода в образце в данный момент времени), возможна зависимость от значения с в определенной точке (например, от значения концентрации на краях), запаздывание (зависимость от с(£ — п, £)) и так далее.
Функционалы типа коэффициента А должны допускать дифференцирование по 4, что более ограничительно: примерами являются функции от 4, х и вектор-функции п(4), где ^п = Д(4, с(-, •)) Л и компоненты Д — функционалы типа коэффициента В. Смысл п — внутреннее состояние системы.
Ограничимся двумя содержательными примерами компонент п. Это температура образца, меняющаяся вследствие тепловыделения поверхностных процессов в зависимости от теплоемкости образца, зависящей от количества содержащегося в нем водорода, в связи с чем в уравнение для температуры входит интеграл от с(£, х) по [0, Ь] [4]. Также при выпрямлении подвижного фронта в задачах с фазовым переходом положение свободной границы — компонента п, причем в уравнении возникает слагаемое с коэффициентом В, зависящим от с(£, 0). Слагаемое типа В возникает, например, при записи уравнения в криволинейных координатах [2].
Граничные условия указанного вида (квазилинейные III рода) описывают нелинейные процессы на границах. Поскольку нормальная производная на границе имеет смысл плотности диффузионного потока, то правая часть О описывает поток. Зависимость от времени связана с переменными условиями, а от предыстории процесса — с теми же причинами, что и зависимость коэффициентов уравнения от предыстории.
Указанные выше примеры применимы и для граничных условий: поток в общем случае зависит от температуры, а выпрямление подвижного фронта приводит к появлению свободной границы в О. Кроме того, О зависит и от значения решения в текущей точке, что имеет физический смысл. Необходимое для полученных ниже результа-
тов свойство монотонности правых частей граничных условий также вполне естественно: повышение концентрации диффузанта (или теплоты) понижает поток.
Следует отметить, что во многих моделях решение предполагается положительным; примеры доставляют абсолютная температура или концентрация. В этих случаях зависимости плотности потока от решения могут описываться функциями, не обладающими монотонностью на Д. Примерами являются закон излучения Стефана — Больцмана &Т4 (& — постоянная Больцмана, Т — абсолютная температура) и квадратичная десорбция 6с2 (6 — константа десорбции, с — концентрация). Применению излагаемых в настоящей работе результатов это не противоречит, поскольку формально можно переопределить правые части для отрицательных значений (соответственно положив, например, sgn(T)^Т4 и sgn(c)6c2). Полученные функции совпадают на положительной полуоси с исходными, также непрерывны вместе с первой производной и монотонны. Если окажется, что решение действительно положительно (а так оно и есть при условии, что начальные значения и граничные условия согласуются с физическим смыслом), то значения потоков при отрицательных аргументах не играют роли, и решение удовлетворяет исходной задаче.
Нелинейный член нулевого порядка Е в уравнении в частных производных обобщает линейный случай; подобные линейные члены возникают в моделях при учете объемных процессов, зависящих от концентрации. Примером является формирование высших гидридов. Нелинейность этого члена незначительно усложняет результаты, полученные для линейного случая, но существенно меняет расчетный алгоритм, поскольку все уравнения разностной схемы (а не только соответствующие граничным условиям) становятся нелинейными. Тем не менее система уравнений имеет особую структуру, допускающую построение эффективного расчетного алгоритма.
Опишем кратко используемый в статье подход к поставленной задаче. Построим неявную разностную схему и докажем равномерную ограниченность сеточного решения и его сеточных производных, применив аналог принципа максимума. Это влечет предкомпактность непрерывных аппроксимаций — интерполяций сеточного решения в пространстве непрерывных функций и слабую предкомпактность в соболевском пространстве функций с первыми обобщенными производными. Предел удовлетворяет некоторому интегральному тождеству, то есть является обобщенным решением задачи. Метод заимствован
из [5]. Констатируем существование обобщенного решения и сходимость к нему разностной схемы.
Разностные аппроксимации. Введем в П равномерную сетку БN с шагами Н и т, 0 < п < N, 0 < * < 1. Узлы этой сетки (п, *)
— это точки (£„,х*), где = пт, х* = *Н. Совокупность узлов при одном п назовем слоем. Введем обозначения с(£„,х*) = с^ и аналогичные. Пусть ст задано в _____1; построим в [0,1] х [0, Ь] непре-
рывную функцию с„_1(£,х) линейной интерполяцией точек с^. Коэффициенты уравнения (2) на слое п аппроксимируем явно: ВП = В(£п-1, х*,2„_1(^, •)) и так же для А и В. Коэффициент Е и правые части (3) — явно-неявным образом: ЕП = Е(4П_1, х*, с^, с„_1(-, •)), ОП = О(4„,х*,сП,с„_1(^, •)) при * = 0,1.
Исключая граничные узлы, получаем множество ОN узлов (п, *) при 1 < п < N, 1 < * < 1 — 1. Границу обозначим N = Б N \ ОN. Введем обозначения:
с__с с^1___с
дтсП = п тп_1, 5ьсП = п Н п, д2сП = д^С1).
В дальнейшем будут использоваться формулы сеточного дифференцирования произведения
дь(фпфп) = дь(ФП)фП + дь(фп)ф дт (ВД) = дт (ФП)ФП + дт (ФП )*П
І+1
п—1 •
Обозначим через Б N и Б N подмножества Б N, на которых определены сеточные производные д^еП и дт сП- Введем также DN и аналогично DN - Заменяя производные сеточными аналогами, получим неявную схему с погрешностью аппроксимации 0(Л. + т):
дт сП = АПд2сП + ВПдьсП - ВП дьС—1 - ЕП = ^П^П, (5)
д^сП—1 = -^П д^сП = С со = ^(х)• (6)
Принцип максимума для сеточной задачи. Докажем аналог принципа максимума для параболических задач и ряд следствий- Экстремумы далее могут быть нестрогими-
Утверждение 1. Пусть сеточная функция с^ определена в Б N, удовлетворяет (5) в ОN и достигает в узле (п*,**) € ОN максимального или минимального в Б N значения с*. Тогда с* ограничено независимо от узлов сетки.
Доказательство. Пусть максимальное значение с* достигнуто в узле (п*,**) € ОN. Применим (5). Левая часть неотрицательна, первые три слагаемых правой части неположительны. Последнее слагаемое ЕП* также неположительно, если с* достаточно велико. Поэтому дтсП * = д^^сП * =0 в частности, и потому с^ *_ 1 = с^ *. Продолжаем рассуждение для узла (п* — 1,**) и далее, пока не получим, что сП * = ^(х^ *). Для отрицательного минимума доказательство такое же. □
Утверждение 2. Пусть сеточная функция с^ определена в О N и удовлетворяет (5)-(6). Тогда с^ ограничена в О N независимо от узлов сетки.
Доказательство. Пусть с^ достигает максимума с* в узле (п*,**) € О N. Если (п*,**) € ОN, то вывод следует из утверждения 1. Проверим граничные узлы дО N. Это узлы вида (0, *), (п, 0) и (п, 1). Пусть в узле вида (п*, 0) достигнут максимум с^ * = с*. Из (6) с учетом условий на О получаем, что д^с^ * > 0, и потому с^ * > с^ * при достаточно больших с*, что противоречит наличию максимума. Аналогично рассуждение для узлов вида (п*, 1). Начальное распределение ^(х) ограничено, поэтому с0 также ограничено. Следовательно, максимальное значение ограничено в О N независимо от узлов сетки. Для нижней оценки рассуждения аналогичны. □
Полученная оценка имеет физический смысл. Отметим, что параметры задачи можно переопределять для слишком больших по абсолютной величине значений с.
Утверждение 3. Пусть с^ удовлетворяет (5)-(6) в О N при достаточно малых Н и т. Тогда |д^с^ | ограничена в О N независимо от узлов сетки.
Доказательство. Произведем сеточное дифференцирование уравнения (5) по Н с учетом равенства д^с^1 = д^с^ + д^д^с^Н:
дт дьсП + дхВ дьсП = АПд2дьсП + (дхА + ВП + джВН)дЛ,дЛ,сП —
— ВПдл.дл.сП 1 + (джВ — дсЕ)дЛ,сП — дхЕсП. (7)
Производные коэффициентов — в точках вида х^ + еН. Введем сеточную функцию времени Мп рекуррентно: Мп = Мп_1/(1 + тдхВп),
Мо = 1. Здесь используем независимость дхВ от х. Функция Мп — сеточный аналог экспоненты: дтМп = —дхВМп. Рассмотрим левую часть (7), обозначая Ц" = д^с"/М„:
д^Мд^сП) + дхВд^с" = дт (М„д^М^ + дх-Вд^с"
М„_1 дт ЦП + ЦП дт М„ + дхВ д^с" =
= М„_1дтЦП - М„иПдхВ + дхВдьсП = М„_1дтЦП.
Поскольку дхВ отрицательна и ограничена, то при достаточно малом шаге т (обеспечивающем положительность знаменателя) функция Мп > 1. Запишем (7) в виде
ММ-!дт ЦП = АПд2 ЦП - В"дь^г1 + (дхВ - дсЕ)ЦП - МЕ+
+ (дхА + ВП + дхВН)дЛ,ЦП.
Пусть — максимальное положительное значение в узле (п, *) € ОN. Тогда левая часть неотрицательна; рассмотрим правую. Первые три слагаемых неположительны. По условиям либо Е = 0, либо дсЕ > 0 и отделено от нуля, а следовательно, при достаточно больших Ц^ сумма первых четырех слагаемых неположительна. Если величина Р = дхА + В^ + дхВН > 0, то последнее слагаемое также неположительно. Противоречие показывает, что положительный максимум функции Цгп, если и достигается в ОN, то ограничен независимо от Н и т, а только от функционала Е.
Рассмотрим случай Р < 0. Заметим, что д^Ц_1 > 0 и Нд£= д^ - дЛЦ‘_1 < дЛЦ‘. Поэтому Рд^ < РНд2и
ММ-1 дтЦ* < (А" + РН)д2Ц* - дхВдьсП + (дхВ - дсЕ)Щ - ^.
Поскольку А > А > 0, а Р ограничено, то при достаточно малом Н множитель в скобках при д^Ц^ положителен; противоречие получается аналогично.
Теперь проверим ограниченность Ц^ в граничных узлах; в силу (6)
Ц0 = дьсО/Мо = дЛс0 = ^(х< + 0), 0 < в < Н
ограничена; а ограниченность G влечет ограниченность производных дл.сП-1 и dhc^ и, следовательно, функций и ЦП.
Итак, мы доказали, что максимальное положительное значение функции Un в D'N ограничено независимо от шагов сетки, если они достаточно малы. Отсюда следует и ограниченность положительного максимума функции д^с^ = МпЦП. В самом деле, при т < то (при некотором то > 0)
T Т_
м„ < (1 - ^т) Т < (1 - ^то) ,
где = max |дхВ|. Здесь использована монотонность функции (1 —
x)-1/x при x € (0,1), а обозначение max означает оценку сверху мно-
жества значений функционала.
Для отрицательного минимума доказательство дословно такое же.
□
Утверждение 4. Пусть сП удовлетворяет (5)-(6) в Dn при достаточно малых h и т. Тогда |дтс^| ограничена в D◦ независимо от узлов сетки.
Доказательство. Пусть максимальное по модулю значение достигается в узле (n, i) € Dn; предположим, что это максимум (для минимума рассуждение такое же) и продифференцируем (5) сеточно по т:
дт дт сП = А^дт сП + В^дт сП — В^дт С-1 — дсЕдт с^+
+ д4АдЬСП + д4ВдЛ,сП — д*ВдЬСП 1 — д*Е <
< д4Ад2сП + д4ВдЬСП — д4ВдЛ,сП 1 — д4Е (8)
Левая часть неотрицательна, поэтому последнее выражение не меньше нуля.
Без ограничения общности считаем, что |д^с^| < 1 в DN. Для упрощения выкладок считаем В = E = 0 (общий случай рассматривается аналогично). Выразим д^с^ из (5) и подставим в полученное неравенство:
д А д А
0 < -t-дтсП — -t-B^n + д;Вд^ < АА
дА
< д4АдтсП-;Г~ВПдЬсП + K1(B) + K2(B)дтсП.
А
Учитывая отрицательность суммы д4А + К2(в), получаем, что при достаточно больших дт сП правая часть отрицательна, что приводит к противоречию. Следовательно, если максимальное в значение функции дтеП достигается в , то оно ограничено независимо от шагов сетки.
Проверим граничные узлы. Для п = 1 рассмотрим (5) с учетом С1 = с0 + тдт :
дт с1 = ^1с0 + т^1_дт с1 — тЕ1.
Если дт с1 > 0 — максимальное значение, то дт с1 < с0 — тЕ1. Однако правая часть ограничена независимо от шагов (это следует из ограниченности коэффициентов уравнения (2) и второй производной начального распределения у>(ж)). Для г = I или г = 0 продифференцируем на сетке соответствующее условие (6):
д^дтеП = д4С + дсС(£,0, и, е)дте^.
С учетом ограниченности частных производных и положительности дсС получаем, что д^дтсП > 0 при достаточно большом положительном максимальном значении дтсП, что противоречит максимальности. Второе граничное условие рассматривается аналогично.
В итоге сеточная производная дт сП ограничена в Б N независимо от шагов сетки. □
Решение системы. Сеточная задача (5)—(6) представляет собой систему нелинейных уравнений особой структуры: подсистемы для каждого п можно решать последовательно (послойно), каждое уравнение содержит не более трех неизвестных (считаем, что ет при т < п уже известны), более того: лишь одна неизвестная входит в каждое уравнение нелинейно. Для решения этой системы применим нелинейный аналог метода прогонки.
Пусть все значения при т < п уже известны; тогда определены коэффициенты АП, ВП и ВП сеточной задачи. Коэффициенты ЕП, СП и СП запишем, явно выделив значения на неизвестном слое: еП(сП) = ЕП, 9П(еП) = СП, 9П(еП) = СП. Выразим из второго уравнения (6) е^:
СП = %0(сП) + СП
Эта зависимость — непрерывная. Уравнения (5) позволяют непрерывно выразить СП+1 через сП и С-1:
е*+1 = ^2дт СП — АПСП~1 + 2АПСП + ^ВПСП + ^-2В^ПдьСПТ1 + ^2еП(еП)
АП +
что, с учетом предыдущей формулы, означает непрерывную зависимость всех сП от с^. Обозначим ее сП = @П(сП). Подставляя эту функцию во второе граничное условие, имеем уравнение
дьвП-1(сП) + 9П (©П(сП))=0. (9)
Покажем, что уравнение (9) имеет решение. В самом деле, левая часть (9) — непрерывная функция от неизвестного е^. Поскольку решение в БП_1 уже построено, оно ограничено константой, не зависящей от узлов сетки: |е^| < Д. Выберем произвольное сП > Д; в силу (6) имеем д^еП > 0 и потому сП > сП > Д. Вычислим еП = ©П(сП) из
(5); если окажется, что Д < сП < сП > еП, то значение сП является максимумом в БП и больше Д, что невозможно в силу утверждения 1. Поэтому еП > сП > Д. Рассуждая аналогично, получим, что ©П(сП) = еП > сП-1 > Д. Это означает, что д^с^-1 > 0 и 9П(@П(сП)) > 0, то есть левая часть уравнения принимает положительное значение. Выбрав произвольное значение сП < —Д, аналогично покажем, что левая часть уравнения принимает отрицательное значение; по теореме Больцано, найдется хотя бы одно значение сП € (—Д, Д), удовлетворяющее уравнению. Далее находим все сП = ©П(сП) (по существу это рекуррентный алгоритм). Полученное решение в попадает в условия доказанных выше утверждений и потому ограничено вместе с первыми сеточными производными.
Сходимость аппроксимаций и обобщенное решение. Итак, сеточные функции сП , д^еП и дт сП равномерно ограничены независимо от шагов сетки. Построим в П функцию с(£, ж) кусочно-линейной интерполяцией сеточной функции еП. По построению с € Н1_(П) (со-болевское пространство функций, обладающих первыми обобщенными производными из Е2(П) [5, 6, 7]) при любых (достаточно малых) шагах Н и т, причем семейство ограничено по норме Н1_(П). Выберем последовательность шагов (т&, Н^), сходящуюся к нулю, и соответствующую последовательность £&(£, ж) (она ограничена в Н_(П)). Пространство Н гильбертово, поэтому из ограниченной последовательности можно выделить слабо сходящуюся подпоследовательность, за которой сохраним прежнее обозначение; ее предел обозначим с(£, ж) € Н1_(П). Эта функция непрерывна, так как равномерно ограниченные непрерывные 2^(£, ж), обладая равномерно ограниченными кусочно-постоянными производными, равностепенно непрерывны, что влечет равномерную сходимость (возможно, после перехода к подпоследовательности — теорема Арцела).
Назовем функцию е(£, ж) € Нх(П) обобщенным решением [5, 6, 7] краевой задачи (2)—(4), если е(£, ж) удовлетворяют интегральному тождеству
при любой непрерывной ^(£, ж) € Нх(П), такой, что «(Т, ж) = 0. Мотивировка определения обычная: гладкое решение краевой задачи тождеству удовлетворяет, а удовлетворяющая тождеству гладкая функция е(£, ж) является решением краевой задачи.
Построенная функция е(£, ж) является решением в смысле данного определения: на аппроксимациях тождества выполняются с погрешностью О (Л., т), а равномерная сходимость аппроксимаций и непрерывность коэффициентов задачи позволяет перейти к пределу. Предложенный сеточный метод решения краевой задачи сходится к обобщенному решению, существование которого тем доказано.
Таким образом, доказано существование обобщенного решения для класса нелинейных параболических краевых задач и построена сходящаяся неявная разностная схема. Класс задач включает, в частности, задачи со свободной границей, возникающие при моделировании взаимодействия газов с твердым телом.
Работа поддержана грантом РФФИ 11-01-00642-а «Разрешимость параболических квазилинейных краевых задач с нелинейными граничными условиями и сходимость разностных схем».
ь
/
/
о
, 0, с(-, •))«(£, 0)^(£, 0, е(£, 0), с(-, •))
о
т
, Ь, е(^, •))«(£, £)^(£, Ь, е(£, Ь), е(^, •))
о
+ / (В — В)«джс^ж^£ — Е(£, ж, с(£, ж), с(-, •))V в,жв;Ь
пп
Resume
We consider the one-dimensional quazi-linear parabolic Neumann boundary-value problem: coefficients of the partial differential equation and right-hand sides of the boundary conditions depend on time, point, and the history of the solution. Convergence of difference approximations to a weak solution to the problem is proved.
Список литературы
[1] Castro F. J., Meyer G. Thermal desorption spectroscopy (TDS) method for hydrogen desorption characterization (I): theoretical aspects // Journal of Alloys and Compounds. 2002. V. 330-332. P. 59-63.
[2] Заика Ю. В., Родченкова Н. И. Диффузионный пик ТДС-спектра дегидрирования: краевая задача с подвижными границами // Математическое моделирование. 2008. Т. 20. № 11. C. 67-79.
[3] Chernov I., Bloch J., Gabis I. Mathematical modelling of UH3 formation // International Journal of Hydrogen Energy. 2008. V. 33. P. 5589-5595.
[4] Чернов И. А. Математическое моделирование экзотермического формирования гидрида // Математическое моделирование. 2010. Т. 22. № 1. С. 3-16.
[5] Ладыженская О. А. Краевые задачи математической физики. М.: Наука, 1973.
[6] Михайлов В. П. Дифференциальные уравнения в частных производных. М.: Наука, 1976.
[7] Evans L. C. Partial differential equations. Providence: Amer. Math. Soc., 1991.
Институт прикладных математических исследований Карельского Научного Центра РАН,
185000, Петрозаводск, ул. Пушкинская, 11,
E-mail: [email protected]