Итерационные методы решения задачи Стокса с переменной вязкостью тема автореферата и диссертации по математике, 01.01.07 ВАК РФ
Гриневич, Петр Петрович
АВТОР
|
||||
кандидата физико-математических наук
УЧЕНАЯ СТЕПЕНЬ
|
||||
Москва
МЕСТО ЗАЩИТЫ
|
||||
2011
ГОД ЗАЩИТЫ
|
|
01.01.07
КОД ВАК РФ
|
||
|
Московский государственный университет имени М.В.Ломоносова Механико-математический ф" "»ьтрт
Гриневич Петр Петрович
Итерационные методы решения задачи Стокса с переменной вязкостью
01.01.07 - Вычислительная математика
АВТОРЕФЕРАТ диссертации на соискание ученой степени кандидата физико-математических наук
0 з г.:АР
Москва 2011
4856569
Работа выполнена на кафедре вычислительной математики Механико-математического факультета Московского государственного университета им. М.В.Ломоносова
Научный руководитель: доктор физико-математических наук, профессор Ольшанский Максим Александрович Официальные оппоненты:
доктор физико-математических наук, профессор Агошков Валерий Иванович
кандидат физико-математических наук Коныиин Игорь Николаевич Ведущая организация: Московский энергетический институт (технический университет)
Защита состоится 9 марта 2011 г. в 16 часов 45 минут на заседании ' диссертационного совета Д. 501.002.16 при Московском государственном университете имени М.В. Ломоносова по адресу: Российская Федерация, 119991, Москва, ГСП-1, Ленинские горы д.1, МГУ, Главное здание, Механико-математический факультет, аудитория Ц-08.
С диссертацией можно ознакомиться в библиотеке Механико-математического факультета МГУ имени М.В. Ломоносова(14 этаж, Главное здание).
Автореферат разослан 9 февраля 2011 г.
Ученый секретарь диссертационного совета, д.ф.-м. н., профессор
Общая характеристика работы
Актуальность работы
Решение многих современных прикладных задач приводит к системам линейных алгебраических уравнений (СЛАУ) с седловой точкой. Характерной особенностью таких систем является знаконеопределенность, В симметричном случае имеются как положительные, так и отрицательные собственные значения. Важной областью, требующей решения задач с седловой точкой, является численное решение линеаризованных уравнений Навье-Стокса, описывающих течение несжимаемой вязкой жидкости. Уравнения Навье-Стокса являются основными уравнениями гидродинамики и, соответственно, играют важную роль в современной науке. В ходе их решения, как правило, возникает необходимость проводить вычисления на мелких сетках, как следствие, решаемые системы имеют большую размерность. Знаконеопределенность, большой размер и зависимость от физических параметров и параметров моделирования делает процедуру выбора метода решения таких систем нетривиальной. Решению таких систем посвящено много работ как отечественных (Г.М. Кобельков, Ю.А. Кузнецов и др.), гак и зарубежных (Р. Гловински, X. Элман, М. Бенци и др.) авторов.
Уравнения Навье-Стокса во многих случаях хорошо описывают поведение жидкостей и газов. Однако, многие вещества в природе описываются моделями с переменным коэффициентом вязкости, зависящим от внешних факторов. Примером могут служить биологические жидкости (например, кровь), нефть, зубная паста, кетчуп, крахмал, разведенный в воде и многие другие вещества. Для моделирования подобных веществ можно рассматривать модифицированные уравнения Навье-Стокса, при этом вязкость является не постоянным параметром среды, а функци-
ей, зависящей от динамических, кинематических или других характеристик среды в данной точке пространства, например, тензора скоростей деформации, давления, температуры, и т.д. Число обусловленности возникающих при дискретизации линейные систем зависит от отношения максимальной вязкости к минимальной. Данное обстоятельство предъявляет дополнительное требование к методам решения СЛАУ, а именно — независимость числа итераций от отношения максимального значения коэффициента вязкости к минимальному и от градиента коэффициента вязкости как функции пространственной переменной. Численные аспекты решения уравнений Навье-Стокса с переменной вязкостью является темой ряда современных исследований: среди них работы Жонга и др.1, Омори и Саито2, Ремана и др.3, Ольшанского и Ройскена4.
В численном анализе методов решения задачи Стокса важную роль играет условие LBB (Ладыженской-Бабушки-Брецци) и его непрерывный аналог, неравенство Нечаса. В диссертации неравенство Нечаса обобщается на случай переменной вязкости и на основе этого обобщения получены оценки эффективности предлагаемого итерационного метода. Некоторые другие обобщения неравенства Нечаса получены в работе Боровикова и Дубинского5.
lSJ. Zhong, Ы. Т. Zuber, L. Moresi, М. Gurnis. Role of temperature-dependent viscosity and surface plates in spherical shell models of mantle convection // Journal of Geophysical Research, Vol. 105, Iss. B5, 2000, pp. 11063-11082.
2/f. Ohmori, N. Saiio. On the convergence of finite element solutions to the interface problem for the Stokes system // Journal of Computational and Applied Mathematics, Vol. 198, Iss. 1, 2007, pp. 116-128.
3M. urRehman, T. Getntn, C. Vuik, G. Segal, S.P. MacLachlan. On iterative methods for the incompressible Stokes problem // International Journal for Numerical Methods in Fluids, D01:10.1002/ftd.2235.
iM.A. Olshamkii, A.Reusken. Analysis of a Stokes interface problem // Numerische Mathematik, Vol. 103, Iss. 1, 2006, pp. 12Э-149.
5 И.А. Боровиков,Ю.А. Дубинский. Некоторые разложения модулей Соболева-Клиффорда и нелинейные вариационные задачи // Труды Математического института им. В. А. Стеклова РАН, Том 260, 2008, с. 57-74.
Одной из задач, где возникают уравнения Навье-Стокса с переменной вязкостью, является моделирование среды Бингама. Среда Бингама является вязко-пластичной средой, которая при напряжениях ниже порогового значения ведет себя как твердое тело, а при превышающих пороговое значение — как вязкая жидкость. Течению среды Бингама посвящено большое количество литературы, среди отечественных работ можно отметить монографию Климова и др.6. Одним из подходов численному решению уравнений Бингама является регуляризация — модель, когда среда Бингама рассматривается как жидкость с переменной вязкостью. Задача Бингама трудно поддается математическому анализу и ее точные решения найдены только для узкого круга модельных задач, например, для течения среды между двумя параллельными пластинами и некоторых других. По этой причине численные методы, зачастую, являются единственным способом анализа многих процессов. Численными методами для решения уравнений Бингама занимаются многие исследователи, среди них Р. Гловински, М. Берковьер, М. Энгельман, Т. Папанастасио и другие.
Задачи с переменной вязкостью появляются и во многих других научных областях, например, в геологии. В мантии Земли температура неоднородна, а вязкость магмы напрямую зависит от температуры.
Цель диссертационной работы
Работа преследует следующие цели.
1. Построение эффективного метода решения систем линейных алгебраических уравнений, возникающих при дискретизации задачи Стокса с переменной вязкостью. Поскольку в реальных приложениях отношение максимальной вязкости к минимальной может быть
6Д.М. Климов, А.Г. Петров, Д.В. Гергиевский. Вязкопластические течения. Динамический хаос, устойчивость, перемешивания. Москва: Наука, 2005.
очень большим, к такому методу решения СЛАУ в настоящей работе предъявляется требование независимости (или слабой зависимости) количества итераций от этого отношения, а также от шага сетки.
2. Теоретический анализ эффективности предлагаемого метода решения систем линейных алгебраических уравнений. В частности, получение оценок скорости сходимости в терминах экстремальных значений коэффициента вязкости и других параметров систем уравнений.
3. Проверка эффективности предлагаемого метода на модельных задачах: регуляризованной задаче моделирования течения среды Бингама в канале и каверне, а также на линейной задаче, возникающей при моделировании всплытия раскаленного пузыря в мантии Земли.
На защиту выносятся следующие основные результаты и положения:
1. Доказано обобщенное неравенство Нечаса для случая переменной вязкости. Доказано обобщение неравенства для случая, когда область представлена в виде объединения непересекающихся подобластей.
2. Предложен переобуславливатель для дополнения по Шуру для дискретной задачи Стокса, учитывающий переменную вязкость. Получена оценка на собственные значения переобусловленного дополнения по Шуру. Получена оценка скорости сходимости метода Узавы-сопряженных градиентов.
3. Предлагаемый переобуславливатель применен для численного решения рсгуляризованной задачи Бингама и линейной задачи , возникающей при моделировании всплытия раскаленного пузыря в магме. Для задачи о течении среды Бингама в канале получены оценки собственных значений.
4. Проведены численные эксперименты с использовании предлагаемого переобуславливателя для задачи течения среды Бингама и линейной задачи, возникающей при моделировании всплытия раскаленного пузыря.
Научная новизна работы
Научная новизна диссертации заключается в доказательстве обобщения неравенства Нечаса для переменной вязкости, построении переобуславливателя для дополнения по Шуру, учитывающего переменную вязкость, а также в анализе его эффективности и эффективности метода Узавы-сопряженных градиентов с применением предлагаемого переобуславливателя. Проведено сравнение с переобуславливанием дополнения по Шуру при помощи матрицы масс.
Практическая ценность работы
В диссертации предложен итерационный метод решения задачи Сток-са с переменной вязкостью с переобуславливателем для дополнения по ШУРУ> учитывающим переменную вязкость. Как теоретические оценки, так и результаты применения к двум модельным задачам показывают, что количество итераций практически не зависит как от шага сетки, так и от отношения максимального значения вязкости к минимальному. Предлагаемый переобуславливатель может быть легко реализован в рамках существующих программных пакетов вычислительной гидродинамики.
Публикации. По теме диссертации опубликовано 6 работ, 3 из них — в изданиях из "Перечня ведущих рецензируемых научных журналов и изданий, в которых должны быть опубликованы основные научные результаты диссертаций на соискание ученых степеней доктора и кандидата наук".
Апробация работы.
1. Международная конференция "Х-я Белорусская математическая конференция". Белоруссия, Минск, 2008.
2. Доклад на семинаре Кафедры вычислительной математики Математического факультета Технического университета Дортмунда под руководством Ш. Турека. Германия, Дортмунд, 2009.
3. XVI Международная конференция студентов, аспирантов и молодых ученых "Ломоносов". Москва, 2009.
4. International Workshop on "Computational Mathematics and Applications". Финляндия, Тампере, 2009.
5. XVII Международная конференция студентов, аспирантов и молодых ученых "Ломоносов". Москва, 2010.
6. Доклад на семинаре Кафедры вычислительной математики, Механико-математический факультет МГУ под руководством Г.М. Ко-белькова. Москва, 2010.
7. Доклад на семинаре Кафедры механики композитов под руководством В.И. Горбачева, Механико-математический факультет МГУ. Москва, 2010.
8. Доклад на семинаре "Технологии математического моделирования
течений со свободной границей" под руководством Ю.В. Василевского и М.А. Ольшанского, ИВМ РАН. Москва, 2010.
9. Доклад на семинаре Кафедры математического моделирования МЭИ(ТУ) под руководством A.A. Амосова и Ю.А. Дубинского. Москва, 2010.
Структура и объем диссертации.
Диссертация состоит из введения, трех глав, заключения и списка литературы, содержащего 70 наименований. Общий объем работы — 105 страниц, работа включает 22 иллюстрации и 18 таблиц.
Содержание работы
Во введении обоснована актуальность диссертационной работы, сформулирована цель и аргументирована научная новизна исследований, показана практическая значимость полученных результатов, представлены выносимые на защиту научные положения.
В первой главе приводится формулировка задачи Стокса с переменной вязкостью, доказывается обобщение неравенства Нечаса на случай весовой нормы. Приводится матричная постановка дискретной задачи, предлагается переобуславливатель для дополнения по Шуру и при помощи обобщенного неравенства Нечаса получаются оценки его эффективности.
В §1.1 дается формулировка задачи Стокса с переменной вязкостью. Основные уравнения, рассматриваемые в первой главе, имеют вид
-divi/(x)Du + Vp = f в П —div u = 0 в il
u = uo на ЗП. 7
Через Du обозначен тензор скоростей деформации |(Vu + VTu). Вязкость v в первой главе считается зависящей только от пространственных координат х.
В §1.2 получен основной теоретический результат первой главы диссертации — неравенство Нечаса в весовой норме. Оно является обобщением хорошо известного неравенства Нечаса 7.
(о, div v) О < Со < inf sup ,, , „_ ,,,
которое, в свою очередь, является непрерывным аналогом условия Jla-дыженской-Бабушки-Брецци (LBB): если константа
CQh = mf sup -—rr——- (1)
^ ^v^JliAllllVVhH ^
равномерно по шагу сетки отделена от нуля, то пара конечноэлементных пространств V/г и Q^ (для скорости и давления, соответственно) приводит к устойчивой дискретной задаче. Под равномерным отделением от нуля понимается существование такой константы с* > 0, что при всех h выполняется со,ь > с*. Результат сформулирован в виде Теоремы 1.
Теорема 1. (Неравенство Нечаса в весовой норме) Пусть область fi связна и имеет кусочно-гладкую липшицеву границу. Предположим, что функция v > 0 достаточно гладкая, чтобы все нормы в (3)-(4) имели смысл. Тогда для всякой функции q 6 для которой
также выполнено (д,^-*) = 0, будет выполнено неравенство
. CrWv'tqW < »Р (2)
veHj |HDv||
Константа cv в двумерном случае задается формулой
c^coti-MMlHMIvHiu/r1 (3)
7 J. NeCas. Les Méthodes Directes en Théeorie des Équations Elliptiques. Paris: Masson, 1067.
с произвольными к > 2 и г > константой с(к,г), зависящей от констант неравенств вложения Нд(П) в Ьь( Г2) сЬ = Ь(к,г) и константой со, зависящей от константы со б неравенстве Нечаса. В трехмерном случае константа си определяется по формуле
с^са^ + срЧ^и-ЦьГ1 (4)
с произвольными к > 3 иг =
Из Теоремы 1 следует обобщение на случай, когда область £1 представлена в виде объединения конечного числа подобластей П,-, которое сформулировано как Теорема 2.
Теорема 2. (Обобщение Теоремы 1 на случай нескольких подобластей) Пусть П = и-^Па, Г24 — связные непересекающиеся подобласти с кусочно-гладкой липшицевой границей; функция V > 0 — кусочно-гладкая на всей области П (т.е. гладкая на каждой из подобластей 0,{). Тогда для любой такой € для которой выполняются условия
(9. = (ЯгУ'^ЩЪ) = 0 будет выполнено неравенство (2) с кон-
стантой с„ = тт1<*<лг где — константы, определяемые
соотношениями (3) или (4), вычисленные в областях
Иногда при помощи удачного разбиения удается добиться того, что пищ Су (П;) > с„(Г2). Ценой этого улучшения являются два условия ортогональности на каждой подобласти: (д, и'1)^^) — (я, ""^¿»(юо = О-
В §1.3 и §1.4 рассмотрены дискретизация задачи и построение СЛАУ. Дискретизация проводилась двумя методами, а именно — методом конечных элементов (¡зоР2-Р1) и методом конечных разностей на разнесенных сетках. Дается описание построения дискретной задачи при помощи обоих методов. Излагается вариационная постановка, дается определение конечноэлементных пространств для метода конечных элементов, вводятся сетки и задаются сеточные аналоги для непрерывных операторов
при использовании метода конечных разностей. Оба метода приводят к решению системы линейных алгебраических уравнений с седловой точкой вида
" :•)(:)-(:) •
В §1.5 описывается метод решения полученной дискретной задачи. Для решения используются методы на подпространствах Крылова (MINRES или GMRES) со специальным переобуславливателем вида
( Á о\ ( Á 0 \
V = . или ?!= 1 . 6
\os) \B-S)
Здесь А — переобуславливатель для блока A, S ~ переобуславли-ватель для дополнения по Шуру S = ВА~1ВТ. Лемма 2 дает оценку эффективности матрицы масс, примененной в качестве переобуславли-вателя для дополнения по Шуру в случае переменной вязкости.
Лемма 2. Если выполняется условие LBB, то есть для константы со,л из (1) выполнено co,/¡ > с* > 0, где с* не зависит от h, то справедливы неравенства
c*2¡/m< S < v~}„M на пространстве Qft.
Следствие.
cond (S'lS) < • если S = M.
с ^min
При большом отношении максимального значения вязкости к минимальному, использование матрицы масс становится неэффективным. Поэтому вместо стандартной матрицы масс для дополнения по Шуру предлагается применять переобуславливатель
при использовании метода конечных элементов (через 1/>; обозначены базисные функции для давления) и
М„ = (Цад^^у)} (8)
при использовании метода конечных разностей. Переобуславливатель М„ строится с учетом функции V, что важно при большом отношении максимального значения вязкости к минимальному.
Эффективность переобуславливателя М„ зависит от констант с„ и С„ в соотношении спектральной эквивалентности
с„М„ < Я < С„М„. (9)
Из Леммы 2 следует оценка константы с„ в соотношении с„Ми < Б, когда используется переобуславливатель М„, задаваемый формулой (7). Лемма 3. Для константы си из соотношения (9) верна оценка
V. «2 "тт
Си > С -,
Цпах
с* — константа из условия ЬВВ (1). Оценку константы С„ дает Теорема 4.
Теорема 3. Пусть и е 1°°(П) и П 6 Е"* и М„ задано в (7). Тогда С„ =
Лемма 3 и Теорема 3 дают оценку эффективности переобуславливателя Ми не лучше, чем дает Теорема 3 для матрицы масс М. Оценка с„ может быть улучшена, если сузить пространство функций Последующие рассуждения верны, если справедлив дискретный аналог Теоремы 1, сформулированный в виде Предположения 1.
Предположение 1. Предположим, что справедливо следующее. Пусть рассматриваемая расчетная область — многоугольник и {т^} — семейство регулярных (то есть таких, где треугольники мо-
гут либо иметь общее ребро, либо общую вершину, либо не пересекаться вовсе) триангуляций. При этом ft = L>KaThK. Пусть также для каждой триангуляции тд определены пространства Уд и Qh и для пары пространств Уд и Qд выполнено LBB-условие (1). Предположим, что для всех qh € Q/, С для которых также выполнено условие
ортогональности (qh= 0, справедливо неравенство
< sup (10)
VhZVk ||I/2Dva||
с некоторой константой с^д, которая выражается в двумерном случае как
Си,/, = со,л(1 +c(fc,r)||i/i||£t||Vi/-i||ir)-1 (И)
с произвольными к > 2 и г > константой c(k,r), зависящей от констант неравенств вложения Hq(£2) в L®(П) с t = t(k,r) и константой со,л, зависящей от константы с* в неравенстве LBB и не зависящей от h. Предположим, что в трехмерном случае константа cv записывается как
c„,h = Mi + (12)
с произвольными к > 3 и г = ^^ ы константой со,д, зависящей от константы в неравенстве LBB и не зависящей от h.
Если Предположение 1 справедливо, то верно его обобщение на случай, когда расчетная область представлена в виде объединения подобластей.
Теорема 4. (Дискретный аналог Теоремы 2.) Пусть расчетная область О, представлена в виде объединения связных подобластей П<, то есть П = Uiljfij. Пусть также границы подобластей проходят только по ребрам треугольников. Тогда неравенство (10) бу-
дет выполняться для всех функций д/, £ <0)/4, для которых справедливы соотношения {ян,^-1)^^) — (яь> ""^¿»(Од = 0, с константой с^н = тщ1<,<лг где — константа, вычисленная па г-й
подобласти.
Если верно Предположение 1 и, как следствие, Теорема 4, то справедлив результат, сформулированный как Теорема 5. Хотя формально она и не дает оценку константы с„ из неравенства (9) (это связано с тем, что на функции ф1 накладываются дополнительные условия ортогональности), но из нее следует оценка на собственные значения переобусловленного дополнения по Шуру, начиная с 2ДГ + 1-го, где N — количество подобластей, на которые разбита расчетная область.
Теорема 5. Пусть функция д/, е (ф/, С удовлетворяет усло-
виям ортогональности (<й, г/_1)х,2(п^ = (%> = 0 на каждой
подобласти. Если справедливо Предположение 1, то выполняется соотношение с2 Л (М¡,д, д) < (Бд, д), где с„,/, — константа из Теоремы 4В §1.6 дается оценка собственных значений для переобусловленной матрицы дополнения по Шуру. Для Л(М~Х5) в диссертации получены оценки
0 = Л1<Л2<Лз---<Лт<Сг, с*2—<л2, с111<х2М+ь (13)
^тах
где <1 — размерность пространства, N — количество подобластей, на которые разбита область Г2, с* — константа из условия ЬВВ (1), а — константа из Теоремы 4.
Далее в этом разделе приводится оценка скорости сходимости метода Узавы-сопряженных градиентов для систем линейных алгебраических уравнений с блочной матрицей вида (5) и переобуславливателем для дополнения по Шуру М„, которая сформулирована в Теореме 6.
Теорема 6. Пусть система с матрицей А решается точно. Тогда
для погрешности на к-й итерации метода Узавы-сопряженных градиентов решения задачи (5) при использовании переобуславливателя Ми выполнены оценка
к-2ЛГ+1
/ j \ ¿i* — i (
к
j \ 2JV-1
а щ л
где с* — константа из LBB-условия (1).
Показано, что влияние возможного присутствия 2N — 1 малого собственного значения приводит к незначительному росту числа итераций метода Узавы-сопряженных градиентов.
Во второй главе описывается применение теоретических результатов
первой главы к решению двух модельных задач: задачи Бингама о течении вязкопластичной жидкости в канале и задаче моделирования мантийной конвекции.
В §2.1 даются уравнения и определяющие соотношения модели Бингама. Соотношения, связывающие тензор скоростей деформации Du и девиатор тензора напряжений г, имеют вид
Di}u
r,J = 2/i£»4u + т5т=г-г, если |Du| > 0;
I I
|т| < Ta, если |Du| = 0,
где т3 — предел пластичности, а ц — коэффициент вязкости. В жидкой зоне коэффициент пропорциональности между девиатором тензора напряжений и тензором скоростей деформаций часто называют эффективной вязкостью.
В диссертации рассматривается подход к численному решению задачи Бингама, основанный на регуляризации, последующей дискретизации на каждом шаге нелинейного процесса и применении нового переобусловленного итерационного метода для решения линейных вспомогательных
задач. Важным свойством предлагаемого в настоящей диссертации метода является отсутствие параметров, требующих тонкой настройки.
В §2.2 дается описание регуляризованной задачи. Регуляризация — подход к решению задачи Бингама, при котором среда рассматривается как жидкая во всей расчетной области, а в жестких зонах эффективная вязкость считается конечной. Величина вязкости в жестких зонах существенным образом определяется параметром регуляризации е. Рассматриваются два варианта регуляризации:
1/£(|Ои|) = 2ц Н--. ' = — Берковьера-Энгельмана;
л/ри2 + £2\
1 _ ехр(-М)
//£(|Ви|) = 2ц + т3-——;—-— — Папанастасио.
|Ои|
В §2.3 определяется итерационный метод, используемый для решения регуляризованной задачи Бингама. Нелинейные уравнения решаются методом Пикара:
и ' ' ^(и"'1)-1 о
Р" / \ Рп~1 '
-:
где через обозначен линейный при заданной функции а оператор
/^(МС И
^ -Фу 0 ]
а через ^(и"-1)-1 — приближенное решение линейной задачи (5), возникающей при дискретизации оператора Р.
В §2.4 выводится оценка собственных значений матрицы М~1Б для задачи о течении среды Бингама в канале. Это одна из немногих задач
и =
вязкопластичности, где известно аналитическое решение u = (u,v): |(1 - 2ts)2, если § - ts < у < \ + т„,
I [(1 - 2та)2 - (1 - 2rs - 2у)2] , если 0 < у < ± - т„ J [(1 - 2ts)2 - (2у - 2т, - I)2] , если 1 > у > \ + т3, v = 0, р = -х.
С помощью разбиения области на три подобласти, соответствующие ядру течения в центре и жидким зонам по краям, непосредственно оценена константа в обобщенном неравенстве Нечаса (2) и, используя соотношения (13), получена оценка на собственные значения переобусловленного матрицей Ми дополнения по Шуру: для любого s 6 (0,1] существует не зависящая от е константа сг > 0, что
cíe < А2, c2es < А7 < • • • < Am < d - 2,
где константа ci не зависит от е.
Последний раздел второй главы посвящен применению результатов первой главы к задаче моделирования мантийной конвекции. Рассматривается линейная задача, возникающая при моделировании всплытия раскаленного пузыря в магме. Считается, что вязкость v в этом случае зависит только от температуры Т, а температура — только от пространственных координат. Приводятся алгоритм разбиения области на подобласти для данной задачи и вычисленные значения констант в обобщенном неравенстве Нечаса (2).
Таким образом, во второй главе описано применение метода, предлагаемого в первой главе, к модельным задачам и получены оценки его эффективности.
В третьей главе представлены результаты численного решения трех модельных задач, включая течение среды Бингама в канале и каверне,
16
а также линейную задачу, возникающую при моделирования всплытия раскаленного пузыря в магме.
В §3.1 и §3.2 представлены и проанализированы результаты решения задач течения среды Бингама в канале и в каверне, соответственно. Приводятся графики, на которых изображены все собственные значения пс-реобусловленного дополнения по Шуру с использованием как матрицы масс М, так и предлагаемого в работе переобуславливателя Ми. Приводится количество нелинейных и линейных итераций для разных конфигураций метода решения СЛАУ. Расчеты проводились с обеими регуля-ризациями и разными значениями параметра £ от 10_х до 10~5. Дано также сравнение блочно - диагонального V и блочно - треугольного (6) переобуславливателей для матрицы Л. Вычисления показали, что для данной модельной задачи численное решение хорошо приближает аналитическое и выбор регуляризации слабо влияет на эффективность метода. При этом, использование блочно - треугольного переобуславливателя Т\ приводит к заметно более высокой скорости сходимости, чем использование блочно - диагонального (6).
Для задачи о каверне представлены графики с приближениями жестких зон для различных значений предела пластичности г5 и рсгуляри-зационного параметра е. Эти результаты (вид и размеры жестких зон) хорошо согласуются с результатами расчетов в других работах. Это свидетельствует о том, что выбранные численные параметры позволяют воспроизводить важные физические эффекты.
Из численных экспериментов по решению регуляризованной задачи Бингама в канале и в каверне можно сделать следующий вывод: предлагаемый в диссертации итерационный метод решения систем линейных алгебраических уравнений показывает практически полное отсутствие зависимости числа итераций как от шага сетки, так и от параметра регу-
ляризации (соответственно, от отношения максимального значения вязкости к минимальному), что хорошо согласуется с теоретическими оценками (в задаче о течении в канале).
В §3.3 обсуждаются результаты численного решения линейной задачи, возникающей при моделировании всплытия раскаленного пузыря в магме. Сравнивается эффективность переобуславливателей М и М„. Численные эксперименты показали, что число линейных итераций в задаче моделирования мантийной конвекции с увеличением отношения максимальной вязкости к минимальной растет незначительно. Хотя теоретические оценки, судя по всему, не являются оптимальными в этом случае.
В заключении сформулированы основные результаты диссертации.
Основные результаты
1. Доказано обобщенное неравенство Нечаса для случая переменной вязкости. Доказано обобщение неравенства для случая, когда область представлена в виде объединения непересекающихся подобластей.
2. Предложен переобуславливатель для дополнения по Шуру для дискретной задачи Стокса, учитывающий переменную вязкость. Получена оценка на собственные значения переобусловленного дополнения по Шуру. Получена оценка скорости сходимости метода Узавы-сопряженных градиентов.
3. Предлагаемый переобуславливатель применен для численного решения регуляризованной задачи Бингама и линейной задачи , возникающей при моделировании всплытия раскаленного пузыря в магме. Для задачи о течении среды Бингама в канале получены
оценки собственных значений переобусловленного дополнения по Шуру и теоретически доказано, что зависимость от параметра регуляризации оценок собственных значений для предлагаемого метода меньше, чем для обычно используемого переобуславливателя на основе матрицы масс.
4. Проведены численные эксперименты, хорошо согласующиеся с теоретическими результатами для задачи течения среды Бингама в канале и показавшие практическую независимость числа линейных итераций как от шага сетки, так и от отношения максимального значения вязкости к минимальному для задачи течения среды Бингама в каверне и линейной задачи, возникающей при моделировании всплытия раскаленного пузыря в магме при использовании предлагаемого переобуславливателя.
Благодарности
Автор выражает глубокую благодарность доктору физико-математических наук, профессору М.А. Ольшанскому за постановку задачи и помощь в работе. Автор выражает благодарность всем сотрудникам кафе-ды вычислительной математики за поддержку и внимание к работе.
Работы автора по теме диссертации
1. P.P. Grinevich, М.А. Olshanskii. An iterative method for the Stokes type problem with variable viscosity // SIAM Journal on Scientific Computing, Vol.31, Iss.5, 2009, pp.3959-3978.
2. П.П. Гриневич, М.А. Ольшанский. Итерационный метод решения
регуляризованной задачи Бингама // Вычислительные методы и программирование, Том 10, 2010, стр. 78-87.
3. П. П. Гриневич. Об итерационном методе решения задачи Стокса с переменной вязкостью // Вестник МГУ, №3, 2010, стр. 38-41.
4. П.П. Гриневич., М.А. Ольшанский Итерационные методы для регуляризованной модели среды Бингама. X Белорусская Математическая Конференция: Тезисы докладов Международной научной конференции. Минск, 3-7 ноября 2008 - Мн.: Институт математики НАН Беларуси, 2008, с. 5-6.
5. P.P. Grinevich, М.А. Olshanskii An iterative method for solving regularized Bingham fluid equations. Тезисы докладов международной конференции "International Workshop on Computational Mathematics and Applications". Department of Mathematics, Tampere Institute of Technology, 2009, c.9.
6. П.П. Гриневич. Об итерационном методе решения задачи Стокса с переменной вязкостью. Материалы Международного молодежного научного форума "ЛОМОНОШВ-2010", М.: МАКС Пресс, 2010.
Отпечатано в отделе оперативной печати Геологического ф-та МГУ Тираж 1С0ЭКЗ. Заказ № 6
Введение
Глава 1. Задача типа Стокса с переменной вязкостью
1.1. Постановка задачи.
1.2. Неравенство Нечаса в весовой норме
1.3. Дискретизация задачи
1.4. Матричная постановка дискретной задачи
1.5. Итерационный метод решения СЛАУ.
1.6. Оценка собственных значений и сходимость итерационных методов.
Глава 2. Задача Бингама и задача мантийной конвекции
2.1. Постановка задачи Бингама
2.2. Регуляризованная задача Бингама.
2.3. Нелинейный итерационный метод для задачи Бингама
2.4. Применение оценок на собственные значения к задаче течения среды Бингама в канале.
2.5. Применение оценок на собственные значения к задаче мантийной конвекции
Глава 3. Численные эксперименты
3.1. Течение среды Бингама в канале
3.2. Течение среды Бингама в каверне.
3.3. Задача моделирования мантийной конвекции
3.4. Выводы из численных экспериментов.
Решение многих современных прикладных задач приводит к системам линейных алгебраических уравнений (СЛАУ) с седловой точкой. Характерной особенностью таких систем является знаконеопределенность. В симметричном случае имеются как положительные, так и отрицательные собственные значения. Подобные системы возникают, например, в теории динамических систем [37], в теории упругости [22, 28], в экономике [16, 32, 42], в теории оптимального управления [17, 18], при моделировании электромагнитных явлений [21, 54] и во многих других областях. Важной областью, требующей решения задач с седловой точкой, является численное решение линеаризованных уравнений Навье-Стокса, описывающих течение несжимаемой вязкой жидкости [38, 39, 55, 62, 63, 67]. Уравнения Навье-Стокса являются основными уравнениями гидродинамики и, соответственно, играют важную роль в современной науке. Отметим, что в ходе их решения, как правило, возникает необходимость проводить вычисления на мелких сетках, а значит решаемые системы имеют большую размерность. Знаконеопределенность делает процедуру выбора метода решения таких систем нетривиальной, так как многие эффективные методы решения СЛАУ, см., например, [58], пригодны только для систем с положительно определенной матрицей.
Уравнения Навье-Стокса во многих случаях хорошо описывают поведение жидкостей и газов. Однако, многие вещества в природе описываются моделями с переменным коэффициентом вязкости, зависящим от внешних факторов. Примером могут служить биологические жидкости (например, кровь), нефть, зубная паста, кетчуп, крахмал, разведенный в воде и многие другие вещества. В некоторых случаях, например, кетчуп, наблюдается уменьшение вязкости при возрастании внутренних напряжений. В других, таких как крахмал, поведение обратное — при возрастании напряжений вязкость резко увеличивается. Для моделирования подобных подобных веществ можно рассматривать модифицированные уравнения Навье-Стокса, при этом вязкость является не постоянным параметром среды, а функцией, зависящей от динамических, кинематических или других характеристик среды в данной точке пространства, например, тензора скоростей деформации, давления, температуры, и т.д. Возникающие при дискретизации линейные системы характеризуются тем, что их число обусловленности зависит от отношения максимальной вязкости к минимальной. Данное обстоятельство предъявляет дополнительное требование к методам решения СЛАУ, а именно — независимость числа итераций от отношения максимального значения коэффициента вязкости к минимальному и от градиента коэффициента вязкости как функции пространственной переменной. Эффективное решение уравнений Навье-Стокса с переменной вязкостью является темой ряда современных исследований. Упомянем некоторые из них. В статье [70] рассматривается моделирование мантийной конвекции в планетарном масштабе. Рассматривается влияние разных моделей вязкости, в том числе — распределенной по слоям, а также вязкости, зависящей от температуры. В [49] получены оценки сходимости конечноэлементного решения задачи Стокса с интерфейсом к решению непрерывной задачи. Вязкость при этом кусочно-постоянна и имеет скачок. Работа [56] посвящена построению эффективных итерационных методов для несжимаемой задачи Стокса с переменной вязкостью, в ней рассматриваются задачи, связанные с моделированием процесса экструзии и с геологией. В статье [68] представлено построение равномерно устойчивых по скачку коэффициента вязкости конечноэлементных аппроксимаций уравнений Дарси-Стокса-Бринкмана. В работе [52] получено тЛ^эир неравенство для задачи Стокса с интерфейсом, которое равномерно по скачку вязкости между подобластями. В [51] авторами предложен метод анализа переобуславливателей для дополнения по Шуру для абстрактных задач с седловой точкой. Применение этого метода к задаче Стокса с интерфейсом приводит к построению переобуславливателя, являющегося устойчивым по отношению к ряду параметров.
В численном анализе методов решения задачи Стокса важную роль играет условие ЬВВ (Ладыженской-Бабушки-Брецци) [23] и его непрерывный аналог, неравенство Нечаса [46]. В диссертации неравенство Нечаса обобщается на случай переменной вязкости и на основе этого обобщения получены оценки эффективности предлагаемого итерационного метода. Некоторые другие обобщения неравенства Нечаса получены в работах [1] и [6].
Одной из задач, где возникают уравнения Навье-Стокса с переменной вязкостью, является моделирование среды Бингама. Среда Бинга-ма является вязко-пластичной средой, которая при напряжениях ниже порогового значения ведет себя как твердое тело, а при превышающих пороговое значение — как вязкая жидкость. Примерами веществ, описываемых моделью Бингама являются зубная паста, различные суспензии, например грязь или незастывший бетон, геоматериалы. Впервые модель независимо была предложена Шведовым [60] и Бингамом [20] для описания движения суспензий в условиях чистого сдвига (одномерная модель). Позднее Генки [2] и Ильюшин [4] предложили пространственное обобщение уравнения состояния Шведова-Бингама и решили ряд задач для случая плоских течений вязко-пластичной среды. Течению среды
Бингама посвящено большое количество литературы, среди отечественной можно отметить труды П.М. Огибалова и А.Х. Мирзаджанадзе [10], П.П. Мосолова и В.П. Мясникова [9], а также Д.М. Климова, А.Г. Петрова, Д.В. Георгиевского. [5]. Задача Бингама является более сложной, чем задача Навье-Стокса. Одним из подходов к ее численному решению является регуляризация — модель, когда среда Бингама рассматривается как жидкость с переменной вязкостью.
Задача Бингама трудно поддается математическому анализу и ее точные решения найдены только для узкого круга модельных задач, например, для течения среды между двумя параллельными пластинами [5]. По этой причине численные методы, зачастую, являются единственным способом анализа многих процессов. Численному решению задачи Бингама посвящены, например, следующие публикации: [19, 29, 30, 33, 44, 48, 53, 57, 65].
Задачи с переменной вязкостью появляются и во многих других научных областях, например, в геологии. В мантии Земли температура неоднородна, а вязкость магмы напрямую зависит от температуры. Численные методы позволяют проводить моделирование внутреннего строения Земли и эта тема рассматривается сейчас многими учеными. Например, в [45, 69] обсуждаются вопросы численного моделирования мантийной конвекции, в [43] — вопросы численного моделирования миграции магмы.
Цель работы.
Работа преследует следующие цели.
1. Построение эффективного метода решения систем линейных алгебраических уравнений, возникающих при дискретизации задачи Стокса с переменной вязкостью. Поскольку в реальных ириложениях отношение максимальной вязкости к минимальной может быть очень большим, к такому методу решения СЛАУ в настоящей работе предъявляется требование независимости (или слабой зависимости) количества итераций от этого отношения, а также от шага сетки.
2. Теоретический анализ эффективности предлагаемого метода решения систем линейных алгебраических уравнений. В частности, получение оценок скорости сходимости в терминах экстремальных значений коэффициента вязкости и других параметров систем уравнений.
3. Проверка эффективности предлагаемого метода на модельных задачах: регуляризованной задаче моделирования течения среды Бингама в канале и каверне, а также на линейной задаче, возникающей при моделировании всплытия раскаленного пузыря в мантии Земли.
Используемые обозначения.
• Обычным шрифтом (например, р) обозначаются скалярные функции, коэффициенты разложения по конечноэлементным базисам и вектора в матричной записи линейных систем.
• Полужирным шрифтом (например, V) обозначаются векторные и тензорные функции и операторы. Используются обозначения
Уи = дих дих дх ду диу диу дх ду дих диг дъ диг дх ду сИУБ = а^Б1 ^ сИУ Г>2 у } где Ю — тензорная функция, действующая из К3 в М3х3, и Т>к =
• Индекс Н у функции или оператора указывает на дискретную (ко-нечноразностную или конечноэлементную) аппроксимацию данной функции или данного оператора.
• Через 1Р обозначается пространство Лебега с индексом суммирования р, через || • \\ьр — норма в I/. Запись || • || означает норму в Ь2. Скалярное произведение в Ь2 и евклидово скалярное произведение обозначаются как (•,•) и (•,■), соответственно.
Используется — пространство всех функций / из таких, что /сЮ, = 0 и — пространство функций / из Ь2{0), для которых выполняется условие /ь>~1сЮ, = 0 с некоторой функцией и > 0 такой, что V"1 € Ь°°.
• Символом Нц обозначено пространство интегрируемых вектор-функций, градиент которых, понимаемый в обобщенном смысле, принадлежит Ь2(Жах<с/ = 2,3 и имеющих нулевой след на границе расчетной области.
Краткое содержание работы.
В первой главе приводится формулировка задачи Стокса с переменной вязкостью, доказывается обобщение неравенства Нечаса на случай весовой нормы. Приводится матричная постановка дискретной задачи, предлагается переобуславливатель для дополнения по Шуру и и при помощи обобщенного неравенства Нечаса получаются оценки его эффективности.
В первом разделе дается формулировка задачи Стокса с переменной вязкостью. Основные уравнения, рассматриваемые в первой главе, имеют вид
-divi/(x)Du + Vp = f в ft —div u = 0 в Q ii = U() на dQ,.
Через Du обозначен тензор скоростей деформации |(Vu + VTu). Вязкость v в первой главе считается зависящей только от пространственных координат х.
Во втором разделе получен основной теоретический результат диссертации — неравенство Нечаса в весовой норме. Оно является обобщением хорошо известного неравенства Нечаса (g,divv) О < с0 < mi sup и „ „ и, ^veHjIMI IIVvH' которое, в свою очередь, является непрерывным аналогом условия Ладыженской - Бабушки - Брецци (LBB): если константа qh, divv/j) со,ft = mf sup т—¡г=—¡7 чьеQhVheYh \Ш\ llVvftll равномерно по шагу сетки отделена от нуля, то пара конечноэлементных пространств Уд и Q/j (для скорости и давления, соответственно) приводит к устойчивой дискретной задаче. Под равномерным отделением от нуля понимается существование такой константы с* > 0, что при всех h выполняется co,ft > с*. Результат сформулирован в виде Теоремы 1, суть утверждения которой заключается в выполнении неравенства для функций д £ Ь2(С1), на которые наложены условия ортогональности
Приводится важное следствие Теоремы 1 — обобщение на случай, когда расчетная область Г2 представлена в виде объединения конечного числа подобластей Г^. В этом случае константа си для всей области О, может быть оценена снизу минимумом из всех констант сДГ^), вычисленных на подобластях.
Иногда при помощи удачного разбиения удается добиться того, что пищ с„(П). Ценой этого улучшения являются два условия ортогональности на каждой подобласти: (д, — (<7> и -1/2)щпг) = о.
Разделы 3 и 4 посвящены дискретизации задачи и построению СЛАУ. Дискретизация проводилась двумя методами, а именно — методом конечных элементов (¡боР2-Р1) и методом конечных разностей (МАС-схема). Дается описание построения дискретной задачи при помощи обоих методов. Излагается вариационная постановка, дается определение конечно-элементных пространств для метода конечных элементов, вводятся сетки и задаются сеточные аналоги для непрерывных операторов при использовании метода конечных разностей. Оба метода приводят к решению системы линейных алгебраических уравнений с седловой точкой с матрицей вида где А = Ат > 0.
В пятом разделе описывается метод решения полученной дискретной 0. Получены оценки на константу с1 V задачи. Основная идея заключается в использовании метода на подпространствах Крылова (MINR.ES или СМНЕЯ) со специальным переобу-славливателем вида
Р=1 Л ° I (3)
О 5 или
7 А О
4)
Здесь А — переобуславливатель для блока А, 5' — переобуславлива-тель для дополнения по Шуру 5 = ВА~гВт. Для дополнения по Шуру предлагается применять переобуславливатель
М„}ц := (5) при использовании метода конечных элементов (через ф{ обозначены базисные функции для давления) и (6) при использовании метода конечных разностей. Переобуславливатель Ми при построени учитывает функцию что важно при большом отношении максимального значения вязкости к минимальному. Сформулированы дискретные аналоги Теорем 1 и 2, обозначенные в тексте как Предположение 1 и Теорема 4, соответственно. Предположение 1 сводится к предположению справедливости дискретного аналога обобщенного неравенства Нечаса для пары конечноэлементных пространств, удовлетворяющих ЬВВ условию устойчивости. Далее доказывается, что
Зд,д) <(1{Мид,д) для всех qh 6 Qh и, если Предположение 1 справедливо, то выполняется соотношение cl,h (Mvq,q) < (Sq,q) для всех qh £ Qh таких, что (qh, v~l)щ^) = (qh, = 0. Здесь d — размерность пространства, т.е. d = 2 или d = 3, через q обозначен вектор коэффициентов разложения по базису функции qh £ Q/г, а с^л — константа из Теоремы 4. Таким образом устанавливается связь между неравенством Нечаса в весовых нормах и эффективностью переобуслав-ливателя Ми для дополнения по Шуру.
В шестом разделе первой главы дается оценка собственных значений для переобусловленной матрицы дополнения по Шуру. Она выводится из представления Куранта-Фишера (Sg, д) ли = max min ----г,
CeVfc 1 qeic± (Mvq, q) где символом Vfci обозначены все (к — 1)-мерные подпространства Rm. Для \(M~lS) в диссертации получены оценки
0 = Ai<A2<A3---<Am<d и clth< А2лг+ь (V где d — размерность пространства, N — количество подобластей, на которые разбита область а с^д — константа из Теоремы 4.
В этом разделе также приводится оценка скорости сходимости метода Узавы-сопряженных градиентов [34] для систем линейных алгебраических уравнений с матрицей (2) и переобуславливателем для дополнения по Шуру Ми и оценка количества дополнительных итераций, вызванных возможным присутствием 2N — 1 малого собственного значения.
Во второй главе описывается применение теоретических результатов первой главы к решению двух модельных задач: задачи Бингама о течении вязкопластичной жидкости в канале и задаче моделирования мантийной конвекции.
В первом разделе даются определяющие соотношения модели Бинга-ма, вводятся понятия жестких зон и эффективной вязкости в жидких зонах. Задача Бингама нелинейна, и ее эффективная вязкость зависит от тензора скоростей деформации, а значит, и от поля скоростей и.
Соотношения, связывающие тензор скоростей деформации Du и де-виатор тензора напряжений т, имеют вид
Dijи rlJ = 2(iDlJu + если |Du| > 0;
M < т8, если |Du| = 0, где т3 — предел пластичности, а ¡1 — коэффициент вязкости. В жидкой зоне коэффициент пропорциональности между девиатором тензора напряжений и тензором скоростей деформаций часто называют эффективной вязкостью. При этом в жестких зонах эффективную вязкость формально можно считать бесконечной.
В диссертации рассматривается подход к численному решению задачи Бингама, основанный на регуляризации, последующей дискретизации на каждом шаге нелинейного процесса и применении нового переобусловленного итерационного метода для решения линейных вспомогательных задач. Важным свойством предлагаемого в настоящей диссертации метода является отсутствие параметров, требующих тонкой настройки.
Во втором разделе дается описание регуляризованной задачи. Регуляризация — подход к решению задачи Бингама, при котором среда рассматривается как жидкая во всей расчетной области, а в жестких зонах эффективная вязкость считается конечной. Величина вязкости в жестких зонах существенным образом определяется параметром регуляризации е. Рассматриваются два варианта регуляризации: г, г/е(|Ви|) = 2[1 Н-- я = — Берковьера-Энгельмана;
Т>и2 + е2 |
1 ехр(Шн1) г/е(|Ви|) = 2/л + т3-|—-— — Папанастасио.
В третьем разделе определяется итерационный метод, используемый для решения регуляризованной задачи Бингама. Нелинейные уравнения решаются методом Пикара:
7 и"-1 \ I { и71 \ ( и"-1 , \Р'" / \ рп~г и"-1) О рп~1 / \ о где через Р обозначен оператор, линейный при фиксированной функции а ч ( -а^(|Ба|)В V 24 а) =
У -сНУ О а через ^(ип-1)-1 — приближенное решение линейной задачи с матрицей Л (2), возникающей при дискретизации оператора Р.
В четвертом разделе выводится оценка собственных значений матрицы М~1Б для задачи о течении среды Бингама в канале. Это одна из немногих задач, где известно аналитическое решение и = (и, у):
1 - 2т3)2, . если \ - т3 < у < \ + та, и = <1 [(1 - 2т3)2 - (1 - 2т3 - 2у)2] , если 0 < у < \ - т„ [(1 - 2т3)2 - (2у - 2та - I)2] , если 1 > у > § + г«,
V = 0, р = —х.
С помощью разбиения области на три подобласти, соответствующие ядру течения в центре и жидким зонам по краям, непосредственно оценена константа в обобщенном неравенстве Нечаса (1) и, используя соотношения (7), получена оценка на собственные значения переобусловленного матрицей Mv дополнения по Шуру: для любого s G (0,1] существует не зависящая от £ константа С2, что
С\£ < C2£s < Л7 < • • • < Am < d = 2, где константа с\ не зависит от е.
Последний раздел второй главы посвящен применению результатов первой главы к задаче моделирования мантийной конвекции. Рассматривается линейная задача, возникающая при моделировании всплытия раскаленного пузыря в магме. В диссертации не ставилась цель решать полную физическую задачу, для простоты считалось, что вязкость является функцией только от пространственных координат. Для данной задачи приводится алгоритм разбиения области на подобласти и вычисленные значения констант в обобщенном неравенстве Нечаса (1).
Таким образом, во второй главе описано применение метода, предлагаемого в первой главе, к модельным задачам и получены оценки его эффективности.
В третьей главе представлены результаты численного решения модельных задач, расмотренных во второй главе и задачи течения среды Бингама в каверне.
В первых двух разделах представлены и проанализированы результаты решения задач течения среды Бингама в канале и в каверне, соответственно. Приводятся графики, на которых изображены все собственные значения переобусловленного дополнения по Шуру с использованием как матрицы масс М, так и предлагаемого в работе нереобуславливателя Mv. Приводится количество нелинейных и линейных итераций для разных конфигураций метода решения СЛАУ. Расчеты проводились с обеими ре-гуляризациями с разными значениями регуляризационного параметра е. Дано также сравнение блочно - диагонального (3) и блочно - треугольного (4) переобуславливателей для матрицы Л. Вычисления показали, что для данной модельной задачи численное решение хорошо приближает аналитическое и выбор регуляризации слабо влияет на эффективность метода. При этом, использование блочно - треугольного переобуславли-вателя (4) приводит к заметно более высокой скорости сходимости, чем использование блочно - диагонального (3).
Для задачи о каверне представлены графики с приближениями жестких зон для различных значений предела пластичности т3 и регуляризационного параметра е. Эти результаты (вид и размеры жестких зон) хорошо согласуются с результатами расчетов в других работах. В [44], например, рассматривается регуляризованная задача, а в работах [30] и [50] — нерегуляризованная. Это свидетельствует о том, что выбранные численные параметры позволяют воспроизводить важные физические эффекты.
Из численных экспериментов по решению регуляризованной задачи Бингама в канале и в каверне можно сделать следующий вывод: предлагаемый в диссертации итерационный метод решения систем линейных алгебраических уравнений показывает практически полное отсутствие зависимости числа итераций как от шага сетки, так и от параметра регуляризации (соответственно, от отношения максимального значения вязкости к минимальному), что хорошо согласуется с теоретическими оценками (в задаче о течении в канале).
В третьем разделе третьей главы обсуждаются результаты численного решения линейной задачи, возникающей при моделировании всплытия раскаленного пузыря в магме. Сравнивается эффективность пере-обуславливателей М и Mv. Численные эксперименты показали, что число линейных итераций в задаче моделирования мантийной конвекции с увеличением отношения максимальной вязкости к минимальной растет незначительно. Более того, предлагаемый итерационный метод оказывается эффективнее, чем предсказывают теоретические оценки.
Апробация работы.
1. Международная конференция "Х-я Белорусская математическая конференция". Белоруссия, Минск, 2008.
2. Доклад на семинаре Кафедры вычислительной математики Математического факультета Технического университета Дортмунда под руководством Ш. Турека. Германия, Дортмунд, 2009.
3. XVI Международная конференция студентов, аспирантов и молодых ученых "Ломоносов". Москва, 2009.
4. International Workshop on "Computational Mathematics and Applications". Финляндия, Тампере, 2009.
5. XVII Международная конференция студентов, аспирантов и молодых ученых "Ломоносов". Москва, 2010.
6. Доклад на семинаре Кафедры вычислительной математики, Механико-математический факультет МГУ под руководством Г.М. Ко-белькова. Москва, 2010.
7. Доклад на семинаре Кафедры механики композитов под руководством В.И. Горбачева, Механико-математический факультет МГУ. Москва, 2010.
8. Доклад на семинаре "Технологии математического моделирования течений со свободной границей" под руководством Ю.В. Василевского и М.А. Ольшанского, ИВМ РАН. Москва, 2010.
9. Доклад на семинаре Кафедры математического моделирования МЭИ(ТУ) под руководством А.А. Амосова и Ю.А. Дубинского. Москва, 2010.
На защиту выносятся следующие основные результаты и положения:
1. Доказано обобщенное неравенство Нечаса для случая переменной вязкости. Доказано обобщение неравенства для случая, когда область представлена в виде объединения непересекающихся подобластей.
2. Предложен иереобуславливатель для дополнения по Шуру для дискретной задачи Стокса, учитывающий переменную вязкости. Получена оценка на собственные значения переобусловленного дополнения по Шуру. Получена оценка скорости сходимости метода Узавы-сопряженных градиентов.
3. Предлагаемый переобуславливатель применен для численного решения регуляризованной задачи Бингама и линейной задачи , возникающей при моделировании всплытия раскаленного пузыря в магме. Для задачи о течении среды Бингама в канале получены оценки собственных значений переобусловленного дополнения по Шуру и теоретически доказано, что зависимость от параметра регуляризации оценок собственных значений для предлагаемого метода меньше, чем для обычно используемого переобуславливателя на основе матрицы масс.
4. Проведены численные эксперименты, хорошо согласующиеся с теоретическими результатами задачи течения среды Бингама в канале и показавшие практическую независимость числа линейных итераций как от шага сетки, так и от отношения максимального значения вязкости к минимальному для задачи течения среды Бингама в каверне и линейной задачи, возникающей при моделировании всплытия раскаленного пузыря в магме при использовании предлагаемого переобуславливателя.
Публикации.
1. P.P. Grinevich, М.А. Olshanskii. An iterative method for the Stokes type problem with variable viscosity // SIAM Journal on Scientific Computing, Vol.31, Iss.5, 2009, pp.3959-3978.
2. U.U. Гриневич, М.А. Ольшанский. Итерационный метод решения регуляризованной задачи Бингама // Вычислительные методы и программирование, Том 10, 2010, стр. 78-87.
3. П.П. Гриневич. Об итерационном методе решения задачи Стокса с переменной вязкостью // Вестник МГУ, №3, 2010, стр. 38-41.
Структура и объем диссертации.
Диссертация состоит из введения, трех глав, заключения и списка литературы, содержащего 70 наименований. Общий объем работы — 105 страниц, работа включает 22 иллюстрации и 18 таблиц.
Основные результаты, полученные в диссертации, могут быть сформулированы следующим образом:
1. Доказано обобщенное неравенство Нечаса для случая переменной вязкости. Доказано обобщение неравенства для случая, когда область представлена в виде объединения непересекающихся подобластей.
2. Предложен переобуславливатель для дополнения по Шуру для дискретной задачи Стокса, учитывающий переменную вязкость. Получена оценка на собственные значения переобусловленного дополнения по Шуру. Получена оценка скорости сходимости метода Узавы-сопряженных градиентов.
3. Предлагаемый переобуславливатель применен для численного решения регуляризованной задачи Бингама и линейной задачи , возникающей при моделировании всплытия раскаленного пузыря в магме. Для задачи о течении среды Бингама в канале получены оценки собственных значений нереобусловленного дополнения по Шуру и теоретически доказано, что зависимость от параметра регуляризации оценок собственных значений для предлагаемого метода меньше, чем для обычно используемого переобуславливателя на основе матрицы масс.
4. Проведены численные эксперименты, хорошо согласующиеся с теоретическими результатами задачи течения среды Бингама в канале и показавшие практическую независимость числа линейных итераций как от шага сетки, так и от отношения максимального значения вязкости к минимальному для задачи течения среды Бинга-ма в каверне и линейной задачи, возникающей при моделировании всплытия раскаленного пузыря в магме при использовании предлагаемого переобуславливателя.
Заключение
1. И.А. Боровиков,Ю.А. Дубинский. Некоторые разложения модулей Соболева-Клиффорда и нелинейные вариационные задачи // Труды Математического института им. В. А. Стеклова РАН, Том260, 2008, с. 57-74.
2. Г. Генки. Пространственная задача упругого и пластического равновесия // Известия Академии наук СССР. Отдел технических наук, №2, 1937.
3. Г.Дюво, Ж.-Л. Лионе. Неравенства в механике и физике. Москва: Наука, 1980.
4. A.A. Ильюшин. Деформация вязко-пластичного тела // Ученые записки МГУ. Механика, №39, 1940, с. 3-81.
5. Д.М. Климов, А.Г. Петров, Д.В. Георгиевский. Вязкопластиче-ские течения. Динамический хаос, устойчивость, перемешивания. Москва: Наука, 2005.
6. A.M. Красногорский. О разложении пространства гармонических функций и некоторых приложениях // Доклады Академии наук, Том411, №4, 2006,с. 452-454.
7. В.И. Лебедев. Разностные аналоги ортогональных разложений, основных дифференциальных операторов и некоторых краевых задач математической физики // Журнал вычислительной математики и математической физики, Том 4, №3, 1964, с. 449-465.
8. C.B. Милютин. Исследование трёхпараметрического итерационного метода, ориентированного на решение двух классов задач с нелинейными седловыми операторами. Диссертация на соискание степени кандидата физико-математических наук. Москва, 2010.
9. П. П. Мосолов, В. П. Мясников. Вариационные методы в теории жёст-ковязкопластической среды. Москва: Наука, 1971.
10. П.М. Огибалов, А.Х. Мирзаджанадзе. Нестационарные движения вязкопластичных сред. Москва: Наука, 1970.
11. М.А. Ольшанский. Лекции и упражнения по многосеточным методам. Москва: ФИЗМАТЛИТ, 2005.
12. М.А. Ольшанский, Е.В. Чижонков. О наилучшей константе в inf-sup-условии для вытянутых прямоугольных областей // Математические заметки, Том 67, №3, 2001, с. 387-396.
13. С. Л. Соболев. Некоторые применения функционального анализа в математической физике, Москва: Наука, 1988.
14. Е.Е. Тыртышников. Методы численного анализа. Издательский центр "Академия", 2007.
15. Е.В. Чижонков. Релаксационные методы решения седловых задач. Москва: ИВМ РАН, 2002.
16. К. Arrow, L. Hurwitz, Н. Uzawa. Studies in Linear and Nonlinear Programming. Stanford, California: Stanford University Press, 1958.
17. A. Battermann, M. Heinkenschloss. Preconditioners for Karush-Kuhn-Tucker systems arising in the optimal control of distributed systems // in Optimal Control of Partial Differential
18. Equations, W.Desch, F.Kappel, and K.Kunisch (Eds.), Vol.126 of International Series of Numerical Mathematics, Birkhäuser Verlag, 1998, pp. 15-32.
19. M. Bercovier, M. Engelman. A finite element method for incompressible non-newtonian flows // Journal of Computational Physics, Vol. 36, Iss. 3, 1980, pp. 313-326.
20. E.C.Bingham. Fluidity and Plasticity. New-York: Mcgraw-Hill Book Company,Inc., 1922.
21. A. Bossavit. "Mixed" systems of algebraic equations in computational electromagnetism // COMPEL, Vol. 17, Iss. 1, 1998, pp. 59-63.
22. D. Braess. Finite Elements. Theory, Fast Solvers, and Applications in Solid Mechanics. Cambridge: Cambridge University Press, 2001.
23. F. Brezzi, M. Fortin. Mixed and hybrid finite element methods. Springer Series In Computational Mathematics, New-York: Springer Vol. 15, 1991.
24. C. Burstedde, O. Ghattas, G.Stadler, T. Tu, L.C. Wilcox. Parallel scalable adjoint-based adaptive solution of variable-viscosity Stokes flow problems // Computer Methods in Applied Mechanics and Engineering, Vol. 198, Iss. 21-26, 2009, pp. 1691-1700.
25. M. Chatzimina, C. Xenophontos, G. C. Georgiou, I. Argyropaidas, E. Mitsoulis. Cessation of annular Poiseuille flows of Bingham plastics // Journal of Non-Newtonian Fluid Mechanics, Vol. 142, Iss. 1-3, 2007, pp. 135-142.
26. E. V. Ghizhonkov, M.A. Olshanskii. On the domain geometry dependence of the LBB condition // Mathematical Modeling and Numerical Analysis, Vol.34, Iss.5, 2000, pp.935-951.
27. U. Christensen, H. Harder. 3-D Convection With Variable Viscosity // Geophysical Journal International, Vol. 104, Iss. 1, 2007, pp. 213-220.
28. P.G.Ciarlet. Mathematical Elasticity, Vol 1: Three Dimensional Elasticity. Amsterdam: North Holland, 1998.
29. E. J. Dean, R. Glowinski. Operator-splitting methods for the simulation of Bingham visco-plastic flow // Chinese Annals of Mathematics, Vol. 23, Iss. 2, 2002, pp. 187-204.
30. E. J. Dean, R. Glowinski, G. Guidoboni. On the numerical simulation of Bingham visco-plastic flow: Old and New results // Journal of Non-Newtonian Fluid Mechanics, Vol. 142, Iss. 1-3, 2007, pp. 36-62.
31. Y. Dimakopoulos, J. Tsamopoulos. Transient displacement of a viscoplastic material by air in straight and suddenly constricted tubes // Journal of Non-Newtonian Fluid Mechanics, Vol.112, Iss. 1, 2003, pp.43-75
32. F. Duchin, D.B. Szyld. Application of sparse matrix techniques to inter-regional input-output analysis // Economics of planning, Vol. 15, No. 2/3, 1979, pp. 142-167.
33. R.N.Elias, M.A.D. Martins, A.L.G.A. Coutinho. Parallel edge-based solution of viscoplastic flows with the SUPG/PSPG formulation // Computational Mechanics, Vol. 38, No. 4/5, 2006, pp. 365-381.
34. H.C. Elman, G.H. Golub. Inexact and preconditioned Uzawa algorithms for saddle point problems // SIAM Journal on Numerical Analysis, Vol.31, Iss.6, 1994, pp. 1645-1661.
35. H.C. Elman, D.J.Silvester, A.J.Wathen. Finite Elements and Fast Iterative Solvers: with Applications in Incompressible Fluid Dynamics. Oxford: Oxford University Press, 2005.
36. B.Fischer, A. Ramage, D.J. Silvester, A.J. Wathen. Minimum residual methods for augmented systems // BIT Numerical Mathematics, Vol. 38, No. 3, 1998, pp. 527-543.
37. R. W. Freund. Model reduction methods based on Krylov subspaces // Acta Numerica, Vol. 12, 2003, pp. 267-319.
38. M. Garbey, Yu. V. Vassilevski. A parallel solver for unsteady incompressible 3D Navier-Stokes equations // Parallel Computing, Vol.27, Iss.4, 2001, pp.363-389.
39. R. Glowinski. Numerical Methods for Nonlinear Variational Problems. New-York: Springer, 1984.
40. P.P. Grinevich, M.A. Olshanskii. An iterative method for the Stokes type problem with variable viscosity // SIAM Journal on Scientific Computing, Vol.31, Iss.5, 2009, pp.3959-3978.
41. M. Gunzburger. Finite Element Methods for Viscous Incompressible
42. Flows: A Guide to Theory, Practice, and Algorithms. Academic Press, 1989.
43. W. Leontief, F. Duchin, D.B. Szyld. New approaches in economic analysis // Science, Vol.228, No.4698, 1985, pp.419-422.
44. D. McKenzie. The generation and compaction of partially molten rock 11 Journal of Petrology, Vol. 25, No. 3, 1984, pp. 713-765.
45. E. Mitsoulis, Th. Zisis. Flow of Bingham plastics in a lid-driven square cavity // Journal of Non-Newtonian Fluid Mechanics, Vol. 101, Iss. 1-3, 2001, pp. 173-180.
46. L.N.Moresi, S.J.Zhong, M.Gurnis. The accuracy of finite element solutions of Stokes' flow with strongly varying viscosity // Physics of the Earth and Planetary Interiors, Vol. 97, Iss. 1-4, 1996, pp. 83-94.
47. J. Necas. Les Méthodes Directes en Théeorie des Equations Elliptiques. Paris: Masson, 1967.
48. R.A. Nicolaides. Analysis and Convergence of the MAC scheme I. The Linear Problem // SIAM Journal on Numerical Analysis, Vol. 29, No. 6, 1992, pp. 1579-1591.
49. E.J. 0'Donovan, R.I. Tanner. Numerical study of the Bingham squeeze-film problem // Journal of Non-Newtonian Fluid Mechanics, Vol 15, Iss. 1, 1984, pp. 75-83.
50. K. Ohmori, N. Saito. On the convergence of finite element solutions to the interface problem for the Stokes system // Journal of Computational and Applied Mathematics, Vol. 198, Iss. 1, 2007, pp. 116-128.
51. M.A. Olshanskii Analysis of semi-staggered finite-difference method with application to Bingham flows // Computer Methods in Applied Mechanics and Engineering, Vol. 198, Iss. 9-12, 2009, pp. 975-985.
52. M.A. Olshanskii, J. Peters, A. Reusken. Uniform preconditioners for a parameter dependent saddle point problem with application to generalized Stokes interface equations // Numerische Mathematik, Vol. 105, Iss. 1, 2006, pp. 159-191.
53. M.A. Olshanskii, A.Reusken. Analysis of a Stokes interface problem // Numerische Mathematik, Vol. 103, Iss. 1, 2006, pp. 129-149.
54. T. C. Papanastasiou. Flows of materials with yield // Journal of Rheology, Vol.31, Iss.5, 1987, pp.385-404.
55. I. Perugia. A field-based mixed formulation for the two-dimensional magnetostatic problem // SIAM Journal on Numerical Analysis, Vol. 34, No. 6, 1997, pp. 2387-2391.
56. A. Quarteroni, A. Valli. Numerical Approximation of Partial Differential Equations. Springer Series in Computational Mathematics, Berlin: Springer, Vol. 23, 1994.
57. M.urRehman, T.Geenen, C.Vuik, G.Segal, S.P. MacLachlan. On iterative methods for the incompressible Stokes problem // International Journal for Numerical Methods in Fluids, To appear in 2010; D01:10.1002/fld.2235.
58. N. Roquet, P. Saramito. An adaptive finite element method for Bingham fluid flows around a cylinder // Computer Methods in Applied Mechanics and Engineering, Vol.192, Iss. 31-32, 2003, pp. 3317-3341.
59. L. Saad. Iterative Methods for Sparse Linear Systems, Boston: PWS, 1996.
60. Y. Saad, M.H. Schultz. GMRES: A generalized minimal residual algorithm for solving nonsymmetric linear systems // SIAM Journal on Scientific and Statistical Computing, Vol. 7, Iss. 3, 1986, pp. 856-869.
61. F.N.Shwedow. La rigidité de liquides // Rapports au Congrès de Physique à Paris, Vol. 1, 1900, pp. 478-486.
62. P.J. Tackley. Effects of strongly variable viscosity on three-dimensional compressible convection in planetary mantles // Journal of Geophysical Research, Vol.101, 1996, pp. 3311-3332.
63. R. Temam. Navier-Stokes Equations: Theory and Numerical Analysis. New-York: North Holland, 1979.
64. S. Turek. Efficient Solvers for incompressible Flow Problems: An Algorithmic and Computational Approach. Lecture Notes in Computational Science and Engineering, Berlin: Springer, Vol. 6, 1999.
65. D. Vola, L. Boscardin, J. C. Latché. Laminar unsteady flows of Bingham fluids: a numerical strategy and some benchmark results // Journal of Computational Physics, Vol. 187, Iss. 2, 2003, pp. 441-456.
66. A.J. Wathen. Realistic eigenvalues bounds for the Galerkin mass matrix 11 IMA Journal of Numerical Analysis, Vol. 7, 1987, pp. 449-457.
67. P. Wesseling. Principles of Computational Fluid Dynamics. Springer Series in Computational Mathematics, Berlin: Springer, Vol.29, 2001.
68. X. Xie, J. Xu, G Xue. Uniformly-stable finite element methods for Darcy-Stokes-Brinkman models // Journal of Computational Mathematics, Vol.26, No3, 2008, pp.437-455.
69. S.J. Zhong D.A. Yuen, L.N. Moresi. Numerical Methods for Mantle Convection // Treaties on Geophysics, Vol.7, 2007, pp.227-252.
70. S.J. Zhong, M.T.Zuber, L. Moresi, M.Gurnis. Role of temperature-dependent viscosity and surface plates in spherical shell models of mantle convection // Journal of Geophysical Research, Vol. 105, Iss. B5, 2000, pp. 11063-11082.