Вестник Сыктывкарского университета. С ер Л. Вып. 13.2011
УДК 004.051+539.3+519.6
ВЫЧИСЛИТЕЛЬНЫЙ КОМПЛЕКС «БИФУРКАЦИОННЫЙ МЕТОД В НЕЛИНЕЙНЫХ МОДЕЛЯХ МЕХАНИКИ»1
Н. А. Беляева, М. Н. Истомина
Вычислительный комплекс объединяет программы по бифуркационным методам в нелинейных моделях механики. В статье рассматривается общая структура комплекса и приводится описание работы входящих в него программ.
Ключевые слова: вычислительный комплекс, бифуркационный метод, куэттовское течение, структурированная жидкость, нелинейная механика, устойчивость, метод прогонки.
1. Общее описание комплекса
Вычислительный комплекс (ВК) «Бифуркационный метод в нелинейных моделях механики», разрабатываемый в рамках проекта «Нелинейные модели и методы механики», объединяет вычислительные программы по нелинейным моделям и методам механики. Программы выполнены на основе математических методов исследования различных критических явлений в нелинейных задачах механики сплошной среды, таких как бифуркационный метод, параметрический анализ, численное решение собственных краевых задач и другие. Целью создания комплекса является объединение в едином программном модуле разработанных на данный момент программ по указанным моделям и методам.
1 Работа выполнена при финансовой поддержке ФЦП "Научные и научно-педагогические кадры инновационной России "на 2009-2013 годы, ГК №02.740.11.0618
© Беляева Н. А., Истомина М. Н., 2011.
Оболочка комплекса реализована в виде web-сайта, выполненного на языке гипертекстовой разметки HTML с использованием каскадных таблиц стилей CSS. Общие правила оформления основной панели и панелей меню прописаны в листе стилей, что значительно упрощает код отдельной веб-страницы и облегчает создание новой страницы. Это немаловажно, т.к. по мере написания новых программ предполагается пополнение комплекса.
ГОУ ВПО "СЫКТЫВКАРСКИЙ ГОСУДАРСТВЕННЫЙ УНИВЕРСИТЕТ" КАФЕДРА МАТЕМАТИЧЕСКОГО МОДЕЛИРОВАНИЯ И КИБЕРНЕТИКИ
ж
ВЫЧИСЛИТЕЛЬНЫЙ КОМПЛЕКС
БИФУРКАЦИОННЫЙ МЕТОД В НЕЛИНЕЙНЫХ ЗАДАЧАХ
СыктГУ
войти кММиК
Работа выполняется при финансовой поддержке ФЦП "Научные и научно-педагогические кадры инновационной России", 2009 - 2013 гг, ГК № 02.740.11.0618
Рис. 1. Главное окно ВК
Комплекс состоит из трёх основных разделов: «О комплексе», «Программы» и «Разработчики».
Раздел «О комплексе» содержит общую информацию и описание структуры комплекса.
Информация, в том числе и контактная, о разработчиках комплекса и авторах входящих в него программ приведена в разделе «Разработчики».
Непосредственный выбор и запуск программ осуществляется в разделе «Программы». Здесь же можно открыть окно с более подробным описанием: авторы, краткая аннотация, постановка задачи, алгоритм, контрольный вариант счёта.
Все программы, входящие в комплекс, были разработаны на кафедре
математического моделирования и кибернетики Сыктывкарского государственного университета. В настоящее время ВК включает в себя следующие программы: «Комбинированный алгоритм перебора вариантов в задачах на собственные значения», «Движение по параметру жесткости в проблеме устойчивости на границе винклеровых сред», «Поиск устойчивого неоднородного решения при течении структурированной псевдопластической жидкости в области сверханомалии», «Бифуркационный анализ куэттовского течения структурированной жидкости в плоском зазоре».
Результаты работы всех программ отличаются наглядностью, что достигается за счёт как численного, так и графического представления результатов вычислений.
Ниже рассматриваются входящие в комплекс программы.
2. Комбинированный алгоритм перебора вариантов в задачах на собственные значения
Проблема устойчивости стержней, пластин и оболочек с односторонними связями относится к конструктивно-нелинейной механике упругих систем. Решение такого рода задач отличается высокой инженерной востребованностью, т.к. полученные результаты могут быть использованы при расчете на устойчивость элементов машин и механизмов.
Задачи на устойчивость одномерных элементов конструкций на границе разномодульных упругих сред винклерова типа при использовании конечно-разностной аппроксимации (на сетке размерностью га) приводятся к решению спектрального уравнения
где А, —ленточные положительно определенные матрицы (га — 1) х (га - 1); й) = • • • wi = & = М, Н = тг/т; Ь =
[&1, • • • ? Ь — вектор формы, такой, что
С = Лс^[с1(А;1&1 + к2( 1 - &х)) ,..., ст_1(£4&ш_1 + к2(1 - Ьт_ 1))] — матрица жесткостей, £4, /^ — безразмерные параметры жесткости; Л —собственное число.
Граничные условия задачи отражены в структурах матриц А и Спектральная задача (2.1) является неопределенной из-за того, что компоненты г Е 1 : га — 1, вектора формы изначально неизвест-
Аи) + Си) = ХС^ш.
(2.1)
(2.2)
ны. Можно использовать следующий алгоритм нахождения первой собственной пары (числа и формы) уравнения (2.1). Задается некоторый вектор формы и решается соответствующая детерминированная спектральная задача. Запоминаются согласующиеся с выбранным вектором формы значения w или — w и отвечающие им собственные числа. После перебора всех возможных значений вектора формы, которые составляют 2m_1 вариантов, выбирается наименьшее собственное число и отвечающая ему собственная форма.
Сформулированный алгоритм полного перебора вариантов (ППВ) всегда приводит к искомому решению. Однако он наталкивается на «проклятие размерности» (curse of dimensionality, R.E. Bellman): при использовании этого алгоритма, например, на сетке размерностью m = 100, потребовалось бы решать свыше Ю30 линейных спектральных задач.
Обойти названное препятствие в реализации алгоритма перебора вариантов позволяет предложенный в работе [3] комбинированный алгоритм «ППВ + ЛПВ» (ЛПВ — локальный перебор вариантов). Сначала на редкой сетке, т.е. такой, чтобы число 2m_1 было не слишком большим, реализуется алгоритм ППВ, и устанавливается отвечающая минимальному собственному числу качественно адекватная собственная форма, т.е. имеющая устойчивый с ростом m вид графика. А затем при последовательном удвоении числа интервалов сетки выполняется перебор вариантов лишь вблизи корней названной формы, т.е. реализуется алгоритм ЛПВ.
В данной программе реализован комбинированный алгоритм перебора вариантов, с помощью которого строится часть собственного спектра конечномерного операторного уравнения (2.1).
Он заключается в следующем:
1. на сетке с относительно небольшим числом узлов применяется к уравнению (2.1) алгоритм полного перебора вариантов возможных форм изгиба и определяется качественно адекватная собственная форма (качественная конфигурация которой остается неизменной);
2. далее последовательно удваивается число узлов сетки путем деления пополам и выполняетсям перебор вариантов лишь вблизи корней искомой собственной формы;
3. процесс продолжается до тех пор, пока соответствующее собственное значение не стабилизируется с требуемой (и достижимой) точностью.
Данная программа реализована на языке С++. На рис. 2 приведен пример работы программы.
Рис. 2. Результаты работы программы «Комбинированный алгоритм перебора вариантов в задачах на собственные значения»
По умолчанию численный результат расчета записывается в файл StergOut.txt. При желании можно задать свое имя выгружаемого файла и путь сохранения. Для вывода графика необходимо, чтобы на компьютере был установлен Borland С++ Builder, придусмотрена также возможность построения графика в Maple.
Подробное описание комбинированного алгоритма поиска части собственного спектра положительного оператора, а также результаты решения задач с использованием данного алгоритма представлены в работах [1]- [3].
3. Движение по параметру жесткости в проблеме устойчивости на границе винклеровых сред
В данной программе реализован алгоритм движения по параметру жесткости одной из винклеровых сред, с помощью которого решается спектральная задача конечномерного операторного уравнения (2.1).
Допустим, что k\ < k2 (для определенности) и к>2 — к\ + fe, í = 0,l,2,..., е — величина шага движения по параметру.
При t — 0 получаем случай однородной упругой среды с параметром жесткости к — к\ — к®. В этом случае решаем задачу на собственные значения детерминированного (все компоненты вектора формы известны) уравнения (2.1) и выбираем минимальное собственное значение,
а соответствующую ему собственную форму рассматриваем как качественно адекватную.
Пусть при к2 = /4_1 собственный спектр найден. Найдем при к2 = кг2. Для этого применим алгоритм локального перебора вариантов, т.е. выполним перебор вариантов вблизи корней качественно адекватной собственной формы (т.е. найденной при к2 — к
- перебираются возможные представления вектора формы Ь = [&ь ...,6Ш_1]Т, где
1,и)1 > 0 0,^0'
- для каждого варианта вектора формы решается задача на собственные значения детерминированного уравнения (2.1);
- запоминается собственная пара (число и форма), для которой форма изгиба ъи^ или —ио^ согласуется с выбранным вектором формы. Найденная собственная форма, используется как качественно адекватная.
Увеличиваем £ до тех пор, пока к2 — к\ не достигнет требуемого значения.
||± Алгоритм движения по параметру :
Параметры задачи
Жесткость среды сверху |до Жесткость среды снизу |-|б
Число точек сетки (п) р[00 Шаг движения по параметру | Граничные условия ¡Шарнирное опирание
Результаты расчета-
Критическая нагрузка
19.72339820540716
Кр.знач. в однор.упр. среде 18.008596360737 Имя выходного Файла
П1е0и( М
|рмы потери устойчивости
' В однор. среде
Вычислить | Прервать | |: График Отмена |
Рис. 3. Результаты работы программы «Движение по параметру жесткости в проблеме устойчивости на границе винклеровых сред»
Данная программа реализована на языке С++. На рис. 3 приведен пример работы программы.
По умолчанию численный результат расчета записывается в файл FileOut.txt. При желании можно задать свое имя выгружаемого файла
и путь сохранения. Для вывода графика необходимо, чтобы на компьютере был установлен Borland С++ Builder, придусмотрена также возможность построения графика в Maple.
Подробное описание алгоритма движения по параметру жесткости в проблеме устойчивости на границе винклеровых сред, а также результаты решения задач с использованием данного алгоритма представлены в работе [4], принятой к печати.
Алгоритм движения по параметру жесткости одной из разномодуль-ных упругих сред является альтернативным комбинированному алгоритму.
Основным отличием алгоритма движения по параметру жесткости от комбинированного алгоритма является возможность получения значений критической нагрузки сразу для заданного интервала. Например, если мы вычисляем критические силы для к\ — 90 и к2 = 16 с шагом 2, то после выполнения программы мы имеем значения критических сред для следующих параметров (£4, к2): (16, 16), (18, 16), (20, 16), ..., (86, 16), (88, 16), (90, 16). Это удобно, если необходимо исследовать поведение системы для большого числа близких друг к другу значений параметров жесткости упругих сред.
Комбинированный алгоритм не дает возможности задать фиксированное значение числа узлов сетки, но удобен, если необходимо вычислить критическое значение сжимающей силы для одного конкретного значения параметров жесткости.
4. Поиск устойчивого неоднородного решения течения структурированной псевдопластической жидкости в области сверханомалии
Неньютоновские жидкости образуют чрезвычайно широкий класс разнообразных материалов, единственными общими свойствами которых являются их текучесть и отклонение от закона трения Ньютона. Обширный класс неньютоновских жидкостей образуют структурированные текучие системы. Структурные превращения в таких системах оказываются ответственными за процессы самоорганизации - спонтанное появление и развитие некоторой структуры в первоначально однородной среде. К таким структурам относятся автоколебания, автоволны, диссипативные структуры и т.д. Под диссипативными структурами понимается качественно новый вид решений, являющихся примером способности неравновесности служить источником упорядоченности. Примером образования диссипативных структур является химическая реакция Белоусова-Жаботинского; диссипативные структуры воз-
никают в процессе клеточного метаболизма и т.д.
Рассматриваемая программа комплекса позволяет проводить исследования куэттовского течения структурированной псевдопластической жидкости в плоском зазоре при потере устойчивости однородного стационарного состояния в области немонотонности реологической кривой (области сверханомалии).
Программа выполнена в среде программирования Delphi. Основной метод расчета — метод прогонки. Результатом работы программы является построение устойчивого неоднородного решения — диссипативной структуры — рис. 4.
Рис. 4. Окно программы «Диссипативная структура в области
сверханомалии»
Пусть структурированная вязкая несжимаемая жидкость заполняет полосу между плоскостями ( = 0 и ( = а ее течение происходит в направлении оси г] и вызывается движением плоскости £ = к. Обозначим через и — вектор скорости жидкости в произвольной точке полосы в момент времени тогда й = (0,г^,0), причем щ — £). Будем предполагать, что жидкость является смесью двух компонент А1 и А2,
которые под действием приложенного механического поля, характеризуемого напряжением а и скоростью деформации 7 = могут взаимно превращаться друг в друга. Суммарная скорость превращения определяется соотношением:
Ф(а, 7) = к2( 1 - а) - а£;0 ехр (р0а + д072) ,
где а — степень структурных превращений (доля А\ в смеси А\ и А2), к2, р, ц — параметры жидкости. Будем предполагать, что жидкость подчиняется реологическому уравнению состояния вида
о- = М7,
где вязкость ¡л = р(а) задается соотношением ц~1 (а) = — а),
здесь /¿1? /¿2 — вязкости компонент Ах и А2 соответственно.
Модель течения жидкости описывается с помощью диффузионно-кинетического уравнения относительно степени структурных превращений а и уравнения движения:
да п<Э2а
т = дё
Ф (а, 7)
д ( , .ди\
Здесь р — плотность, И — коэффициент диффузии. Граничные и начальные условия рассматриваемой задачи имеют вид:
да
Ж
е=о
да
Ж
= 0, щ=о = 0, щ > 0,
{=А
£ — 0 . ^^о^^^с/?/ — ^ ^и
= 0,
= Ио-
Введем безразмерные переменные и параметры:
- 1 - й - ЕЕ. - рк2Ь2 _
— I'7" - 7 9 1 Р — ' ^ — ' X — 7 ;
п р1г1 ¡л2 Ц 2 к2
ц{а) Ро1л2Щ «о?о тл и
¡/(а) =-= —-—= ——, V = —.
П П Щ
Получим безразмерную модель: /
(4.3)
начальные условия:
г = 0 : = 0. alo^i = «о
(4.4)
граничные условия:
да да
(4.5)
Для численного решения системы (3) - (5) применяется метод прогонки на простраственно-временной сетке (х,т).
Формула прогонки для расчета степени структурирования среды:
Ei, Fi- прогоночные коэффициенты. Прямой ход прогонки состоит в нахождении по рекуррентным формулам прогоночных коэффициентов, обратный — вычисление степени структурирования.
Аналогичный метод применяется для нахождения скорости течения.
Таким образом вычисляется скорость течения и степень структурных превращений во всех узлах построенной сетки.
Программа, выполненная в среде программирования Delphi, используется для численного построения устойчивого неоднородного решения — диссипативной структуры, возникающей в результате потери устойчивости однородного стационарного решения задачи в области сверханомалии. Впервые указанная структура была получена в работе [5]. Это — монотонно-возрастающее или монотонно-убывающее решение а = а(х), удовлетворяющее стабилизационному соотношению.
5. Бифуркационный анализ куэттовского течения структурированной жидкости в плоском зазоре
Основной целью данной программы является получение диссипативной структуры в области сверханомалии на основе бифуркационных исследований, сопоставление полученного результата с результатами численных исследований, представленных в работе [5].
ai-i,j — Ei_iaij + Fi-ij
Применение процедуры Ляпунова - Шмидта позволяет определить вид аналитического решения стационарной задачи куэттовского течения структурированной жидкости в плоском зазоре. Полученное решение качественно совпадает с указанным выше численным решением. Метод и результаты аналитических исследований представлены в работах [6] - [8]. В силу громоздкости необходимых математических преобразований они проведены с использованием математического пакета МАРЬЕ. На рис. 5 приведен график соответствующего аналитического стационарного неоднородного решения, качественно совпадающего с дисссипативной структурой [5].
О 0,2 0.4 „ 0.6 0.8 1
Рис. 5. График аналитического стационарного неоднородного решения, полученного посредством бифуркационного метода
Для определения типа бифуркационной точки решаем задачу (3) -(5), рассматривая при этом стационарные решения, т.е. систему приводим к виду:
лч/ \
Ф (а) = к (1 — а — ах ехр [у + 5у2 (1 + Аа)2]) .
(5.6)
Выберем в качестве бифуркационного параметра /3, остальные параметры считаем фиксированными. Существуют критические значения /30, при которых возможны бифуркации неоднородных структур из однородного распределения а^. Из условия обращения в ноль собственного значения оператора для задачи (6) с граничными условиями (4), линеаризованного на постоянном решении ао, находим значение /30 :
Ь(а0,(3)г = /3
дг2
(5.7)
с1Ф
где ^ = — аа
. Получаем:
а=а о
А) = 7—42, П = 1,2,...
При /3, близком к /30, происходит рождение неоднородного стационарного решения.
Бифуркационный анализ состоит в представлении решения вблизи критического значения параметра в виде асимптотического ряда по степеням некоторого параметра характеризующего отклонение бифуркационного параметра от его критического значения, и нахождения коэффициентов указанного разложения. Известный метод Ляпунова-Шмидта определяет рекуррентную процедуру для получения коэффициентов разложения, на каждом шаге которой необходимо решать линейные уравнения с оператором Ь (7) и правыми частями, определяемыми разложением нелинейных членов в уравнении (6) и зависящими от решений, полученных на предыдущих шагах алгоритма. В итоге получаем аналитическое решение в виде асимптотического ряда.
Проведенное исследование позволило сделать вывод о докритиче-ском типе бифуркации [9].
6. Заключение
Таким образом, в ВК собраны программы по нелинейным моделям и методам механики, разработанные на кафедре математического моделирования Сыктывкарского государственного университета. Комплекс позволяет демонстрировать работу программ для решения ряда задач без дополнительных затрат на разработку и написание новых программ. Результаты вычислений отличаются наглядностью представления, возможностью сравнения результатов решения одной задачи разными методами. Вычислительный комплекс дает возможность изучения алгоритма каждой из представленных программ. Все это позволяет использовать данный ВК в учебном процессе. В частности, он может быть использован в рамках таких спецкурсов как «Теория пластин и оболочек», «Устойчивость упругих систем», «Математические модели теории упругости», «Математические модели гидромеханики», «Математические модели термовязкоупругости», «Механика сплошной среды», «Основы гидродинамики» и др. Помимо активного внедрения комплекса в учебный процесс планируется его размещение на сайте Сыктывкарского университета. Предполагается дальнейшее применение комплекса в задачах исследования на устойчивость критических точек в нелинейных
задачах механики, при написании курсовых и дипломных работ, подготовке магистерских и кандидатских диссертаций. В настоящее время работа над комплексом продолжается. По мере появления новых программ комплекс будет пополняться.
Литература
1. Михайловский Е.И., Тулубенская Е.В. Программная реализация комбинированного алгоритма перебора вариантов в задачах на собственные значения. Свид. о госуд. регистрации программы для ЭВМ № 2010615791. 7 сентября 2010 г.
2. Тулубенская Е.В. Устойчивость оболочек и пластин конструктивно-нелинейной механики. Автореферат дис. на соискание ученой степени к. ф. - м. н. СПб: СПбГУ, 2008.
3. Михайловский Е.И., Тулубенская Е.В. Алгоритм локального перебора вариантов в одной существенно нелинейной спектральной задаче //РАН. ПММ. 2010. Т. 14. Вып. 2. С. 299 - 310.
4. Михайловский Е.И., Тулубенская Е.В. Алгоритм движения по параметру жесткости в проблеме устойчивости на границе вин-клеровых сред // Вестн. Сам. гос. техн. ун-та. Сер. Физ.-мат. науки
5. Беляева H.A. Неоднородное течение структурированной жидкости // Математическое моделирование. 2006. Т. 18. № 6. С. 3 -14.
6. Беляева H.A., Кузнецов К.П. Бифуркационный анализ области сверханомалии куэттовского течения структурированной жидкости в плоском зазоре // Материалы Всероссийской научно-методической конференции «Инновации и традиции Науки и Образования■», Часть Сыктывкар: Сыктывкарский государственный университет, 2010. С. 3 - 7.
7. Кузнецов К.П., Беляева H.A. Параметры области сверханомалии куэттовского течения структурированной жидкости // Аналитические и численные методы моделирования естественнонаучных и социальных проблем: Сб. статей V между нар. науч.-техн. конф. 25-28 октября 2010 г. / Под ред. И.В. Бойкова. Пенза, 2010. С. 205 - 209.
8. Беляева H.A., Кузнецов К.П. Бифуркационный анализ куэттовского течения структурированной жидкости в области сверханомалии // Материалы Второй международной конференции
«Математическая физика и ее приложения», Самара, 29 августа - 4 сентября 2010 г. Самара: Изд-во «Книга», 2010. С. 51 -52.
9. Кузнецов К.П. Бифуркация стационарных неоднородных структур в области сверханомалии куэттовского течения структурированной жидкости в плоском зазоре // Материалы Международного научного форума «Ломоносов-2011»/ Отв. ред. А.И. Андреев, А.В. Андриянов, Е.А. Антипов, М.В. Чистякова. [Электронный ресурс] - М.: МАКСПресс, 2011.
Summary
Belyaeva N. A., Istomina М. N. Computing Systems «Bifurcation method in nonlinear models Mechanics»
Computing system includes programs for the branching method in nonlinear mechanics models. The article discusses the general structure of the complex and a description of its constituent programs. Keywords: computing system, bifurcation method, couette flow, structured liquid, non-linear mechanics, stability, sweep method.
Сыктывкарский государственный университет Поступила 20.04-2011