CC BY
42
УДК 532.517,571.938
А. Н. Нуриев
ИСПОЛЬЗОВАНИЕ МЕТОДОВ БИФУРКАЦИОННОГО АНАЛИЗА ПРИ ИССЛЕДОВАНИИ СИСТЕМЫ УРАВНЕНИЙ НАВЬЕ-СТОКСА ДЛЯ ПРИЛОЖЕНИЯ В ЗАДАЧАХ ГИДРОМЕХАНИКИ И ХИМИЧЕСКОЙ ТЕХНОЛОГИИ
Ключевые слова: уравнения Навье-Стокса, бифуркационный анализ, метод продолжения по параметру, метод
гомотопии.
Представлена схема для бифуркационного анализа нелинейных систем большой размерности, возникающих, в частности, при дискретизации стационарных уравнений Навье-Стокса для несжимаемой жидкости. Схема включает методы локализации различных решений системы, построения ветвей и идентификации локальных бифуркаций системы: бифуркаций Андронова-Хопфа и бифуркаций складки.
Keywords: Navier-Stokes equations, bifurcation analysis, continuation method, homotopy method.
The scheme for bifurcation analysis of the large systems of nonlinear equations arising from the discretization of the stationary incompressible Navier-Stokes equations is presented. The scheme includes methods for location of different solutions, continuation of solution branches and identification of local bifurcations of the system: Andronov-Hopf bifurcations and fold bifurcations.
Современные методы бифуркационного анализа нелинейных систем успешно используются при решении задач во многих областях науки. Однако действительно широкое применение они получили для случаев систем небольшой размерности. По сложившимся представлениям, результативность и эффективность большинства бифуркационных методов достигается только в таком классе задач. Это существенно ограничивает их применение для задач гидродинамики и задач химической технологии, которые используют математические модели, основанные на системе уравнений Навье-Стокса. Численные подходы к решению этих задач приводят к нелинейным системам, содержащим сотни тысяч неизвестных. В данной работе представлены методы бифуркационного анализа систем большой размерности, которые возникают при дискретизации стационарных уравнений Навье-Стокса.
Бифуркационный анализ в работе проводиться в рамках классического подхода (описанного например в [1]) для анализа однопараметрических нелинейных систем. Он состоит из следующих этапов:
1. Выбор исследуемого диапазона изменения параметра.
2. Локализация решений при фиксированном значении параметра
3. Построение ветвей решения.
4. Анализ устойчивости и исследовании точек бифуркации.
Представим дискретизованную систему уравнений Навье-Стокса в следующем виде:
F(y) = 0, (1)
где F - нелинейный оператор задачи, y = (x, Re) - расширенный вектор, состоящий из вектора неизвестных задачи x (в зависимости от постановки это дискретные давление и скорости, или завихренность и функция тока) и параметра задачи Re (числа Рейнольдса).
На первом этапе выбор исследуемого диапазона изменения параметра системы Re полностью зависит от конкретной задачи и целей исследования. Следует отметить, что даже небольшой диапазон в области умеренных чисел Рейнольдса может содержать большое количество различных по структуре ветвей решения. Локализация всех этих ветвей важна для практических задач, так как каждая из них может описывать различные варианты развития физических процессов.
Локализация ветвей заключается в поиске возможных вещественных решений системы
(1) при фиксированном значении параметра Re = Re0 из исследуемого диапазона. Для этих целей использовался FPN (Fixed-Point Newton) метод гомотопии [2]. Он объединяет в себе преимущества методов гомотопии Ньютона и гомотопии фиксированной точки и подходит для задач большой размерности. Функция гомотопии H(y, t) строится, как линейная комбинация двух вещественных функций: H(y, t) = t R(y) + (1 - t) G(y). Первая функция, R(y), содержит решения исходной задачи и записывается как R(y) = F(y)(y - y0), вторая задается в виде суммы G(y) = Gfp(y) + Gn(y), где Gfp(y) = (y - y0) - функция, представляющая гомотопию фиксированной точки, Gn(y) = (R(y) - R(y0)) - функция, представляющая гомотопию Ньютона. Здесь t является параметром гомотопии, y0 — некоторое решение задачи. После упрощения задача гомотопии представляется как H(t, y) = (F(y) — 1 + t)(y0 - y) = 0. Эта задача имеет одно тривиальное решение y = y0, не зависящее от параметра гомотопии, и нетривиальные, удовлетворяющие системе уравнений F(y) = 1 - t. При t = 1 решения такой системы являются решениями исходной задачи (1). Таким образом, во-первых, задача локализации сводиться к задаче продолжения по параметру гомотопии t, во-вторых, множество решений исходной задачи может быть найдено с одно начального приближения. В качестве такого приближения в данной работе выбиралось решение главной ветви, являющейся продолжением стоксовского решения задачи при Re = 0. Хотя для общего случая FPN метод не может гарантировать нахождение всех вещественных решений задачи, на практике он позволяет найти много различных решений.
Локализованные при фиксированном числе Рейнольдса ветви необходимо построить во всей исследуемой области. Для этого используется метод продолжения по параметру [1]. Суть метода заключается в построении последовательности решений {y0, y2, ..., ym}, где yk = (xk, Rek) (k = 1, ., m), аппроксимирующих ветвь, которой они принадлежат, в выбранном диапазоне изменения параметра. Для этих целей рассматривается расширенная задача, в которой наряду с вектором xk, неизвестным также считается число Рейнольдса Rek. Основными компонентами метода продолжения являются предиктор и корректор. В качестве предиктора в настоящей работе использовался метод касательной. Начальное приближение для k-ого решения последовательности строилось как yk0 = yk-1 + hkdyk-10, где hk = Rek0 - Rek-1 - шаг по числу Рейнольдса, dyk-1° = (dXk-1°, dRek-10) - вектор касательной, который находится из системы уравнений:
J(yk-1) dyk-10 = 0, dRek-10 = 1. (2)
Здесь J(yk-1) - Якобиан расширенной системы для (k —1) решения последовательности. Коррекция полученного приближения проводилась итерационно с помощью метода Мура-Пенроуза[1], который является модификацией известного метода Ньютона. На каждой итерации метода решалась задача вида
J(yk ”'1) ykn = 0, (ykn - ykn-1, dYk"-1) = 0. (3)
Скобки в последнем уравнении обозначают скалярное произведение, dyk " - вектор
касательной для yk n1, который находился из решения линейной задачи (4).
J(yk-1 n-1) dyk-1 n-1 = 0, dRek-1 n-1 = 1. (4)
Описанный метод обладает высокой скоростью сходимости: суперлинейной в окрестностях точек бифуркаций и квадратичной в других областях.
Для нахождения каждого решения аппроксимирующей последовательности необходимо решать одну линейную задачу (2) для предиктора и две линейные задачи (3),(4) на каждом шаге корректора. Учитывая, что для достижения необходимой точности в среднем требуется 4 итерации метода Мура-Пенроуза, то в сумме для нахождения одного решения нелинейной задачи решается 9 линейных задач. Все эти задачи имеют большую размерность, осложнены отсутствием диагонального преобладания матриц линейных систем. Для их решения использовался метод обобщенных минимальных невязок (GMRES) с предобуславливателем [3],[4]. В качестве предобуславливателя для GMRES использовался многосеточный метод, описанный в [5].
Для идентификации точек бифуркации и анализа устойчивости на четвертом этапе
рассматривалась задача на собственные числа. Устойчивым решениям задачи (1) отвечает отсутствие в спектральной задаче собственных чисел с положительной вещественной частью. В общем случае, потеря устойчивости может произойти в результате одой из следующих ситуаций: либо одно из вещественных собственных чисел перейдет через ноль (бифуркация складки), либо вещественная часть у пары комплексно-сопряженных собственных чисел перейдет через ноль (бифуркация Андронова-Хопфа). Таким образом, анализ устойчивости и идентификация точек бифуркации сводятся к нахождению нескольких собственных чисел задачи с наибольшей вещественной частью. Для этого в данной работе использовался метод ортогональных проекций Арнольди, реализованный в пакете ARPACK [6].
Описанный подход бифуркационного анализа был успешно применен для исследования классической задачи гидромеханики о циркуляционном течении вязкой жидкости квадратной каверне и для исследования стационарной составляющей решения в задаче о гармонических колебаниях цилиндра в вязкой жидкости. В результате были найдены и изучены новые (раннее не опубликованные) решения и построены бифуркационные диаграммы для данных задач.
Литература
1. Kuznetsov, Y.A. Elements of applied bifurcation Theory / Y.A Kuznetsov. - Berlin: Springer, 1995. -591 p.
2. Rahimian, S. K. A new homotopy for seeking all real roots of a nonlinear equation / Saeed Khaleghi Rahimian, Farhang Jalali, J. D. Sea Saeed. // Computers & Chemical Engineering. - 2011. - V. 35, Issue 3. -P. 403-411.
3. Barrett, R. Templates for the solution of linear systems: building blocks for iterative methods / Richard Barrett et al. - Society for Industrial Mathematics, 1987. - 143 p.
4. Ju, J. Predictor-corrector preconditioned Newton-Krylov method for cavity flow / J. Ju, G. Lapenta // Lecture Notes in Computer Science - 2005. - V. 3514/2005. - P.601-666.
5. Егоров, А. Г. Неединственность стационарного течения вязкой жидкости в квадратной каверне / А.Г. Егоров, А.Н. Нуриев // Учeн. зап. Казан. гос. ун-та. Сер. Физ.-матем. науки. - 2009. - Т. 151, № 3. -С.130-143.
6. Lehoucq, R.B. ARPACK users' guide - solution of large-scale eigenvalue problems with implicitly restarted Arnoldi methods. / Richard B. Lehoucq, Danny C. Sorensen, Chao Yang. - SIAM. Software, environments, tools. - 1998. - 142 p.
© А. Н. Нуриев - инж. каф. информатики и прикладной математики КНИТУ, асп. каф.
аэрогидромеханики КФУ, Artem.Nuriev@ksu.ru.
Все статьи номера поступили в редакцию журнала в период с 28.06.11 по 15.08.11.
