Димова стефка николаева численное исследование нестационарных тепловых структур
На правах рукописи
УДК 519.6, 517.9 Димова Стефка Николаева ЧИСЛЕННОЕ ИССЛЕДОВАНИЕ НЕСТАЦИОНАРНЫХ ТЕПЛОВЫХ СТРУКТУР Специальность 05.13.18 – математическое моделирование, численные методы и комплексы программ
АВТОРЕФЕРАТ
диссертации на соискание ученой степени доктора физико-математических наук
Дубна – 2004
Работа выполнена в Лаборатории информационных технологий Объединенного института ядерных исследований
Научный консультант:
доктор физ.-мат. наук, чл.корр. РАН, проф. Курдюмов Сергей Павлович.
Официальные оппоненты:
• доктор физ.-мат. наук, проф. Абрамов Александр Александрович;
• доктор физ.-мат. наук, проф. Вабищевич Птр Николаевич;
е • доктор физ.-мат. наук, проф. Гулин Алексей Владимирович.
Ведущая организация: Институт прикладной математики им. М.В.Келдыша РАН, г. Москва.
Защита диссертации состоится 22 октября 2004 г. в 14:00 ч. на заседании Дис сертационного совета Д720.001.04 в Объединенном институте ядерных исследований (Лаборатория информационных технологий), г. Дубна, Московской области.
С диссертацией можно ознакомиться в библиотеке Объединенного института ядерных исследований.
Автореферат разослан 21 сентября 2004 г.
Ученый секретарь Совета:
кандидат физико-математических наук Иванченко З.М.
1
ОБЩАЯ ХАРАКТЕРИСТИКА РАБОТЫ
1.1 Актуальность работы Новые тенденции в современной науке наиболее ярко проявляются в синерге тике [1],[2]. Одной из главных задач синергетики – выявление общих принципов, лежащих в основе процессов самоорганизации в системах самой разной приро ды: физических, технических, биологических, социальных. При всем различии таких систем, все они характеризуются несколькими общими признаками, среди которых определяющими являются открытость, диссипативность и нелинейность [2].
Открытость системы означает наличие в ней обмена веществом и энергией с окружающей средой не только через границы системы, а также через объемные ис точники и стоки.
Диссипативность системы означает наличие в ней рассеивающего, размываю щего неоднородности фактора (теплопроводность, диффузия, дисперсия и т.д.).
Нелинейность системы выражается нелинейной зависимостью диссипирующего фактора от состояния среды. В такой среде действие нелинейных объемных источ ников и стоков создает и может усиливать ее неоднородность. Взаимодействие этих двух нелинейных факторов (с одной стороны сглаживание неоднородностей благо даря диссипации, а с другой – их усиление объемными источниками и стоками) мо жет приводить к локализации процессов на определенных участках среды – т.е., к возникновению структур, к возникновению самоорганизации. Чтобы отразить влияние диссипации на формирование структур, Пригожин ввл понятие диссипа е тивной структуры [3]. Действие стоков приводит к возникновению стационарных диссипативных структур (стоячие волны, ячейки Бенара), а действие источников – к нестационарным, эволюирующим диссипативным структурам (тепловые струк туры, в частности Т-слой [4], бегущие и спиральные волны, вихри).
Благодаря нелинейности в одной и той же среде без изменения е параметров е могут возникать разные структуры. Но не произвольные структуры, а свойственные только этой среде. Одна из важнейших и актуальных задач синергетики – опреде ление спектра структур, которые могут возникать и самоподдерживаться в открытых нелинейных системах.
Важный класс диссипативных структур образуют сильно нестационарные структуры [5]. Они возникают в результате сверхбыстрых процессов в нелинейных открытых системах с объемными источниками (так называемые системы с поло жительной нелинейной обратной связью). Характерные величины в таких про цессах (температура, концентрация, энергия) возрастают на несколько порядков за конечное время. Это время принято называть временем обострения, а сам процесс – режимом с обострением (blow-up).
Проблема режимов с обострением поставлена в 1940-х и 50-х годах в связи с тео рией цепных реакций Семенова и с теорией горения (Зельдович, Баренблатт, Либро вич). Их интенсивное изучение в 70-х годах началось благодаря предложенному на основе численных экспериментов [6] процессу сверхсжатия центральных частей кап ли из дейтерия и трития путем ее облучения профилированным во времени лазерным импульсом. Режимы с обострением рассматриваются в области физики плазмы [7], лазерного термоядерного синтеза [8], магнитной гидродинамики [9], астрофизики.
Методология решения “задач на обострение” позволяет с нетрадиционной точки зре ния рассмотреть ряд классических задач. Она открывает новые подходы к задачам коллапса (быстрое сжатие вещества), химической кинетики, метеорологии (катастро фические явления в атмосфере Земли), экологии (рост и вымирание биологических популяций), эпидемиологии (вспышка инфекционных заболеваний), экономики (фе номен бурного экономического развития) и т.д.. Поэтому разработка надежных и эффективных методов численного решения задач на обострение является актуальной проблемой, выходящей за рамки конкретного класса задач.
Механизмы самоорганизации можно проследить, изучая сравнительно простые на первый взгляд математические модели. Вся синергетика работает с несколькими такими моделями – это нелинейные уравнения (или системы уравнений) типа реакция-диффузия с нелинейными источниками и стоками.
Одна из самых богатых моделей – модель тепловых структур. В общем виде она выглядит так:
N (ki (u)uxi )xi + Q(u), t 0, x RN (1) ut = i= Эта модель описывает процессы распространения тепла и горения в среде с коэффи циентами теплопроводности ki (u) 0 и источником Q(u) 0, которые являются нелинейными функциями температуры u(t, x) 0.
Модели вида (1) в различном контексте изучались многими исследователями.
Большее количество работ посвящено полулинейным уравнениям: ki (u) 1, Q(u) = eu (уравнение Франка-Каменецкого), Q(u) = u, 1. После пионерской работы H. Fujita (1966 г.), они и их обобщения изучались интенсивно многими авторами. Сре ди них: J. Bebernes, A. Bressan, H. Brezis, D. Eberly, A. Friedman, В.А. Галактионов, И.М. Гельфанд, M.A. Herrero, R. Kohn, Л.А. Лепин, С.А. Посашков, А.А. Самарский, J.L. Vazquez, L.J.L. Velazquez. Монография [28] содержит обширную библиографию и отражает часть этих исследований.
Квазилинейное уравнение изучалось в работах D.G. Aronson, A. Friedman, H.A.
Levine, S. Kaplan, L.A. Pelitier и др.. Большой вклад ученых русской школы. Необыч ный эффект локализации граничных режимов с обострением обнаружен численным экспериментом в работе А. А. Самарского и М. И. Соболя 1963 года [10]. Проблема локализации в квазилинейных уравнениях с источником поставлена С.П. Курдю мовым в 1974 г.. Работы И.М. Гельфанда, А.С. Калашникова, ученых школы А.А.
Самарского и С. П. Курдюмова (их работы цитированы ниже) посвящены исследова нию интереснейших проблем физического и математического характера, связанных с моделью (1) или ее обобщениями. Среди них: локализация (строгая и эффективная) процесса горения в пространстве, развитие разных типов режимов с обострением, возникновение структур – бегущие и стоячие волны, сложные структуры с различ ной степенью симметрии. Секретом успеха этих исследований явилось сочетание вычислительного эксперимента с применением и развитием качественных и аналитических методов теории обыкновенных и в частных производных дифференциальных уравнений, теорий групп Ли и Ли-Бэклунда.
Особое место в этих исследованиях занимает изучение разных типов автомодель ных и инвариантно-групповых решений уравнения (1) со степенными нелинейно стями:
ki (u) = ui, Q(u) = u (2) Такой выбор не случаен. Зависимости вида (2) от температуры встречаются в мно гих реальных процессах [11]. Например, в случае i = = 2.5, 5.2, уравнение (1) моделирует термоядерное горение с электронной теплопроводностью в плазме [7];
значения = 0, 2 3 соответствуют моделям автокаталлитических про цессов с диффузией в химических реакторах;
6.5 соответствует радиационной теплопроводности высокотемпературной плазмы в звездах и т.д..
Установлено [31],[17], что на классе степенных функций симметрия уравнения (1) в определенном смысле максимальна, оно допускает богатый набор инвариантно групповых решений. Нестационарные диссипативные структуры режимов с обостре нием, рассматриваемые в предлагаемой работе, связаны со степенными автомодель ными решениями или с более общими инвариантными решениями, когда по радиусу автомодельность степенная, а по углу – бегущая волна. Исследования в области дис сипативных структур дают основание считать, что именно инвариантные решения описывают “аттракторы” эволюции диссипативных структур и тем самым характе ризуют важные внутренние свойства нелинейной диссипативной среды.
Этот богатый запас инвариантных решений уравнения (1) со степенными нели нейностями необходим для успешного применения методов исследования того же уравнения в случае более общих зависимостей k(u), Q(u). С помощью методов опе раторного сравнения и стационарных состояний удается исследовать свойства реше ний (такие как локализация, неограниченность, асимптотическое поведение) целых классов более общих нелинейных уравнений.
Наконец, в случае степенных коэффициентов есть нужные соотношения между диссипацией и источником, при которых происходит их согласование. В результате этого возникают сложные структуры, более того, возникает спектр структур, кото рые горят согласованно, с одним моментом обострения. Чтобы установить в системе эту сложную организацию и вывести ее на один из возможных путей эволюции, надо знать что заложено в ней и спровоцировать этот путь заданием начальных данных из области притяжения аттракторов-структур. Поэтому нахождение аттракторов процессов в системах различной природы является очень актуальной про блемой теории самоорганизации.
1.2 Цели работы Основная цель работы – исследование условий возникновения и эволюции силь но нестационарных тепловых структур в среде, описываемой уравнением нелинейной теплопроводности со степенными нелинейностями и некоторыми его модификация ми. В частности, это выявление и исследование возможных способов сложной ор ганизации нелинейной диссипативной среды. И так как исследования проводятся в основном численным экспериментом, то важной задачей является также создание и апробирование эффективных вычислительных алгоритмов для рассматриваемых классов задач:
– нелинейных параболических задач с сингулярными по времени решениями – за дач Коши и краевых задач для квазилинейных параболических уравнений, системы таких уравнений и полулинейных уравнений;
– автомодельных нелинейных задач – краевых задач для квазилинейных обыкно венных дифференциальных уравнений (ОДУ), системы таких уравнений и квазили нейных эллиптических уравнений.
1.3 Научная новизна работы Впервые реализованы численно инвариантные решения, которые описывают на правленное распространение тепла и горения в двумерной нелинейной анизотропной среде со степенными коэффициентами теплопроводности и степенными источника ми. Так конструктивным образом доказано их существование.
Найдены два новых типа решений автомодельной задачи (так называемых соб ственных функций горения нелинейной среды) – “с перетяжкой” и с полостью около центра симметрии, в цилиндрически- и сферически-симметричном случаях при па раметрах среды + 1. На основе детального численного исследования предло жен сценарий поведения сложных собственных функций (четных и нечетных) при + 1 + 0.
Открыт новый мир структур при параметрах среды + 1. Это удалось благодаря:
– изменения идеологии – искать структуры, “живущие” и развивающиеся не на нулевом фоне, а на гомотермическом фоне (т.е. изменения математической модели тепловых структур постановкой новых условий на бесконечности);
– нахождения подходящего (комплекснозначного) разделения переменных в лине аризованной автомодельной задаче и тем самым нахождения линейных приближений к спиральным с.ф.;
– детального аналитического и численного исследования линейных приближений, которое показало, что их асимптотики автомодельные.
Все это дало возможность замкнуть автомодельную задачу в случае +1 под ходящими краевыми условиями, конструировать начальные приближения к различ ным решениям и решить ее численно. Таким образом впервые показано, что усложне ние организации среды может происходить не только за счет локализованных во вне и сходящихся к центру структур, но и за счет расходящихся волн горения сложной структуры, в том числе раскручивающихся спиральных волн.
1.4 Достоверность результатов Проведено детальное исследование поведения решений при изменении как физиче ских, так и вычислительных параметров для каждого из рассмотренных классов задач. В частности, расчеты проводились на последовательностях вложенных сеток и таким образом исследовался их порядок точности;
когда это было возможно, чис ленные результаты сравнивались с известными аналитическими или качественными результатами.
Новые, неожиданные результаты подтверждены разными методами.
Разработанные методы и алгоритмы для параболических и для эллиптических задач имеют внутренний критерий точности:
– восстановление в параболической задаче времени обострения, заложенном в со ответствующей автомодельной задаче, если начальные данные автомодельные;
– сохранение автомодельности вплоть до момента обострения при использовании точных автомодельных начальных данных.
Полученные новые типы решений автомодельной задачи в радиально-симметрич ном случае инициировали исследования других авторов [12]–[14] другими методами (методом динамической аналогии, бифуркационным анализом). Цитированные ис следования подтвердили существование этих новых решений.
1.5 Практическая ценность работы Математическая модель тепловых структур возникла в связи с решением акту альных технологических проблем физики плазмы, лазерного термоядерного синтеза, магнитной гидродинамики и др.. Как было отмечено выше, режимы с обострением встречаются также во многих других областях. Разработанная здесь и тщательно апробированная методика решения задач на обострение применима к моделям фи зики плазмы, моделям автокаталлитических реакций, в теории горения (модели вос пламенения твердого и жидкого горючего). Она легко может быть перенесена и на другие классы задач. Отметим особо метод построения адаптивных сеток по про странственным переменным, согласованных с известными асимптотиками процессов в среде.
1.6 Апробация результатов диссертации Основные результаты диссертации докладывались и обсуждались на семинарах в ЛИТ (ЛВТА), ОИЯИ, Дубна;
в ИПМ им. "М.В. Келдыша", РАН, Москва;
в ВЦ РАН, Москва;
в ФМИ Софийского университета “Св. Кл. Охридского”;
в Институте математики и информатики Болгарской АН;
в Университете в Линфилд, Орегон, США;
в Департаменте по Прикладной математике Университета в Коллорадо, Бо улдер, США;
на коллоквиуме Факультета математики Университета в Кайзерслау терне, Германия;
в “Макс Планк” институте по физике плазмы, Мюнхен, Германия;
в Техническом Университете в Гетеборге, Щвеция, и на следующих международных конференциях:
Numerical Methods and Applications, 1988, 1994, 1998, София;
2002, Боровец, Бол гария;
Qualitative Theory of Dierential Equations, март 1989, Рига, Латвия;
International IMACS Conference on Mathematical Modelling and Applied Mathema tics, июнь 1990, Москва;
International Colloquium on Dierential Equations, август 1991, Пловдив, Болгария;
Elliptic Boundary Value Problems, March 1992, Росток, Германия;
International Colloquium on Numerical Analysis, август 1992, 1996, 1997, Пловдив, Болгария;
International Conference on Programming and Mathematical Methods for Solving Physical Problems, июнь 1993, Дубна;
International Conference “Criteria of Self-organization in Physical, Chemical and Biological Systems, июнь 1995, Москва-Суздаль;
IMA Minisymposium “Mathematical Investigations of Models in Combustion”, ноябрь 14-17, 1999, Минеаполис, США;
International Congress on Mathematical Modelling, октябрь 2002, Дубна.
Исследования проводились в соответствии с тематическими планами ЛИТ ОИЯИ, ИМИ Болгарской АН, ФМИ Софийского Университета и были поддержаны гранта ми Министерства образования, науки и технологий Болгарии, и Научного фонда Софийского университета.
1.7 Публикации По материалам диссертации опубликовано 22 работы, из которых 5 работ в рос сийских ведущих научных журналах.
1.8 Структура и объём диссертации Диссертация состоит из введения, четырех глав и заключения. Полный объем диссертации – 211 страниц, включая 65 рисунков, 22 таблицы и список литературы, содержащий 138 наименований.
1.9 Личный вклад автора Автор диссертации, работая в коллективе соавторов из Института Прикладной Математики, РАН, Москва, из Института математики и информатики БАН, София, из Факультета математики и информатики Софийского университета, в большин стве случаев был инициатором данных исследований. Двое из соавторов – Даниела Василева и Милена Колева – были аспирантками автора диссертации (вторая ру ководимая совместно с проф. М. Касчиевым) и защитили успешно кандидатские диссертации. Все аналитические исследования, приведенные в диссертации, являют ся личным вкладом автора. Автор непосредственно участвовала в математической постановке задач, в разработке численных методов и алгоритмов, в их программной реализации, а также в численном решении поставленных задач, в обработке, анализе и интерпретации результатов численных экспериментов.
2 ОСНОВНОЕ СОДЕРЖАНИЕ РАБОТЫ
Во Введении oбоснована актуальность проблем, исследованных в диссертацион ной работе, и сформулированы цели работы. Введены необходимые понятия, пере числены основные задачи и сделан обзор известных автору результатов, связанных с ними. Коротко изложены и прокомментированы полученные в диссертации резуль таты.
Основные используемые понятия введены на самой простой модели вида (1), когда предполагается радиальная симметрия. В этом случае задача Коши с финитными начальными данными формулируется так:
ut = x1N (xN 1 u ux )x + u, (3) x R+, t 0, 0, 1, (4) ux (t, 0) = 0, u(0, x) = u0 (x) 0, 0 x l, u0 (x) 0, x l.
Если начальные данные удовлетворяют дополнительным условиям: u0 (x) C(R+ ), (u u0 )(0) = 0, то существует единственное локальное (по времени) обобщенное ре шение u = u(t, x) задачи (3)-(4), которое является неотрицательной непрерывной функцией в R+ (0, T0 ), где T0 (0, ] – конечное или бесконечное время суще ствования решения (см. список литературы в обзоре [15]). При этом u(t, x) является классическим решением в окрестности любой точки (t, x), где u(t, x) является стро го положительным. В точках вырождения оно может не иметь необходимой глад кости, но тепловой поток xN 1 u ux должен быть непрерывным. Это означает, что u ux = 0 всюду, где u = 0. Когда время T0 конечное и limtT0 supxRN u(t, x) = +, T0 называется временем обострения неограниченного (сингулярного по времени) решения.
Неограниченное решение называется локализованным (в строгом смысле), если множество L = {x RN : |u(T0, x) := limtT0 u(t, x) 0} ограничено в RN.
Неограниченное решение называется эффективно локализованным, если мно жество L = {x RN : |u(T0, x) = } ограничено.
Уравнение (3) допускает автомодельное решение (а.р.) [16]:
ua (t, x) = (t)() = (1 t/T0 ) 1 (), (5) m = x/(t) = x (1 t/T0 ) 1, m = ( 1)/2. (6) Оно соответствует начальным данным ua (0, x) = (x). Функцию () 0, опреде ляющую пространственно-временную структуру автомодельного решения, согласно принятой терминологии [18] называют собственной функцией (с.ф.) горения нелинейной среды, описываемой уравнением (3). Далее для краткости исполь зуется название собственная функция. Она удовлетворяет в R+ вырождающемуся ОДУ:
1 1 N 1 ( N 1 ) + = 0 (7) + 2( 1)T0 ( 1)T и краевым условиям:
(0) = 0, () = 0, = 0, если = 0. (8) Уравнение (7) имеет два решения, являющимися константами: () 0 и () H = (T0 ( 1))1/(1). Они играют существенную роль при анализе различных ре шений уравнения (7). Если задача (7), (8) имеет решение для некоторого T0 = T0,1, то она имеет решение и при любом другом T0 = T0,2, при этом решения связаны преобразованием подобия [19]. Это дает возможность без ограничения общности по ложить T0 = 1/( 1), и тем самым, H 1.
Тогда уравнение (7) принимает вид:
L() 1N ( N 1 ) + m + = 0. (9) Анализ решений задачи (9), (8), (см. [17], [16], Глава IV, и ссылки в них), дает сле дующие результаты:
При любых 1 + 1 существует финитное решение () 0.
При + 1, N 1 и = + 1, N 1 задача не имеет немонотонных решений. Единственность решения доказана для случая + 1, N = 1.
При + 1, N 1 задача не имеет финитных решений.
Если + 1 s = ( + 1)(N + 2)/(N 2)+, то задача имеет по крайней мере одно решение () 0 в R+, которое строго монотонно убывает по и имеет асимптотику () = Ca 2/(1) [1 + ()], (10) () 0,, Ca = Ca (,, N ) – постоянная. Позже [27] интервал по был расширен.
При N = 1, + 1 задача имеет не менее чем K = [a] 1, [.] целая часть, a = ( 1)/( 1) 1 (11) решений, различающихся по числу экстремумов на полуоси [0, ) [20]. Обо значим их через a,i (), i = 1, 2,... K. На основе линейного анализа и некоторых численных результатов в работе [21] было высказано предположение, что число раз личных решений a,i () при + 1 и N 1 равняется K + 1. В случае N = этот результат был уточнен недавно [14] бифуркационным анализом: число реше ний дается формулой K = [a], если a – нецелое, и K = a 1, если a – целое.
Для N = 2, 3 бифуркационный анализ дает ту же оценку числа решений, но при + 1, + 1 она нарушается (см. раздел 2.2 и [14]).
На основе этих результатов устанавливаются базисные режимы горения сре ды, описываемыми а.р. (5),(6). Для их характеризации вводятся также следующие понятия: полуширина xs = xs (t), которая для монотонных по x решений с един ственным максимумом в точке x = 0 определяется уравнением u(t, xs ) = u(t, 0)/2, и точка фронта xf : u(t, xf ) = 0, u ux (t, xf ) = 0.
HS -режим с обострением, 1 + 1. Формируется тепловая волна, которая за время T0 охватывает все пространство. Процесс с обострением не лока лизован: mes L = mes L =, xs, xf, t T0.
S -режим с обострением, = +1. Диффузия тепла и интенсивность нагрева согласуются так, что приводят к локализацию процесса в области диаметром Ls, на зываемой “фундаментальной длиной” S -режима: L = L = {|x| Ls /2}. Внутри L среда нагревается до бесконечной температуры за время T0 ;
xs = const, xf = Ls /2.
LS -режим с обострением, + 1 f = + 1 + N. Интенсивность источника сильнее, чем диффузия. Среда нагревается до бесконечной температуры за время T0 только в одной точке: mes L = 0, xs 0, t T0. В соответствии с различными решениями a,i (), i = 1, 2,..., среда горит в виде простых (i = 1) и сложных структур (i 1) с одним и тем же моментом обострения.
Чтобы показать значимость автомодельных решений как аттракторов широких классов других решений того же уравнения, введем еще несколько понятий.
В случае произвольных финитных начальных данных u0 (x) (4) вводится в рас смотрении автомодельное представление решения u(t, x) задачи (3), (4), определя емое в каждый момент времени в соответствии с видом а.р. (5), (6):
m = 1 (t)u(t, (t)) (12) (t, ) = (1 t/T0 ) 1 u t, (1 t/T0 ) Автомодельное решение ua (t, x) называется асимптотически устойчивым, ес ли существует достаточно широкий класс решений u(t, x) задачи (3), (4) с началь ными данными u0 (x), автомодельные представления (t, ) которых стремятся в некоторой норме к () когда t T0 :
(13) (t, ) () 0, t T0.
В определении автомодельного представления (12) участвует неизвестное вре мя обострения T0 и оно непригодно для численного исследования асимптотической устойчивости а.р. В работе [19] был предложен и численно реализован (при N = 1 ) другой подход, позволяющий исследовать “структурную” устойчивость неогра ниченных а.р. в специальной “автомодельной” норме, согласованной при каж дом t с геометрической формой структуры и не использующей в явном виде време ни обострения T0. Вводится новое автомодельное представление в соответствии со структурой а.р. ua (t, x) :
maxx u(t, x) (t, ) = u(t, ((t))m )/(t), (14) (t) =.
max () Автомодельное решение ua (t, x) называется структурно устойчивым, если сходи мость (13) имеет места для (t, ), заданной через (14).
Для автомодельных решений со сложной пространственно-временной структу рой, отвечающих немонотонным собственным функциям горения нелинейной сре ды, вводится понятие метаустойчивости. Автомодельное решение ua (t, x) называет ся метаустойчивым, если для каждого 0 существует класс начальных данных u0 (x) (x) и время T, T0 T T0 таких, что для автомодельных представлений соответствующих решений выполнялось (15) (t, ) () C(), 0 t T.
Глава 1 посвящена численным методам решения основных задач – нелинейных автомодельных задач (раздел 1.1) и нелинейных параболических задач с сингуляр ными по времени (blow-up) решениями (раздел 1.2). Детальное описание и иссле дование методов сделано для случая радиальной симметрии.
Общими трудностями как для автомодельных, так и для параболических за дач, являются: нелинейность;
зависимость от нескольких параметров (как правило – трех);
недостаточная гладкость решений на поверхностях вырождения эллиптиче ского оператора, где решения обращаются в ноль. В случае радиальной симметрии при N 1 и в полярных координатах в двумерном случае к ним добавляется син гулярность в точке x = 0 ( = 0).
Что касается автомодельных задач, то самая существенная трудность – это неединственность решения, что существенно отличает их от классических задач с единственным решением. Возникают следующие проблемы:
найти “хорошее” приближение к каждому решению;
сконструировать итерационный метод, который: сходится всегда к искомому решению (соответствующему начальному приближению);
сходится быстро;
обеспе чивает достаточную точность;
автоматизировать процесс вычислений так, чтобы единообразно и быстро нахо дить все различные решения при данных параметрах задачи (,, N );
определить априори куда перенести условия из бесконечности, чтобы асимпто тика (10) выполнялась (отметим, что чем больше номер с.ф., тем больше константа Ca в асимптотике).
Для преодоления указанных проблем здесь использованы:
– подход “линеаризации” [19] автомодельного уравнения относительно нетри виального пространственно-однородного решения H и сшивание решений линеа ризованного уравнения с известной асимптотикой, для нахождения приближений к различным решениям;
– итерационные схемы, основанные на непрерывном аналоге метода Ньюто на (НАМН) для решения нелинейных функциональных уравнений [35], который имеет большую, по сравнением с классическим методом Ньютона, область сходимо сти;
– метод конечных элементов (МКЭ) разного порядка точности на квази равномерных сетках, для решения линейного уравнения НАМН на каждом шаге итерационного процесса;
– МКЭ, основанный на несимметричном методе Галеркина [22], для сингу лярных в центре симметрии задач при N 3.
В разделе 1.1.1 поставлена радиально-симметричная автомодельная задача:
Найти функцию () 0, удовлетворяющую уравнению (9) и краевым услови ям (8), которые с учетом перечисленных ее свойств и асимптотику (10), уточнены следующим образом:
для + 1, (16) (0) = 0, (l) = 0, где l выбирается так, чтобы краевое условие не влияло на решение, т.е., l было больше носителя 0 решения;
для + 1, (17) (0) = 0, (l) + (p/l)(l) = 0, p = 1/m = 2/( 1), где l выбирается так, чтобы с хорошей точностью выполнялась асимптотика (10).
Детальное изложение численного метода сделано на примере этой задачи. В нью тоновской итерационной схеме (раздел 1.1.2) (18) L (k )vk = L(k ), k+1 = k + k vk, 0 k 1, k = 0, 1,..., k = k () = (, tk ), vk = vk () = v(, tk ), (19) 0 () – начальное приближение, уравнение (18) имеет вид 1 N 1 k vk N 1 k k vk + mvk + N 1 N (20) 1 N 1 k k +(1 k )vk = + mk + k k.
N Итерационные поправки vk () удовлетворяют однородным граничным условиям vk (l) = 0 для + 1, (21) vk (0) = 0, vk (l) + (p/l)vk (l) = 0 для + 1. (22) vk (0) = 0, В разделе 1.1.3 проведена дискретизация по МКЭ задач (20), (21) и (20), (22) в слабой форме:
Найти функцию vk () H 1 (0, l), удовлетворяющую интегральному тождеству l N 1 v()w()d (23) (L (k )vk, w) = (L(k ), w), w H (0, l), (v, w) = и граничным условиям (21) или (22), где k () – заданная функция, k D = {() : +1, d+1 /d L2 (0, l)}, H (0, l) = {w : (N 1)/2 w, (N 1)/2 w L2 (0, l), (1 )w(l) = 0}.
Значение = 0 отвечает условию (21), а = 1 – условию (22).
Аппроксимация интегрального тождества (23) проводится на основе МКЭ с ис пользованием Лагранжевых (линейных и квадратичных) конечных элементов. Так для вектора V с компонентами значения итерационных поправок в узлах конечно элементной сетки получена линейная система уравнений (24) A()V = B() с несимметричными ленточными матрицами A и B. Для решения системы (24) использовано LU разложение матрицы A, которое не увеличивает ее профиля.
Шаг k в (19) определяется экстраполяционной формулой [23]:
k1 k, если k k1,, если k k1, k = min 1, k1 k = max 0, k k k где k = max h |B(k )k |. Этот выбор k при уменьшении k обеспечивает стрем ление k к 1, и скорость сходимости (18), (19) становится квадратичной.
Выбор начальных приближений обсуждается в Разделе 1.1.5. Для + предложено (N ) (N ) A(N ) cos2 (/20 ), 0, N (25) 0 () = (N ) 0, 0 l, где A(N ) 1 – константы. Приближения (25) основаны на известном [5] точном решении для S -режима при N = 1 :
2( + 1) (1) (1) cos2 (1), 0 = LS /2 = + 1/, (26) () =, ( + 2) (1) 0, на полученных в [21] приближенных оценках для 0 при N = 2, 3 :
(N ) (N ) LS /2 = 0 + 1 + 1.84(N 1)//, (N ) (N ) и на установленном нами в численном эксперименте соотношении 0,HS 0.
Для + 1 “линейные приближения” искались в виде:
1 + i y(), 0 i (27) 0,i () =, Ci p, i l.
где y() = F (( 1)/( 1), N/2;
( 1) 2 /4) – вырожденная гипергеомет рическая функция. Она является решением линеаризованного относительно H автомодельного уравнения ( () = 1 + y(), | y() | 1 ):
1N ( N 1 y ) + m y ( 1)y = 0. (28) Параметры i, i, Ci определялись из условий 0,i () C 2 [0, l] ( C 2 – приближения) 1 или 0,i () C [0, l] ( C – приближения). В обоих случаях Ci = ip+1 y (i )/[p y(i ) + i y (i )]. (29) i = p/[p y(i ) + i y (i )], Для C 1 -приближений точка i является корнем уравнения y() = 0 (всего K + корней), а для C -приближений – корнем уравнения (p + 1)y () + y () = 0.
На основе изложенной методики создан комплекс программ на языке F ORT RAN.
Он дает возможность вычислять решения автомодельной задачи для всех видов ре жимов с обострением. В случае LS -режима вычисление “линейных приближений” и сшивание с асимптотикой реализовано программно, так что процесс полностью ав томатизирован – задается только номер искомой с.ф.. Это сделало возможным тща тельное исследование сходимости и точности численного метода (Раздел 1.1.6). Оно проводилось посредством сравнения с известным точным решением (26) S режима при N = 1 и методом Рунге на последовательностях вложенных сеток для HS и LS режимов. Исследование показало: в той части области, где решение имеет необходимую гладкость, при использовании линейных элементов имеет место оп тимальный порядок точности O(h2 ), a при использовании квадратичных эле ментов – сверхсходимость порядка O(h4 ) в узлах разбиения.
Многочисленные эксперименты, проведенные в широком диапазоне изменения значений параметров и при N = 1, 2, 3, указывают на быструю сходимость итерационного процесса (18), (19). Замечены однако следующие особенности:
– при + 1 + 0 гипотеза о малых колебаниях с.ф. около гомотермического решения, на которой подход линеаризации основан, не выполняется. При N = 2, были обнаружены с.ф., значения которых равны нулю в кольцевой области около центра симметрии. Это обусловило необходимость детального исследования случая + 1 + 0.
– при N = 2, 3 в сингулярной точке = 0 порядок точности ниже ( 1.7 для линейных и 3.7 для квадратичных элементов).
Для преодоления потери точности в окрестности = 0 была использована моди фикация метода Галеркина на сингулярный случай. Она предложена и теоретически исследована для линейных сингулярных при x = 0 задач в работе [22] и названа несимметричным методом Галеркина. Идея несимметричного метода – приме нить метод Галеркина к специальной, несимметричной форме оригинальной (самосо пряженной) задачи. Применительно к автомодельной задаче он изложен в Разделе 1.1.7. Вычислительная схема построена на основе альтернативной формы запи си уравнения (9):
L() = ( ) (N 2) + m 1+ + (1 1 ) = 0 (30) где = 0 при N = 1, = 1 при N = 2, 3. При N = 1, 2 она совпадает с (9).
Многочисленные эксперименты на последовательностях вложенных сеток, прове денные для S, HS и LS режимов, N 3, показывают оптимальный порядок точности этой вычислительной схемы вплоть до центра симметрии.
Раздел 1.2 посвящен численным методам решения параболических задач.
Основная трудность численного решения параболических задач - неограничен ность решений, при этом неограниченность в изолированных точках, в конечной об ласти или во всем пространстве. И еще связанные с ней: подвижная граница, на которой решение не во всех случаях достаточно гладкое;
неустойчивость неогра ниченных решений.
Неограниченные решения режимов с обострением неустойчивы в том смысле, что малым возмущениям начальных данных соответствует малое (априорно неизвестное) изменение времени обострения, и следовательно – сколь угодно большие изменения в решении вблизи момента обострения. Именно для устранения проблем с неустой чивостью вблизи момента обострения введены понятия структурной устойчивости и метаустойчивости. Они отражают сущность режимов с обострением. Простые струк туры устойчивы в этом смысле вплоть до момента обострения, сложные структуры вырождаются в простые вблизи момента обострения. При этом вырождение может происходить разным образом - сложная структура может выродиться в одну или несколько простых. Эта неустойчивость стохастическая.
Одна из задач этой работы – исследование структурной устойчивости и мета устойчивости различных инвариантных решений – накладывает существенные тре бования на методы решения параболических задач – реализовать методы, вычис лительная неустойчивость которых существенно “меньше” неустойчивости решения.
Как показано дальше, разработанные методы дают возможность, во-первых, вычис лять автомодельные решения вплоть до момента обострения, когда максимум реше ния нарастает на 3-15 порядков (в зависимости от параметров задачи) с сохранением автомодельности, и во-вторых, восстанавливать время обострения T0, заложенное в автомодельной задаче (T0 = 1/( 1)), с точностью до 105 107. Это дает возможность прогнозировать с большой точностью возникновение критических ре жимов в среде.
Для преодоления указанных проблем и для удовлетворения поставленных требо ваний, здесь использованы:
– трансформация Кирхгофа для линеаризации эллиптической части диффе ренциального оператора;
– МКЭ для дискретизации по пространственным переменным с интерполирова нием нелинейных коэффициентов по базису конечномерного пространства и концен трацией матрицы массы;
– динамические адаптивные сетки в LS -режиме (сгущающиеся сетки) и в HS -режиме (разреживающиеся сетки с сохранением числа узлов), согласованные с автомодельными (и приближенными автомодельными в полулинейном случае) зако нами;
– несимметричный метод Галеркина в радиально-симметричном случае при размерностей пространства N 3 ;
– явный метод типа Рунге-Кутта второго порядка точности с расширенной областью устойчивости и с автоматическим выбором шага, обеспечивающим выпол нение слабого принципа максимума и, в случае глобальных решений, достижение заданной точности в конце интервала интегрирования.
В разделе 1.2.1 поставлена радиально-симметричная параболическая задача:
(rN 1 u ur )r + u, (31) ut = 0 r R, 0 t T0, rN (32) ur (t, 0) = 0, 0 t T0, u(t, R) = 0, 0 t T0, (33) u(0, r) = u0 (r) 0, 0 r R, Величина R выбирается так, чтобы граничное условие в (32) не влияло на решение u(t, r). При + 1 это возможно из-за локализации процесса, при + величина R = R(t) меняется автоматически в процессе счета.
Симметричный метод Галеркина изложен в разделе 1.2.2. Применение транс формации Кирхгофа u u+ w dw = (34) p(u) =, + является существенным для последующего интерполирования нелинейных коэффи циентов и оптимизации процесса вычислений.
Задача (31)-(34) в слабой форме формулируется так:
Найти функцию u(t, r) D, D = {u : u, u+1 /r L2 }, которая удовлетворяет интегральному тождеству (35) (ut, v) = A(t;
u, v), v H0 (0, R), 0 t T0, начальному условию (33) и граничным условиям (32). Здесь R R p(u) v N rN 1 q(u) = u, (u, v) = r u(r)v(r)dr, A(t;
u, v) = + q(u)v dr, r r 0 H0 (0, R) = {v : r (N 1)/2 v, r(N 1)/2 v L2 (0, R), v(R) = 0}.
Для дискретизации (35), (33) использованы Лагранжевы (линейные и квадратич ные) конечные элементы.
Пусть {i }n – стандартный узловой базис Лагранжа пространства Sh. Интер i= полируем u(t, r), u0 (r) и нелинейные функции p(u), q(u) по базису Sh :
n n (36) uh (t, r) = ui (t)i (r), u0 (r) = u0 (ri )i (r), i=1 i= n n (37) p(u) pI = p(ui )i (r), q(u) qI = q(ui )i (r).
i=1 i= Использование аппроксимаций (36),(37) и концентрация матрицы массы дают сле дующую систему ОДУ (записанной в нормальной форме):
u = M1 Kp(u) + q(u), (38) u(0) = u0.
Здесь u = [u1 (t),..., un (t)]T, p(u) = [p(u1 ),..., p(un )]T, q(u) = [q(u1 ),..., q(un )]T, n R R rN 1 i j dr, K = {kij }, kij = rN 1 i j dr.
M = diag{mii }, mii = mij, mij = 0 j= Отметим следующие два преимущества этого подхода: матрица жесткости K по стоянная;
матрица массы M диагональная, так что произведение M1 K имеет ту же разреженную структуру, как и матрица K.
Несимметричный метод Галеркина обобщен для нелинейной параболической за дачи (31)–(33) в разделе 1.2.3. После трансформации Кирхгофа (34) уравнение (31) записывается в следующей несимметричной форме:
r ut = (r pr )r + (N 2)pr + r q(u), = 0, N = 1;
= 1, N 1.
на основе которой делается дискретизация. Применяя все следующие шаги как в симметричном методе, получаем систему ОДУ:
u = M1 (Kp(u) + (N 2)Bu) + q(u), (39) u(0) = u0, R R R r i j dr, r i j dr, B = {bij }n, r i j dr.
mij = kij = bij = i,j= 0 0 Матрица B несимметричная.
Модификация явного метода Рунге-Куттa второго порядка точности с расширен ной областью устойчивости, использованная для решения систем ОДУ (38) и (39), приведена в разделе 1.2.4. Реализован автоматический выбор шага так, чтобы од новременно выполнялись условия устойчивости и достижения заданной точности в конце интервала интегрирования. Для сингулярных по времени решений исполь зован стоп-критерий 1016 и тогда T0 – время, достигнутое в счете.
Выбор явного метода решения системы ОДУ обусловлен следующими причинами:
– исследование чисто неявной нелинейной разностной схемы для решения задачи (31)-(33) при N = 1, проведенное в [16], §5, Глава V II, показало, что ее одно значная разрешимость на временном слое накладывает зависимость между шагами по пространству и по времени такова порядка, как и условия устойчивости явной схемы;
– вблизи момента обострения, когда значения решения на временном слое разли чаются на несколько порядков, связывать их в системе уравнений необоснованно;
– при одинаковых ограничениях на шаги по времени, явные методы намного эко номичнее неявных;
они допускают и распараллеливание;
– явные методы лучше приспособлены к локальному сгущению пространствен ной сетки, которое абсолютно необходимо в случае обостряющихся в изолированных точках решений.
Идеология динамических адаптивных сеток, согласованных с автомодельными (и приближенными автомодельными) законами, изложена в самом общем случае в разделе 1.2.5.
Если дифференциальное уравнение ut = Lu допускает автомодельное решение ua (t, x) = (t)(), = x/(t), то связь между x и дает идею как адаптировать пространственную сетку. Для задачи (31)-(33) эта связь дает:
maxr u(t, r) = r(t)m, = r(t)m, m = ( 1)/2. (40) (t) =, maxr u0 (r) На основе соотношений (40) выработана следующая стратегия адаптации. Пусть r(k) – шаг по пространству (длина конечного элемента) в момент времени t = tk.
В LS режиме, m 0 шаг r (k) выбирается так, чтобы шаг (k) был огра ничен сверху: (k) = r(k) (t)m h0. Когда условие нарушается, делается следу ющее:
– каждый элемент делится пополам и значения решения в новых точках находятся линейной или квадратичной интерполяцией по старым значениям;
– чтобы не увеличивать излишне число точек, сгущение сетки делается только в той части области, где решение меняется;
точки сетки, где оно уже установилось с заданной точностью u (в расчетах обычно u = 107 ), отбрасываются.
Отметим особо, что условие для (k) проверяется на каждом временном шаге, но сетка сгущается только тогда, когда оно нарушается. Тогда производится и проверка для отбрасывания точек сетки из-за установления решения.
В HS режиме, m 0 шаг r (k) выбирается так, чтобы шаг (k) был огра ничен снизу: h0 / (k) = r(k) (t)m. Когда условие нарушается, делается следу ющее:
– точки сетки, начиная с второй, отбрасываются через одну;
– длина R(t) интервала интегрирования увеличивается вдвое добавлением но вых элементов: R(tk ) = 2R(tk1 ). Адаптация сетки в этом случае производится с сохранением числа элементов.
Эта идея – вложение структурных свойств (геометрии, разных видов симмет рий, законов сохранения) непрерывных задач при разработке дискретных методов – лежит в основе нового важного направления в вычислительной математике – гео метрическое интегрирование, которому посвящено много работ и монографий (см.
например [25],[26]). Достоинство предложенного подхода по сравнению с применяе мыми обычно методами подвижных сеток [26] в том, что для сетки здесь не решается дополнительная дифференциальная задача.
На основе изложенного метода создан комплекс программ на языке F ORT RAN для решения задачи (31)-(33) на множестве всех допустимых физических парамет ров. Создана также компьютерная система для визуализации особенностей нелиней ных процессов. Комплекс программ и система визуализации используются в процессе обучения по курсам “Нелинейные математические модели” и “Математические моде ли и вычислительный эксперимент” в ФМИ Софийского университета.
В разделе 1.2.6 приведены результаты численного исследования сходимости и точности метода решения задачи (31)–(33) как с автомодельными, так и с неавтомо дельными начальными данными.
Пример 1. Это случай = + 1, N = 1, когда известно точное решение (с.ф.) (26) автомодельной задачи. Оно использовано в качестве начальных данных для параболической задачи при разных вычислительных параметрах и анализиро вались его автомодельные представления. Как было сказано, для обостряющихся решений сходимость и оценка ошибки в автомодельной норме имеет смысл вплоть до момента обострения. Зависимость этой ошибки от дискретизации по про странству и по времени показана в Таблице 1. Там содержатся также: вычисленное время обострения T0, число шагов по времени, максимум решения за это время, и относительные ошибки e(t, ) = |(t, ) ()|/max () автомодельных представ лений решения ua,h (t, r) для t = T0 в трех точках по r. Как показывают эти и все проведенные эксперименты, восстановление (в процессе решения параболической за дачи) времени обострения, заложенном в автомодельной задаче, как и сохранение автомодельности, зависят больше всего от дискретизации по пространству.
Табл. 1. S режим: = 2, = 3, N = 1, T0 = 0. h T0 N ST EP u(T0, 0) e(T0.1.8) e(T0, 2.7) e(T0, 2.8) 0.1 102 0.507910 1262 0.13E+8 0.11E - 2 0.33E - 2 0.18E - 0.05 102 0.507077 5266 0.62E+7 0.28E - 3 0.44E - 3 0.11E - 0.05 104 0.499976 5870 0.54E+7 0.28E - 3 0.44E - 3 0.11E - 0.05 106 0.499958 9057 0.50E+7 0.28E - 3 0.44E - 3 0.11E - 0.025 104 0.499998 21737 0.31E+7 0.70E - 4 0.11E - 3 0. Пример 2. В Таблице 2 показаны некоторые результаты численного исследова ния эволюции второй с.ф. 2 () при = 2, = 3.25, N = 2, T0 = 1/( 1) = 0.(4).
Они показывают, что и для немонотонных с.ф. при уменьшении шага h время t a, до которого автомодельность сохраняется, увеличивается.
Табл. 2 LS режим, 2 () : = 2, = 3.25, N = 2, T0 = 0.(4).
h T0 ua,m ta ta / T 0.4414948 99.85405% 0.25 0.4421401 102 0.4438662 99.99799% 0.125 0. 103 0.4443100 99.99996% 0.0625 0. 0.03125 0.4444104 2 104 0.4444104 99.99999% Глава 2 посвящена исследованию автомодельных решений и их структурной устойчивости в радиально-симметричном случае.
В разделе 2.1 проведено исследование “линейных приближений” к собственным функциям при + 1 + 0. Получены асимптотики решений yN () = F (a, N/2, ( 1) 2 /4), a = ( 1)/( 1) линеаризованного уравнения (28) ( N – размерность пространства):
1 1 1 Y1 () = lim y1 () = lim F a, ;
= cos( ) 2a a a 1 1 Y2 () = lim y2 () = lim F a, 1.0;
= J0 ( ) a a a 3 1 1 2 sin( ) = Y3 () = lim y3 () = lim F a, ;
.
2a 4 a a где J0 (x) – функция Бесселя. Отсюда для расстояний между C 1 -точками сшивания i (27) (а в пределе + 1 + 0 и между C 2 -точками сшивания) получено i+1 i /, + 1 + 0.
Исследовано также поведение точек пересечения с.ф. i с пространственно-однород ным решением () 1 и поведение констант i в (27),(29) при + 1 + 0. На основе результатов даны рекомендации по использованию “линейных приближений”.
В разделе 2.2 исследована и уточнена структура с.ф. LS -режима при переходе к S -режиму.
1. Показано, что структура с.ф. при N 1 и +1, +1 существенно от личается от структуры с.ф. при N = 1. При N = 1 максимумы с.ф. увеличиваются с удалением от центра симметрии [20], a при N 1 вычислительным эксперимен том было обнаружено, что при + 1 + 0 первый, центральный максимум с.ф.
становится больше других (такие с.ф. были названы потом [12] с.ф. “с перетяжкой”).
2. Различен и переход LS -режима к S -режиму. При N = 1 переход +1+ “непрерывен” – с.ф. k () LS -режима стремится к с.ф. S -режима, составленной из k элементарных решений (26) (Рис. 1а). При N 1 такого “непрерывного” перехода нет.
2.1. Для фиксированного и некоторого j = j (, N ) центральный (первый) (N ) минимум с.ф. с четными номерами 2j, j = 1, 2,..., становится равным нулю и N для j (, N ) все с.ф. 2j () аннулируются в некоторой области около начала координат (Рис. 1б):
(N ) (N ) (N ) (2j ) (2j ) (j ) = 0, 2j () = 0, 0 j = j (,, N ), j = 1, 2,..., Рис. 1. Поведение с.ф. при + 1 + (N ) при этом j растут при +1+0. Все максимумы 2j () стремятся к максимуму (N ) с.ф. S -режима при тем же и N = 1. Таким образом с.ф. 2j (), “уходя на бесконечность” при + 1 + 0, стремится к с.ф. S -режима при N = 1 и того же, состоящей из j элементарных решений (26) (Рис. 1е).
2.2. Для фиксированного существует такое значение j = j (, N ), что при (N ) + 1 j с.ф. 2j+1 (), j = 1, 2,... распадается на две части: центральная, которая стремится к с.ф. S -режима для соответствующего значения N ( (N ) () ), (N ) и вторая, которая совпадает с с.ф. 2j (), “уходящей на бесконечность” при + 1 + 0. Для фиксированных, N всегда j j (Рис. 1в-1д).
Согласно описанному “сценарию” перехода LS к S -режиму, при + 1 + остается только первая с.ф. LS -режима, которая стремится к единственной с.ф.
S -режима. Этот “сценарий” является результатом детального вычислительного экс перимента.
Так в процессе этого исследования найдены с.ф. новой структуры – с.ф. “с пе ретяжкой” и с.ф. с левым фронтом.
Существование с.ф. с левым фронтом было подтверждено асимптотическим ана лизом – получена асимптотика с.ф. в окрестности левого фронта + 0 :
() C( ), = 1/, C = (( 1) /2)1/.
3. Собственные функции с левым фронтом представляют собой тела вращения с “пустотой” в центре. Численные эксперименты показывают, что при + 1 + (N ) радиусы “пустот” увеличиваются и четные с.ф. 2j, j = 1, 2,... “выстраиваются” так, что разности этих радиусов стремятся к одной и той же величине. Гипотеза автора о ней: = () = 2 + 1/ = LS / (Рис. 1ж – 1и).
Полученные новые типы решений инициировали исследования других авторов [12]-[14] другими методами (методом динамической аналогии, бифуркационным ана лизом). Их исследования подтвердили существование этих новых решений, и в част ности [14] описанного выше сценария.
В разделе 2.2.3 исследован аналитически предельный случай 0, =+1.
Показано, что при 0+ автомодельное решение S-режима стремится к авто модельному решению с конечной энергией соответствующего линейного уравнения ut = uxx + u, x R+, t 0.
В разделе 2.3 исследовано существование с.ф. и структурная устойчивость а.р.
LS -режима при переходе параметра через некоторые пороговые значения, полу чившие название “критические экспоненты” (см. [27] и ссылки там), когда свойства решений как автомодельной, так и параболической задачи, меняются существенно.
До сих пор обнаружены следующие критические экспоненты:
f = + 1 + 2/N (экспонента Фуджита);
s = ( + 1)(N + 2)/(N 2), N 3 (экспонента Соболева);
u = ( + 1)(1 + 4/(N 4 2 N 1)), N 11;
u =, N 11 ;
3( + 1) + ( 2 (N 10)2 + 2(5 + 1)(N 10) + 9( + 1)2 )1/ p = 1 +, N 11.
N В разделе 2.3.1 исследовано асимптотическое поведение неограниченных реше ний с конечной энергией для параметров f = + 1 + 2/N, когда автомодельное решение (5) ua (t, r) L1 (RN ). В этом случае качественная теория нестационарно го осреднения “амплитуда – полуширина” предсказывает автомодельное поведение амплитуды решения задачи (31)–(33) и возможное неавтомодельное поведение полу ширины [16]. Там был поставлен вопрос: какое инвариантное или ПАР описывает асимптотическую (t T0 ) стадию процесса? Численное исследование стало воз можным благодаря адаптации сетки. Показано, что а.р. (5), соответствующее 1 (), структурно устойчиво: во всех экспериментах с финитными начальными данными, обеспечивающими неограниченность решений, на асимптотической стадии они стре мятся к автомодельному.
В разделе 2.3.2 конструированы численно решения а.з. при соотношениях па раметров за критическими экспонентами s, N 3;
u N 11;
p, N 11. Таким образом конструктивно установлено существование ре шений а.з. при u N 11, что является открытой проблемой [27]. Показано тоже, что соответствующие им автомодельные решения структурно устойчивые, что подтверждает гипотезу из [27].
В разделе разработанные для радиально-симметричных задач (9), (16),(17) и (31)-(33) методы обобщены на случай системы нелинейных уравнений с источниками, которая описывает процессы теплопроводности и горения двухкомпонентной среды:
u1t = x1N (xN 1 u1 u1x )x + u1 u2, x R+, N = 1, 2, 3, 1 (41) u2t = x1N (xN 1 u2 u2x )x + u1 u2, i 0, i 1, i 0, i = 1, 2.
2 Система допускает неограниченные автомодельные решения ua = (u1a, u2a ) u1a = g1 (t)1 () = (1 t/T0 )m1 1 (), = x/(1 t/T0 )n, (42) u2a = g2 (t)2 () = (1 t/T0 )m2 2 (), mi 0, i = 1, 2, mi = i /p, i = i + 1 i, i = 1, 2, p = (1 1)(2 1) 1 2, n = (m1 1 +1)/2 = (m2 2 +1)/2, 1 (2 +12 ) = 2 (1 +11 ), 0 N n mi, i = 1, 2.
Функции 1 (), 2 () удовлетворяют системе автомодельных уравнений L1 (1, 2 ) 1N ( N 1 1 1 1 ) + n1 m1 1 1 1 2 2 = 0, (43) L2 (1, 2 ) 1N ( N 1 2 2 2 ) + n2 m2 2 1 1 2 2 = и краевым условиям (44) fi () = (mi /n) (fi ()/), =l 1, i = 1, В уравнениях (43) без ограничения общности положено T0 = 1. Поведение автомо дельных решений при t T0 зависит от знака параметра n : n 0 – HS режим с обострением;
n = 0 – S режим с обострением;
n 0 – LS режим с обострением.
Системы (41) и (43) называются системами с сильной обратной связью, если p 0, и системами с слабой обратной связью, если p 0.
В случае N = 1 автомодельное решение (42) и задача (43), (44) рассматривались в работах Курдюмова С.П., Куркиной Е.С., Малинецкого Г.Г. и Тельковской О.В.. При численном построении с.ф. для конкретных параметров использованы разностные методы и метод стрельбы.
Для построения с.ф. и детального исследования структурной устойчивости авто модельных решений при N = 1, 2, 3 здесь использованы обобщения методов, изло женных в Главе 1. Сформулируем основные результаты:
– структурно устойчивыми являются только а.р. с двумя компонентами простой структуры систем с сильной обратной связью;
– все другие а.р. метаустойчивы.
Таким образом показана принципиальная возможность применения разработан ных методов для исследования процессов самоорганизации в широких классах нели нейных диссипативных сред, описываемых системами типа реакция–диффузия.
В Главе 3 исследуется асимптотическое поведение сингулярных по времени ре шений трех классов полулинейных параболических уравнений с источниками ut = r1N (rN 1 ur )r + Q(u), r R+, t 0, N = 1, 2, 3 (45) Q(u) = (1 + u) ln (1 + u), 1, Q(u) = u, 1, Q(u) = eu.
Условия существования неограниченных решений таких уравнений изучаются с 1966 г. и этой проблеме посвящено большое количество работ (см. обзор в книге [28]). Случай полулинейных уравнений, где из-за недостаточности диффузии не воз никает сложных структур, оказался более трудный для исследования, чем квазили нейный случай. Причина этому – несуществование точных инвариантных решений (кроме тривиальных гомотермических) вообще или для некоторых размерностей N пространства. В этих случаях асимптотическое поведение неограниченных решений описывается приближенными автомодельными решениями (ПАР), которые являются точными а.р. других уравнений – уравнений типа Гамильтона-Якоби. Это означает, что на асимптотической стадии, при t T0, параболическое уравнение вырождается в уравнении первого порядка. Строгие результаты о вырождении на асимптотической стадии уравнения (45) получены в основном для N = 1 и при ограничениях на начальные возмущения. Поэтому численное исследование, осо бенно при N 1, является актуальной задачей.
В разделе 3.1 исследуется уравнение с источником Q(u) = (1+u)ln (1+u), кото рое не имеет нетривиальных инвариантных решений с обострением [16,17]. Это един ственное полулинейное уравнение, для которого до сих пор найдены аналоги трех режимов с обострением: HS режим: 1 2;
S режим: = 2;
LS режим:
2. В начале [16] качественными методами было показано, что асимптотическое поведение его неограниченных решений описывается автомодельными решениями v(t, r) = exp{(T0 t)1/(1) ()} 1, = r(T0 t)m/(1), m = ( 2)/2, (46) m ( )2 + = 0, () = cos2 (/2), 0, при = 2 (47) 1 уравнения первого порядка типа Гамильтона-Якоби vt = vr /(1 + v)+(1+v) ln (1+v).
Уже имеются строгие результаты [29] об асимптотической устойчивости ПАР (46).
При достаточно ограничительных предположениях [29] на начальное возмущение до казано, что при = 2, N = 1 область эффективной локализации LSL = [, ].
Для = 2 и N 1 нет такого точного результата – доказано только, что решение растет неограничено в сфере B {r }. Нами показана численно (раздел 3.1.3) структурная устойчивость ПАР (46) при практически произвольных начальных дан ных во всех трех режимах. В S режиме ( = 2 ) все численные эксперименты подтверждают сходимость автомодельных представлений к точному решению (47) для произвольных начальных данных при N = 1 и для достаточно больших при N = 2, 3 (чтобы обеспечить неограниченность решения, см.[16]). Показано также, что при = 2, N 1 неограниченное решение эффективно локализовано в сфере радиуса : L = {r }. Кроме того показано: если начальное возмущение нецентральное, достаточно далеко от центра симметрии и с достаточно большой энергией, то тепловая энергия не расплывается и горение ло кализовано эффективно в кольцевой области толщиной 2, что является новым результатом (Рис. 2).
Методом конечных элементов с кон центрацией матрицы массы, интерполи рованием нелинейных коэффициентов и динамической адаптацией сетки реша лись численно как оригинальная задача (для переменной u ), так и трансформи рованная, для U (t, r) = ln(1 + u(t, r)) (раздел 3.1.2).
В разделе 3.2 исследовано асимпто тическое поведение решений уравнения с источником Q(u) = u. При 1 1+2/N, u0 0, задача Коши имеет толь ко неограниченные решения, но их асимп тотика существенно отличается от той, которую следует ожидать по аналогии с квазилинейным случаем. Автомодельная задача, которая получается по этой ана логии, не имеет решения при 1 s = (N + 2)/(N 2)+.
Рис. 2. = 2, N = Качественный анализ [16] показал, что пространственно-временная структура u(t, x) при t T0 описывается ПАР r u (t, r) = (T0 t) 1 (), = (48), a (T0 t)| ln(T0 t)| () = { 1 + (( 1))2 /(4)} 1 (49) которое является точным а.р. уравнения vt + xvx {2(T0 t)| ln(T0 t)|}1 = v.
Строгое доказательство этого результата дано в [30] для N = 1 и при ограниче ниях на начальные данные. Для N 1 строгого доказательства нет. Численным экс периментом нами показана структурная устойчивость ПАР (48),(49) для N = 1, 2, при практически произвольных начальных данных.
В разделе 3.3 исследовано асимптотическое поведение решений уравнения с ис точником Q(u) = eu (модель воспламенения твердого горючего). Задача Коши имеет только неограниченные решения, но уравнение (45) при N = 1, 2 не имеет точных инвариантных неограниченных решений, кроме гомотермического uH = ln(T0 t).
Качественный анализ [16] показал, что асимптотическое поведение неограниченных решений описывается ПАР 2 r u (t, r), () = ln((1 + ) ) (50) = ln(T0 t) + (), = a (T0 t)| ln(T0 t)| которое является точным а.р. уравнения vt + xvx {2(T0 t)| ln(T0 t)|}1 = ev.
Строгое доказательство для N = 1 и при ограничениях на начальные данные дано в [30];
для N = 2 нет строго доказательства. Численным экспериментом на ми показана структурная устойчивость ПАР (50) для N = 1, 2 при практически произвольных начальных данных.
Отметим, что ПАР (48), (50) существенно отличаются от а.р., рассматриваемых до сих пор. Это отличие, во-первых, в зависимости автомодельной переменной от произведения (T0 t)| ln(T0 t)|, и, во-вторых, в структуре ПАР (50), которое является суммой гомотермического решения и функции (). Разработанные чис ленные методы (отметим особо адаптацию сетки, согласованную с ПАР) работают достаточно хорошо и при наличии этих особенностей и дают возможность устано вить структурную устойчивость ПАР (48), (50) на достаточно больших компактах || Глава 4 посвящена исследованию новых способов организации в двумерной нелинейной теплопроводной среде.
В разделе 4.1 построены численно инвариантные решения, которые описывают направленное распространение тепла и горения в нелинейной анизотроп ной среде со степенными коэффициентами теплопроводности ut = (u1 ux1 )x1 + (u2 ux2 )x2 + u, x = (x1, x2 ) R2, (51) 1 0, 2 0, 1.
В работе [31], (см. тоже [17]) методами группового анализа установлено, что уравне ние (51) допускает неограниченное инвариантное решение ua (t, x1, x2 ) = (1 t/T0 ) 1 (), = (1, 2 ) R2, i = xi (1 t/T0 )mi /(1) mi = ( i 1)/2, i = 1, 2. Функция () удовлетворяет нелинейному эллиптиче скому уравнению (в нем положено T0 = 1/( 1)) :
i i + = 0. (52) L() + i i i 2 i i= Характер распространения тепла определяется соотношением между величинами, (1 + 1), (2 + 1). По аналогию с радиально-симметричным случаем в каждом направлении i рассматриваем три типа соотношений:
1 i + 1 – HS -режим с обострением;
1 = i + 1 – S -режим с обострением;
i + 1 – LS -режим с обострением.
В разделе 4.1.1 поставлена краевая задача для уравнения (52), приведена ите рационная схема на основе НАМН и конструированы начальные приближения.
Для LS S ( LS HS ) режимов использовались начальные приближения (1) i (1, 2 ) = i,0 (1 )0 (2 ), i = 1, 2,..., k1, которые являются произведением одномерных приближений (27) для LS и (25) для S ( HS ) режимов.
Метод решения параболической задачи обсуждается в разделе 4.1.2. И здесь основной момент при реализации МКЭ – это трансформации Кирхгофа:
u () =, = 1, 2.
p (u) = () d, Они, вместе с последующей интерполяцией нелинейных функций p (u) и q(u) = u позволяют перенести нелинейность задачи на векторы p U и q U в системе ОДУ dU = M 1 K1 p1 U + K2 p2 U + q (U ), (53) U (0) = U dt относительно U (t) = {U1 (t),..., Un (t)}T, что очень важно для эффективности мето да в двумерном случае.
Метод реализован в виде комплекса программ, который представляет собой раз витие Пакета прикладных программ TERMO для решения тепловых задач методом конечных элементов, созданного в ИМИ БАН и предназначенного для решения дву мерных стационарных и нестационарных задач теплопроводности (электро- и маг нитостатики, диффузии, фильтрации) для кусочно-неоднородных сред.
Рис. 3. HS S режим: 1 = 3, 2 = 2, = 3.
Рис. 4. HS LS режим: 1 = 3, 2 = 1, = 3.
В разделе 4.1.3 приведены три типа решений, отвечающие одним и тем же на чальным данным:
1, (x1, x2 ) u(0, x1, x2 ) = u0 (x1, x2 ) =, 1 = {0 x1 b1, 0 x2 b2 }.
0, (x1, x2 ) / – HS S режим: 1 = 3, 2 = 2, = 2 + 1 = 3 (Рис. 3). По направлению x1 тепловой фронт за конечное время T0 достигает бесконечности ( S режим), а в направлении x2 тепло распространяется лишь до некоторой глубины LS, после чего происходит остановка фронта ( S режим). Носитель решения представляет выпук лую область.
– HS LS режим: 1 = 3, 2 = 1, = 3 (Рис. 4). По направлению x1 тепловой фронт за конечное время T0 достигает бесконечности ( S режим), а в направле нии x2 происходит сокращение эффективной глубины прогрева ( LS режим). На финальной стадии решение представляет собой бесконечно прогретую линию, вне которой температура конечна.
Для сравнения там приведен также S режим: 1 = 2 = 1 = 2.
В разделе 4.2 исследуются автомодельные решения двумерной модели тепловых структур в полярных координатах:
1 ut = (ru ur )r + 2 (u u ) + u, t 0, 0 r, 0 2. (54) r r Методом инвариантно-группового анализа С.Р. Свирщевским было установлено, (см.
тоже [17]), что оно допускает инвариантные решения следующего вида:
t 1 t 1 C0 t 1 m ln(1 ), (55) ua (t, r, ) = (1 ) (, ), = r(1 ), =+ T0 T0 1 T Здесь C0 – параметр семейства решений. Функция (, ) 0 удовлетворяет нели нейному эллиптическому уравнению (как и раньше, положено T0 = 1/( 1) ):
1 1 ( ) 2 ( ) + m + = 0. (56) L() C Из (55) при C0 = 0 непосредственно следует es = res = const, s = ( 1)/2C0 = m/C0. (57) Это означает, что траектории неоднородностей в среде (например, локальных макси мумов) будут логарифмическими спиралями при = + 1 или окружностями при = + 1. Направление движения для фиксированного C0, например C0 0, зави сит от соотношения и : при + 1 – к центру (скручивающиеся спирали), при + 1 – от центра (раскручивающиеся спирали).
Случай C0 = 0 рассматривался в работах [32],[33]. Там впервые конструирова ны радиально-несимметричные с.ф. горения среды в LS режиме. В HS режиме до сих пор были получены только радиально-симметричные с.ф.. Причина этому – решения автомодельной задачи всегда искались при условии lim (, ) = H 0.
Изменение идеологии – искать решения а.з., которые на бесконечности стремятся не к тривиальному решению 0, а к гомотермическому H 1 :
lim (, ) = H 1.
дало возможность получить два новых класса с.ф. в HS режиме – радиаль но несимметричные с.ф. сложной симметрии при C0 = 0 и спиральные с.ф. при C0 = 0. Изменение идеологии было подсказано в процессе тщательного исследова ния (аналитического и численного) решений линеаризованного относительно H автомодельного уравнения и их эволюции во времени.
Раздел 4.2.1 посвящен радиально-несимметричным с.ф. в LS режиме при C0 = 0. На основе метода линеаризации и сшивания с асимптотикой в рабо тах [32],[33] предложен эффективный способ определения начальных приближений ( EiM m ) для решения автомодельной задачи итерационным методом. Проведенные там эксперименты показали: структура с.ф. не всегда следует структуре начальных приближений;
итерационный процесс сходится не при всех начальных приближений.
Там был поставлен вопрос о создании более эффективного итерационного метода и это обусловило наше исследование.
Рис. 5. = 2, = 3.25. слева: E1M m, m = 2,..., 7;
справа: E2M m, m = 2,..., Используя обобщенные на случай полярных координат методы из Раздела 4. и способ получения линейных приближений из [33], нами были сконструированы с.ф. для разных значений параметров и и разной степени симметрии. В неко торых случаях однако итерационный процесс сходился к радиально-симметричным с.ф. Впоследствии в диссертации М. Колевой были предложены и реализованы на чальные приближения, число и поведение которых в области немонотонности тоже самое, как у предложенных в [33], но различающиеся от них процедурой сшивания с асимптотикой – сшивание проводилось по окружности, а не по лучам. С их помощью удалось получить радиально-несимметричные с.ф. в тех случаях, когда прежними приближениями это не удавалось. На Рис. 5 слева показаны изолинии с.ф. E1 Mm для m = 2, 3, 4, 5, 6, 7, а справа – изолинии с.ф. E2 Mm, m = 2, 3, 4, 5, 6, 7 для значений параметров = 2, = 3.25.
Отметим следующее:
– как и в работе [33], при m 1 наблюдалось раздвоение некоторых максимумов начальных приближений, при этом для данного значения m оно происходило для двух видов начальных приближений (на Рис. 5. это E2M 2 и E2M 3 );
– при одном и том же наборе параметров нам удалось получить большее количе ство с.ф., чем в работе [33].
– исследование эволюции во времени радиально-несимметричных с.ф. LS режи ма показало, что соответствующие им автомодельные решения метаустойчивые.
Они сохраняют свою структуру (в смысле (15)) до времени T, очень близкое к вре мени обострения T0. После этого сложная тепловая структура вырождается в одну или несколько простых структур, соответствующих самой простой (в общем случае радиально-симметричной) с.ф. с теми же параметрами,.
Разделы 4.2.2-4.2.7 посвящены конструированию радиально-несимметричных с.ф. в HS режиме.
Задача о численной реализации “спиральных” инвариантных решений была по ставлена в 1984 г.. Как отмечалось в [34], на пути к е решению были существенные e трудности. Во первых, считалось, что линеаризация автомодельного уравнения здесь не дает желаемого результата, потому что в полученном линейном уравнении нельзя сделать разделения переменных. Во вторых, не была известна асимптотика решений автомодельного уравнения (56) на бесконечности.
Для преодоления этих проблем сделано следующее.
1. Предложено (Раздел 4.2.3) подходящее (комплекснозначное) разделение пе ременных Yk (, ) = Rk ()eik, k = 0, k целое (для периодичности) (58) в линеаризованном относительно гомотермического решения автомодельном уравне нии 1 1 (59) (y ) + 2 y y + C0 y + ( 1)y = 0.
Для функции Rk () получено уравнение:
k 1 (60) Rk + Rk + 2 + C0 ki + 1 Rk = 2 2. Рассмотрены ограниченные в т. = 0 решения уравнения (60) при k 0 :
z = ( + C0 ki)1/2, = + 1 : Rk () = Jk (z), Jk (z) – функция Бесселя первого рода порядка k, 1 + C0 ki k 1 = +1 : Rk () = k 1 F1 (a, b;
z), a = +, b = 1+k, z =, 1 2 – вырожденная гипергеометрическая функция.
1 F1 (a, b;
z) 3. На основании свойств функций Jk (z) и 1 F1 (a, b;
z) установлено:
– можно рассматривать только случай k 0, C0 0;
– если = + 1 (и C0 k = 0) или + 1, то |Rk ()| когда ;
– если + 1, то |Rk ()| 0, когда.
4. Обосновано, что без ограничения общности можно использовать только реаль ную часть yk (, ) = |Rk (r)| cos(arg(Rk () + k)) частных решений (58) линейного уравнения (59).
5. Исследована геометрическая структура решений yk (, ) в зависимости от па раметров,, C0 и k. Самый существенный результат: линии уровня функций yk (, ) при + 1 ( + 1) и при больших близки к логарифми ческим спиралям с параметром s(s), определенным в (57). Этот факт дал возможность вывести асимптотику с.ф. (, ) из асимптотики решений yk (, ) :
yk (, ) c (1)/m cos k + (61) ln,, c = c(k,,, C0 ).
s 6. Предложены “линейные приближения” к с.ф. HS режима:
1 + yk (r, ), 0 r r0, 0 (62) k (r, ) =, 1, r r0, 0, r0 – константы. В экспериментах r0 выбиралось так, чтобы круг Vr0 (0) был телом спирали для = 0.01.
7. Чтобы установить как близки “линейные приближения” к с.ф., была исследо вана их эволюция во времени и асимптотика при t T0 (Раздел 4.2.4). Использо вано автомодельное представление (t,, ) решения u(t, r, ), соответствую щего начальным данным u0 (r, ), согласованное с автомодельным законом (55) при C0 = 0 :
(t,, ) = u(t, r, )/(t), = r(t)m, = C0 ln((t)), (63) (t) = max u(t, r, )/max u0 (r, ), = [0, ) [0, 2) Численные эксперименты, проведенные для широкого набора физических и вычис лительных параметров, показали, что при 1 приближения (62) близки с.ф.: они сохраняют свою структуру до времени, близкое к T0, при этом T0 очень близко к времени обострения а.р., отвечающих тем же параметрам,. И это вполне есте ственно – линеаризация “уничтожает” S и LS -режимы (линеаризованное уравнение имеет линейную диффузию и линейный источник, т.е. для его решений полуширина будет увеличиваться, как это имеет место в HS -режиме).
8. Асимптотика (61) предсказывает следующую асимптотику с.ф. (, ) = k (, ), k = 1, 2,... :
k (, ) 1 + (1)/m cos k + (64) ln,, = c s и такую аппроксимацию:
k (, ) = 1 + yk (, ), = const. (65) 9. Из (64) выведено граничное условие третьего рода при = l 1:
k k 1 k m = m/( 1). (66) = (m 1)/m sin k + ln, m s s При решении задачи о спиральных с.ф. возникла идея поиска в случае C0 = сложных немонотонных волн HS -режима, выходящих на гомотермический фон при. Для единообразного рассмотрения автомодельная задача ставится в обла сти = [0, l] [0, ], где = /k для C0 = 0 и = 2/k для C0 = 0, с условиями симметрии при C0 = 0 :
k k (67) (, 0) = (, ) = 0, 0 l, k с условиями периодичности при C0 = 2 k k (68) k (, 0) = k (, ), (, 0) = (, ), 0 l, k k и с условием симметрии в нуле. Задача замыкается краевым условием (66) при C0 = и краевым условием k k, m = m/( 1) = m при C0 = 0. В Разделе 4.2.5 изложен численный метод решения автомодельной задачи, в Разделе 4.2.6 – параболической задачи.
Раздел 4.2.7 содержит результаты численного исследования асимптотического поведения радиально несимметричных волн в HS режиме. Одновременно проде монстрированы надежность и хорошая точность методов решения как эллиптиче ской, так и параболической задачи.
Все проведенные численные эксперименты демонстрируют высокую точность вос становления времени обострения T0 и метаустойчивость этих волн. Эволюция во времени пяти с.ф., отвечающих одним и тем же параметрам = 2 и = 2.4, пока зана и анализирована. Первые две с.ф. отвечают параметру C0 = 0, остальные три – C0 = 1. Точное время обострения равняется T0 = 1/( 1) = 0.(714285).
Здесь приведена эволюция двух с.ф.. На Рис. 6, (t = 0), показаны поверхность и линии уровня с.ф. с одним максимумом, одним минимумом и с одной осью симмет рии, вычисленной из начального приближения (65), где k = 1, = 1. При t 0 по казаны линии уровня решения параболической задачи на четырех временных слоях.
Время обострения, найденное в счете, T0 = 0.714475. Для 0 t 0.7108 = 99.49%T движение максимума следует автомодельному закону m r(t) = r(0)(1 t/T0 ) 1 r(0)(t)m.
t=0 t=0. (, ) 10 5 0 -5 -5 -10 - -10 -5 0 5 10 -10 -5 0 5 t=0.6979 t=0.7108 t=0. 10 20 5 10 0 0 0 46 20 16 -5 -10 - -10 -20 - -10 -5 0 5 10 -20 -10 0 10 20 -120 -60 0 60 Рис. 6. Эволюция сложной волны HS режима: = 2, = 2.4, C0 = 0, k = 1.
(, ) Рис. 7. Эволюция двухрукавной спиральной волны: = 2, = 2.4, C 0 = 1, k = 2.
На асимптотической стадии фон решения замирает, сложная волна вырождается в радиально симметричную, развивающуюся на этом фоне. Автомодельные пред ставления (t,, ) решения u(t,, ) стремятся к радиально-симметричной с.ф. с конечным носителем, отвечающей тем же параметрам и. Отметим, что это про исходит при t 99.49%T0. Следовательно, сложная волна, соответствующая этой с.ф., является метаустойчивой.
Для двухрукавной спиральной волны (Рис. 7) T0 = 0.71428579 еще ближе к точ ному T0. Движение максимума и минимума следует автомодельному закону до вре мени t = 0.7065 = 98.91%T0. Вырождаясь на асимптотической стадии, волна перехо дит через стадию волны сложной симметрии, перед тем как выродится в радиально симметричную.
Отметим, что существование континуума решений а.з. в радиально-симметричном случае и HS режима, которые выходят на гомотермическое решение H при, было отмечено в [16], но этому результату не было уделено достойное внима ние. Оказывается, именно такие решения а.з. определяют найденные здесь спираль ные волны.
Вопрос о существовании спиральных структур в LS -режиме (скручивающиеся спирали) пока остается открытым. Хотя “хребты” линейных приближений к с.ф.
LS -режима при стремятся тоже к автомодельным (см. п.5. выше), но их ам плитуда стремится к бесконечности. Остается неясным с какой асимптотикой надо сшивать линейные приближения. Автор считает, что в задаче нахождения спираль ных структур в LS -режиме надо отказаться от условия ограниченности с.ф. при = 0. На эту идею наводят другие решения линеаризованного уравнения (логариф мическое решение вырожденного гипергеометрического уравнения с особенностью в нуле), а также и некоторые физические соображения.
3 ОСНОВНЫЕ РЕЗУЛЬТАТЫ ДИССЕРТАЦИИ 1. Предложена, исследована и реализована общая методика решения нелинейных ав томодельных задач, основанная на непрерывном аналоге метода Ньютона и методе конечных элементов. Она дает возможность автоматизировать процесс вычислений (включая случай неединственности решений) и единообразно решать краевые зада чи для нелинейных ОДУ второго порядка, систем таких уравнений и нелинейных эллиптических уравнений второго порядка.
Предложена, исследована и реализована общая методика решения нелинейных па раболических задач с сингулярными по времени решениями, основанная на методе конечных элементов с концентрацией матрицы массы и интерполированием нелиней ных коэффициентов. Реализованы адаптивные динамические сетки по пространству, согласованные с автомодельным или приближенным автомодельным законом.
2. В радиально-симметричном случае исследовано поведение собственных функций горения нелинейной среды при стремлении LS -режима к S -режиму, когда органи зация среды существенно усложняется. Показано существенное различие плоского и радиально-симметричного случаев. Найдены два новых типа собственных функций – “с перетяжкой” и с левым фронтом. На основе детального численного исследования предложен сценарий поведения сложных собственных функций (четных и нечетных) при стремлении к S -режиму.
Исследовано асимптотическое поведение локализованных обостряющихся реше ний параболической задачи в LS режиме, когда автомодельное решение не при надлежит L1 (RN ). Численно показано, что и в этом случае монотонное по r авто модельное решение структурно устойчиво.
Конструированы численно решения автомодельной задачи при соотношениях па раметров за критическими экспонентами. Показано, что соответствующие им моно тонные по r автомодельные решения структурно устойчивые, подтверждая таким образом гипотезу Галактионова и Вазкеса.
3. Конструированы численно собственные функции горения двухкомпонентной нели нейной среды в радиально-симметричном случае и исследована структурная устой чивость соответствующих автомодельных решений. Показано, что в LS режиме структурно устойчивыми являются только автомодельные решения простой струк туры (с монотонными по r компонентами) систем с сильной обратной связью.
4. Численно исследовано асимптотическое поведение обостряющихся решений полу линейных уравнений с различными нелинейными источниками в радиально-симме тричном случае. Подтверждены результаты качественного анализа о том, что асимп тотическая стадия эволюции обостряющихся решений описывается приближенными автомодельными решениями (точными автомодельными решениями уравнений пер вого порядка типа Гамильтона-Якоби). В случае источника (1 + u)ln2 (1 + u) показа но, что обостряющиеся решения эффективно локализованы в сфере радиуса или в кольцевой области толщиной 2.
5. Численно реализованы инвариантные решения, которые описывают направленное распространение тепла и горения в двумерной нелинейной анизотропной среде со степенными коэффициентами теплопроводности. Так конструктивным образом под тверждено их существование.
6. Расширено понятие собственной функции горения нелинейной среды введением в рассмотрение инвариантных решений, которые стремятся к нетривиальному го мотермическому решению, а не к нулю на бесконечности. Это дало возможность получить радиально-несимметричные сложные волны в HS -режиме и численно ре ализовать инвариантные решения, которые описывают распространение спиральных волн в нелинейной изотропной среде. Таким образом показано, что усложнение орга низации среды может происходить не только в LS -режиме за счет локализованных во вне и сходящихся к центру структур, но и в HS -режиме за счет метаустойчи вых расходящихся волн горения сложной структуры, в том числе раскручивающихся спиральных волн.
Публикации по теме диссертации:
[Д01] М.И. Бакирова, С.Н. Боршукова, В.А. Дородницын, С.Р. Свирщев ский. О направленном распространении тепла в нелинейной анизотроп ной среде //Препринт ИПМ АН СССР, No 182, 1985, 12 стр.
[Д02] S.N. Dimova-Borshukowa. Numerical analysis of the directed heat diusion in a nonlinear anysotropic medium. Proc. ICCM, Tokyo, 1986, VIII, 78-82.
[Д03] Бакирова М.И., Димова С.Н., Дородницын В.А., Курдюмов С.П., Са марский А.А., Свирщевский С.Р., Инвариантные решения уравнений теплопроводности, описывающие направленное распространение горе ния и спиральные волны в нелинейной среде, ДАН СССР, 299 (2), 1988, 346–350.
[Д04] С.Н. Димова, М.С. Касчиев, С.П. Курдюмов. Численный анализ одно мерных собственных функций горения нелинейной среды. I. Численный метод и эксперименты. Препринт ОИЯИ, Дубна, Р11-88-473, 1988,16 стр.
[Д05] С.Н. Димова, М.С. Касчиев, С.П. Курдюмов. Численный анализ одно мерных собственных функций горения нелинейной среды. II. Некоторые предельные случаи. Препринт ОИЯИ, Дубна, Р11-88-831, 1988, 15 стр.
[Д06] С.Н. Димова, М.С. Касчиев. Численный анализ двумерных собствен ных функций горения нелинейной анизотропной среды. Сообщения ОИ ЯИ, Дубна, Р11-88-876, 1988, 9 стр.
[Д07] С.Н. Димова, М.С. Касчиев. Численный анализ двумерных собствен ных функций горения нелинейной среды. Proc. Intern. Conf. Numerical Methods and Applications, Soa, 1988, Publ. House BASci, Soa, 1989, 121-125.