CC BY
11
УДК 517.956.226
Т. О. Капустина1
О КРАЕВОЙ ЗАДАЧЕ ДЛЯ УРАВНЕНИЯ СМЕШАННОГО ТИПА: АСИМПТОТИЧЕСКИЙ АНАЛИЗ И ЭФФЕКТИВНЫЙ ЧИСЛЕННЫЙ АЛГОРИТМ
Рассматривается сингулярно возмущенная краевая задача для уравнения смешанного эллиптико-параболического типа. Первая часть работы посвящена асимптотическому исследованию решения задачи. Используется модификация метода пограничных функций для уравнений смешанного типа с малыми параметрами при старших производных. Вторая часть посвящена созданию численного метода, учитывающего структуру решения при малых значениях параметра. Применяется идея приближенного разложения эллиптического оператора на произведение двух параболических. Цель работы — на основе асимптотического представления решения разработать эффективный численный алгоритм.
Ключевые слова: сингулярные возмущения, уравнения смешанного типа, метод пограничных функций, факторизация оператора.
1. Введение. Уравнения смешанного типа возникают как математические модели физических процессов, где происходит взаимодействие двух сред с разными свойствами. В частности, модель стационарного распределения температуры в системе, состоящей из нагретого твердого тела, помещенного в движущуюся холодную жидкость, приводит к уравнению смешанного эллиптико-параболического типа [1]. Независимыми переменными в этом уравнении будут координаты рассматриваемой точки, а искомой функцией — стационарная температура в этой точке. Температура в твердом теле описывается эллиптическим уравнением, а в жидкости — параболическим, причем параболическое уравнение содержит первую производную по той из пространственных координат, которая отвечает за расстояние от рассматриваемой точки до границы раздела сред. Далее, естественно потребовать, чтобы температура и тепловой поток были непрерывны на границе между твердым телом и жидкостью. Это требование ведет к условиям склейки — непрерывности решения и его нормальной производной.
В данной работе рассматривается краевая задача для эллиптико-параболического уравнения с малыми параметрами при старших производных, с условиями Дирихле на соответствующей части границы области, а также двумя условиями склейки на линии изменения типа.
При численном решении задачи стандартным методом конечных разностей параболическая задача требует значительно меньше арифметических операций, чем эллиптическая. Цель работы — используя малость параметра, на основе асимптотического представления решения построить численный алгоритм, позволяющий найти приближенное решение быстрее и с меньшими затратами компьютерных ресурсов, чем стандартным численным методом.
Для создания эффективного численного алгоритма воспользуемся идеей факторизации эллиптического оператора из [2, 3]. При малых значениях параметра приближенно разложим эллиптический оператор на произведение двух параболических. Далее, заменим затратную для численного решения эллиптическую задачу двумя последовательными более экономными параболическими задачами с разным направлением "времени". Направление выбирается в зависимости от знаков коэффициентов уравнения так, чтобы обе параболические задачи были корректными. Для того чтобы начать вычисления по такому алгоритму, нужно знать значение решения и его производных на линии изменения типа. Поскольку эти значения неизвестны, заменим их асимптотическими приближениями необходимой точности относительно малого параметра. Асимптотическое представление построим с помощью метода пограничных функций [4], а также используем идеи применения асимптотических методов к уравнениям смешанного типа из [5-7].
1 Механико-математический факультет МГУ, доц., к.ф.-м.н., e-mail: kapustina-tatianaQyandex.ru
4 ВМУ, вычислительная математика и кибернетика, № 4
Таким образом, в данной работе предлагается идея сочетания асимптотического и численного подходов для приближенного решения сингулярно возмущенной задачи. Преимущество предложенного алгоритма по сравнению со стандартным состоит в значительном сокращении числа вычислительных операций. Статья является продолжением работ [8, 9], выполненных совместно с Ж.-П. Лоэаком.
2. Постановка задачи. Рассмотрим прямоугольник П = (0,1) х (—1,1) на плоскости (х,у). Введем обозначения: П+ = (0,1) х (0,1) — верхняя часть прямоугольника, П_ = (0,1) х (—1, 0) — нижняя, Б = (0,1) х {0} — линия, разделяющая П+ и П_. Через Г± обозначим боковую границу И : . а через 7 — нижнюю границу П_:
Г+ = {0} х (0,1) и {1} х (0,1), Г_ = {0} х (—1,0)и{1} х (—1,0), 7 = (0,1) х {-1}.
Рассмотрим уравнение смешанного типа, эллиптическое в П_ и параболическое в П+. Зададим граничное условие Дирихле на части внешней границы Г = Г_ и Г+ и7, а также условия склейки, т. е. условия непрерывности решения и его нормальной производной, на линии изменения типа Б:
Ми = аи + Ьиу — еихх = /, (ж, у) € 0+,
Ьи = си +(1иу - е(ихх+ иуу) = д, (ж,у)еП_,
и\г = 0, ^
V,-= И4. (%)-= (%) + |_5"
Здесь е — малый положительный параметр, а, Ь, /, с, й, д — гладкие функции переменных (ж, у); и± обозначает функцию и в областях В книге [1] доказано, что при положительных коэффициентах а, бис решение задачи (1) существует и единственно.
Как структура асимптотики решения, так и численный алгоритм существенно зависят от знака коэффициента й. Мы рассмотрим два случая: положительного и отрицательного й.
3. Асимптотическое представление решения. Построим асимптотическое представление решения задачи (1) с помощью метода пограничных функций [4] и воспользуемся идеями применения асимптотических методов к сингулярно возмущенным уравнениям смешанного типа [5-7].
Асимптотическое представление й решения и задачи (1) состоит из регулярной части и нескольких пограничных функций. Главное слагаемое и0 регулярной части удовлетворяет вырожденному уравнению (полученному из (1) при е = 0) с дополнительными условиями, выбор которых зависит от знака коэффициента й.
3.1. Случай отрицательного коэффициента <1. При ё < 0 вырожденная задача имеет вид
Х° + Ц| = /, (ж,у)еП+,
(®,у)€ 0_, (2)
(11у )_ =
Взяв первые два уравнения (2) на линии Б, получим систему линейных алгебраических уравнений относительно функций и°(х, 0) и иу(х, 0). При положительных коэффициентах а, Ь, с и отрицательном й определитель этой системы не равен нулю, следовательно, система однозначно разрешима. Решая два обыкновенных дифференциальных уравнения (2), первое в области П+, второе в области П_, с общим начальным условием и°(х, 0), находим функции г^ в областях соответственно. Решение вырожденной задачи и0 удовлетворяет двум условиям склейки на линии 5 и не удовлетворяет граничному условию на Г из (1). Следовательно, решение и(х,у,е) исходной задачи (1) имеет пограничные слои в окрестности части границы прямоугольника Г. Главное слагаемое асимптотики имеет вид (см. [10]):
v > \и-(х,у,£), (х,у) еО_, ^
й+(х,у,е) = и\{х,у) + д+ + '
й-(х,у,е) = и°_(х,у) + q- ( ) + ( ) + (х У + 1
£ J \ у/£ / \ £
\у/£ £ / \£ £ J \ Vе е } \ е е
где q, z, v — пограничные функции, a w, р, т, к — угловые пограничные функции, зависящие от соответствующих быстрых переменных.
Особенностью данной задачи является разный масштаб пограничных слоев. Слои в окрестности боковых границ Г1*1 имеют ширину (т. е. масштаб быстрой переменной) порядка у/е, а в окрестности нижней границы 7 — ширину порядка е. В вершинах прямоугольника (0, —1) и (1,-1) происходит пересечение этих слоев. Поэтому в окрестностях данных точек нужны угловые пограничные функции, по две в каждой вершине из-за разного масштаба взаимодействующих слоев [11].
Все пограничные функции либо могут быть вычислены в явном виде, либо для них могут быть получены двусторонние оценки. Все пограничные функции экспоненциально убывают при удалении от своего пограничного слоя или вершины. Например, оценки функций q и w имеют вид
где С; V обозначают быстрые переменные, С, а = const > 0. Оценки пограничных функций z и v аналогичны оценке для q, а оценки угловых функций р, т, к аналогичны оценке w с заменой V на соответствующие быстрые переменные.
Теорема. При всех достаточно малых £ функция й(х,у,е) является равномерным асимптотическим представлением решения и(х, у, е) задачи (1), т.е. справедлива оценка
тах_ (ж,з/)еп
и(х,у,е) - й(х,у,е) < В у/е, (4)
где постоянная В не зависит от £.
Доказательство теоремы основано на оценках пограничных функций и на принципе максимума для параболических и эллиптических уравнений [4, 5, 11].
При необходимости более высокой точности по е можно стандартным способом продолжить процесс построения асимптотики. Без дополнительных условий на данные задачи асимптотику можно построить до второго порядка [10].
3.2. Случай положительного коэффициента <1. При ё > 0 вырожденная задача ставится
так
(си° + ¿и® = (ж,у) е П_, «^1^ = 0,
+ = (ж, у) е П_|_, |5 = и°_ |5.
Сначала находится функция и°_ в области П_ с начальным условием на 7. Затем с использованием в силу первого условия склейки (1) уже найденного значения г^ в качестве начального условия на Б вычисляется в П+. В этом случае решение вырожденной задачи и0 непрерывно, имеет разрыв нормальной производной на 5, удовлетворяет граничному условию на 7 и не удовлетворяет на Г±. Решение и(х,у,е) исходной задачи (1) имеет пограничные слои ширины у/е в окрестности Г±, а также внутренний слой ширины е на 5. В окрестности точек (0,0) и (1,0) происходит взаимодействие пограничного и внутреннего слоев разного масштаба. Поэтому в окрестности этих точек асимптотика будет содержать угловые функции двух типов, аналогичные функциям и), р, т, к из первого случая. Все пограничные функции экспоненциально убывают при удалении от соответствующих линий или вершин. Для асимптотического представления справедлива оценка, аналогичная приведенной в теореме из п. 3.1.
4. Стандартный численный алгоритм. Нестандартными для стандартного численного метода будут условия склейки на Б. Чтобы применить метод конечных разностей, заменим (1) двумя последовательными задачами. Сначала рассмотрим эллиптическую задачу
(си + йиу - £{ихх + иуу) = д, (ж,у)еп_, 1Нг_и7 = 0' (аи + Ьиу - £ихх) |5 = /,
5 ВМУ, вычислительная математика и кибернетика, № 4
где в качестве граничного значения на S, в силу двух условий склейки (1), используем параболическое уравнение из П+. Затем рассмотрим параболическую задачу
Г(ш + buy - еихх = /, (ж, у) G П+,
^ I Г - ^ 1 ^ | (Ч - W - | {У Ч
где в качестве начального условия на S, в силу первого условия склейки (1), берем значение и-из эллиптической задачи (5).
Численное решение уравнения смешанного типа (1) получаем, последовательно применяя стандартный метод конечных разностей к задачам (5) и (6). Известно, что при решении линейных систем, полученных дискретизацией задач (5) и (6), число арифметических действий для параболической задачи значительно меньше, чем для эллиптической (см., например, [12]).
5. Факторизация эллиптического оператора. Поскольку численное решение параболической задачи требует меньших затрат, чем решение эллиптической задачи, то для создания эффективного численного алгоритма при малых е приближенно разложим эллиптический оператор Ь задачи (5) на произведение двух параболических. Воспользуемся идеей факторизации эллиптического оператора и оценками из [2, 3]. Факторизация и численное решение рассматриваются для случая постоянных коэффициентов с, й.
Запишем разложение решения и задачи (5) в ряд Фурье по переменной ж:
сю
и(х,у) = (у) зт(кттх).
к=1
Применение к и эллиптического оператора Ь дает для каждого к коэффициент
(^¿з +й^7у + (с + Ф^)2)) (Заполучим обыкновенный дифференциальный оператор по у. Обозначим £ = ки и запишем символ этого оператора как А(А) = —еХ2 + й\+ (с + е£2). Далее, разложим корни многочлена А до первого порядка по степеням е:
л+=-Нс -4-+ее)+0(е2)' =-К-т - *+4 - +о(е2)'
тогда А(А) = —е(Л — А+)(А — А~) + 0(е2). Теперь применим обратное преобразование Фурье (при этом —д2/дх2) и получим следующее приближенное разложение эллиптического оператора Ь на произведение параболических операторов:
¡^£й~2Р-Р+ + 0(е2), <¿>0, ..
\^ей-2Р+Р- + 0{е2), й<0,
/ с2\ д д2 ( (I2 с2 \ д д2 /' ' I С Г — + <1----7. Р =---г + Г — + <1 —
й2) ду ° дх2 ' " Vе ^ ' ^ (Р) ' " ду ' £ дх2 '
Выбор между формулами (7) зависит от знака коэффициента й.
Если й > 0, то Р+ — это прямой параболический оператор, задача для него будет корректна при у, возрастающем от — 1 до 0; Р~ — обратный параболический оператор, задача для которого корректна при у, убывающем от 0 до — 1. Если й < 0, то Р+ и Р~ меняются ролями.
Таким образом, для повышения эффективности численного алгоритма мы заменяем эллиптическую задачу (5) двумя последовательными параболическими задачами: первая — для у, убывающего от 0 до — 1, вторая — для у, возрастающего от — 1 до 0.
6. Эффективный численный алгоритм. Эффективный алгоритм основан на факторизации эллиптического оператора и на асимптотическом представлении решения. Обозначим через V приближенное решение задачи (5), (6). Рассмотрим два случая в зависимости от знака коэффициента (1.
6.1. Случай отрицательного коэффициента <!. В случае й < 0 используем вариант Ь и —£<1~2Р+Р~ разложения (7) эллиптического оператора Ь на произведение двух параболических операторов. Заменим уравнение Ьи = д уравнением —ей~2Р+ {Р~у) = д. Полагая ад = Р~у, получим два параболических уравнения
—ей~2Р+,ш = д, Р~у = ад,
которые надо решить последовательно.
Для того чтобы начать вычисление, нужно задать начальное значение функции ад на линии Б и граничное значение на Г_. Граничное условие следует из (1). Поскольку г>|г = 0, то
уу\г = %|/|р = 0. Рассматривая уравнение Ьу = д на Г_, получим
Чг_ = (Р~*)\г- = _ = = (8)
Теперь найдем начальное условие на Б. Из (5) имеем еухх\3 = (ау + Ьуу — /)|5, следовательно,
и}\я = (Р г>)|„ =
с12 с2
— с + а + е-^ ) V + (й + Ь)уу — /
т.е. для вычисления начального условия потребуются значения и|5 и Так как они неизвест-
ны, заменим их точные значения асимптотическими приближениями, полученными в п. 3. При постоянных коэффициентах с, ё пограничные функции д и г могут быть вычислены в явном виде. Используем значение главного слагаемого асимптотики й, заданного формулой (3), на линии Б. Согласно оценке (4) имеем
и\
= и°(х, 0) + С1! ехр + С12 ехр ^ +0(^),
где С\у2 и «1)2 — известные константы, а1)2 > 0. Нормальная производная йу |5 определяется аналогичной формулой. При необходимости более высокой точности вычислим столько слагаемых асимптотического представления, сколько нужно для достижения требуемой точности. Итак, в качестве начального условия для ад на линии Б вместо Р~у берем
ад
= (Р~й)\3. (9)
Следующий шаг — численное решение двух последовательных параболических задач, заменяющих эллиптическую часть (5). В первой задаче вычисляем промежуточную функцию ад в нижней области П_. Эта функция удовлетворяет параболическому уравнению с оператором Р+, где переменная у убывает от 0 до —1, с граничным условием (8) и начальным условием (9):
Р+п) = 1д, (х,у) б П_, у б (0, -1),
ад1г =—д, ги|5 = (Р~й) |5. ^ '
Вторая задача — это вычисление функции V в нижней области П_. Функция V удовлетворяет параболическому уравнению с оператором Р~, где переменная у возрастает от — 1 до 0, с правой частью ад, уже найденной в первой параболической задаче (10), и с заданными в (1) начальным и граничным условиями:
'р-у = гц, (ж,у) € 0_, у € (-1,0),
^|г_=о, «|7 = о. ^^
Осталось вычислить функцию V в верхней области П+. Функция V удовлетворяет параболическому уравнению с оператором Л /. где переменная у возрастает от 0 до 1, с заданным в (1) граничным условием и начальным условием на Б — функцией у_, найденной во второй параболической задаче (11):
'Му = /, (ж,у) € 0_|_, у € (0,1), Чг+=0> 45 = 456 ВМУ, вычислительная математика и кибернетика, № 4
6.2. Случай положительного коэффициента 6,. При 6, > 0 действуем аналогично, меняя ролями операторы Р+ и Р~. Прежде всего, с помощью асимптотического представления (3) вычисляем начальное значение функции т — Р+ь на Я, а из (1) выводим граничное условие на Г_: "ш|г — 9- Далее рассматриваем две параболические задачи, заменяющие эллиптическую в нижней области
Сначала решаем параболическое уравнение для у, убывающего от 0 до —1:
р-ю = -сРе~хд, (х,у) у е (0,-1),
— (I, 'Ш| „ = (р + м) |„.
(13)
Затем решаем параболическое уравнение для у, возрастающего от —1 до 0:
р+у = и>, (ж,г/)еп_, у е (-1,0),
-у I _ = 0, г>| = 0. 11 _ 17
(14)
И наконец, решаем параболическую задачу в верхней области для у, возрастающего от 0 до 1:
= (х,у) е 0+, У е (0,1),
I П I I (15)
Заметим, что при построении численного алгоритма мы дважды воспользовались малостью параметра е: один раз — при приближенной замене эллиптического оператора произведением параболических операторов, и второй раз — при замене начального условия его асимптотическим представлением. Обе замены дают хорошее приближение при достаточно малых значениях е.
7. Результаты численного эксперимента. Здесь представлены графики решения V задачи (10)—(12) для случая отрицательного 6, и е = 0.05. Вычисления проведены методом конечных разностей на сетке с постоянным шагом с помощью эффективного алгоритма, построенного в п. 6.1. На рис. 1 показан график ■и(х1у) при у — —0.3. На рис. 2-4 представлены графики ■и(х1у) для нескольких значений х: х — 0.1 (рис. 2), х — 0.2 (рис. 3), х — 0.5 (рис. 4).
Заметим, что структура решения, полученного в результате численного эксперимента, расположение и масштаб пограничных слоев полностью согласуются с асимптотическим представлением (3). Таким образом, предложенная в работе идея позволяет объединить асимптотический анализ и численные методы и создать эффективный вычислительный алгоритм для сингулярно возмущенной задачи.
СПИСОК ЛИТЕРАТУРЫ
1. Джура о в Т.Д. Краевые: задачи для уравнений смешанного и смешанно-составного типов. Ташкент: Фан. 1979.
2. Loheac J.-P., Nataf F.. Schatzman M. Parabolic approximations of the convection diffusion equation // Mathematics of Computations. 1993. 60. N 202. P. 515 530.
3. Loheac J.-P. Some numerical applications of pseudo-differential computation // Proceedings of International Miniconference "Qualitative Theory of Differential Equations and Applications". M.: MESI. 2010. P. 132 140.
4. Васильева А.Б.. Бутузов В.Ф. Асимптотические методы в теории сингулярных возмущений. М.: Высшая школа. 1990.
5. Сушко В.Г. Асимптотические разложения решений бисингулярных задач: Дис. д.ф.-м.н. М.. МГУ. 1997.
7 ВМУ, вычислительная математика и кибернетика, № 4
6. Сушко В. Г. Асимптотические решения некоторых сингулярно возмущенных уравнений смешанного типа // Фундаментальная и прикладная математика. 1997. 3. № 2. С. 579-586.
7. Розов Н.Х., Капустина Т.О. Одна краевая задача для эллиптико-параболического уравнения // Дифференц. уравн. 2001. 37. № 6. С. 847-848.
8. Loheac J.-P., Kapustina Т. Asymptotic and numerical analysis of mixed type equation // Proceedings of International Miniconference "Qualitative Theory of Differential Equations and Applications". M.: MESI, 2013. P. 179-189.
9. Лоэак Ж.-П., Капустина Т. О. Асимптотическое и численное исследование эллиптико-параболического уравнения // Таврический вестник информатики и математики. 2015. 26. № 1. С. 50-61.
10. Капустина Т. О. Об асимптотическом решении сингулярно возмущенной краевой задачи для уравнения смешанного типа // Ученые записки Таврического Национального ун-та им. В. И. Вернадского. Сер. Математика. Механика. Информатика и кибернетика. 2006. 19(58). № 2. С. 31-34.
11. Бутузов В.Ф. Угловой погранслой в сингулярно возмущенных задачах с частными производными // Дифференц. уравн. 1979. 15. № 10. С. 1848-1862.
12. Тихонов А.Н., Самарский A.A. Уравнения математической физики. М.: Наука, 1966.
Поступила в редакцию 21.06.17
УДК 519.63:535.42
Ю. А. Еремин, И. В. Лопушенко2
УЧЕТ ЭФФЕКТА НЕЛОКАЛЬНОСТИ ПРИ РАССЕЯНИИ СВЕТА НА ПЛАЗМОННЫХ НАНОЧАСТИЦАХ В ГИБРИДНОЙ СХЕМЕ МЕТОДА ДИСКРЕТНЫХ ИСТОЧНИКОВ
Рассмотрена трехмерная задача дифракции поля плоской волны на плазмонной на-ночастице с учетом эффекта нелокальности. Решение построено на основе модифицированной вычислительной схемы метода дискретных источников. Проведено численное исследование влияния эффекта нелокальности на рассеивающие свойства сферических наночастиц при их деформации.
Ключевые слова: метод дискретных источников, нанооптика, дифракция, эффект нелокальности.
1. Введение. В настоящее время ведутся активные исследования в области проектирования и оптимизации наноплазмонных структур с использованием математического моделирования, основанного на классической теории Максвелла. Если размер рассматриваемой структуры составляет менее 10 нм, то классическая электродинамическая теория не применима для описания возникающих физических эффектов. В этом случае необходим учет так называемого эффекта нелокальности (ЭНЛ).
В силу выдающихся успехов наноплазмоники и ее широкого внедрения во многие области человеческой деятельности, начиная от диагностики и лечения заболеваний до жидкокристаллических дисплеев и солнечных элементов, проблема учета ЭНЛ становится все более актуальной [1]. Как было установлено, ЭНЛ может существенно исказить рассеивающие свойства плазмонных структур. Задача исследователя, как правило, состоит в необходимости точно определить положение максимума и интенсивность плазмонного резонанса, происходящего в наноструктурах. Это положение зависит от материала их составляющих, от формы частиц и их размеров, от свойств окружающей среды и поляризации внешнего возбуждения. Учет ЭНЛ может существенно изменить положение и интенсивность пика плазмонного резонанса (ПР), а также структуру ближнего поля, которая может включать в себя неизлучающие компоненты полей [2, 3].
1 Факультет ВМК МГУ, вед. науч. сотр., д.ф.-м.н., e-mail: ereminQcs.msu.su
2 Физический факультет МГУ, асп., e-mail: lopushenko.ivanQphysics.msu.ru
* Работа выполнена при финансовой поддержке G-RISC, программа DAAD, проект № М-2017а-6.
