Повышение эффективности явной схемы «Крест» численного решения разномасштабных волновых задач механики сплошных сред Текст научной статьи по специальности «Физика»

МЕХАНИКА

УДК 519.6:539.3

ПОВЫШЕНИЕ ЭФФЕКТИВНОСТИ ЯВНОЙ СХЕМЫ «КРЕСТ» ЧИСЛЕННОГО РЕШЕНИЯ РАЗНОМАСШТАБНЫХ ВОЛНОВЫХ ЗАДАЧ МЕХАНИКИ СПЛОШНЫХ СРЕД

© 2010 г. К.А. Кастальская, Д.Т. Чекмарев

Нижегородский госуниверситет им. Н.И. Лобачевского 4ekm@mm.unn.ru

Поступила в редакцию 29.06.2010

Рассматривается метод повышения эффективности явной численной схемы типа «крест» путем введения неявного стабилизирующего оператора простой структуры. Метод применяется в случае одновременного использования крупных и мелких разностных сеток. Рассмотрена его реализация для плоских динамических задач теории упругости. Приводятся примеры численных расчетов, показывающих эффективность метода.

Ключевые слова: явная схема, регуляризация, метод расщепления, разномасштабные процессы.

Устойчивость явных разностных схем типа «крест» решения волновых задач механики сплошных сред определяется условием Куранта - Фридрихса - Леви или максимальной собственной частотой рассматриваемой задачи

[1]. Возникающие при этом ограничения на временной шаг разностной схемы могут быть весьма обременительными и существенно снижать ее эффективность. Они могут быть сняты путем дополнения схемы неявным оператором простой структуры, который, с одной стороны, увеличивает область зависимости сеточной задачи, а с другой, - понижает высшие собственные частоты сеточной задачи [1]. При этом введенный оператор должен легко обращаться и мало влиять на существенные для данной задачи низшие частоты. Данный прием эффективен при наличии в задаче разномасштабных процессов - мелкомасштабного квазистатического и крупномасштабного динамического. В целом, применение стабилизирующего оператора близко по идеологии к методам расщепления по физическим процессам [2].

Предположим, имеется явная разностная схема типа «крест» вида

Чи + Аги = 0,

(1)

где Ьь - самосопряженный положительно определенный сеточный оператор, действующий в пределах одного временного слоя, - опера-

тор второй разностной производной по времени, Бпи = -у (и(1 + т) - 2и(1) + и(1 -т)) (здесь / т

- время, т - шаг по времени).

Предлагается ввести неявный стабилизирующий оператор простой структуры А, после чего схема (1) преобразуется к виду

Ь^и + АБии = 0. (2)

Введение неявного оператора может быть эффективным в ситуациях, когда в процессе деформирования совместно протекают несколько процессов: крупномасштабный волновой и мелкомасштабный «квазистатический», требующий мелкой сетки в локальной области. Это может иметь место, например, при деформировании тонкостенных конструкций или при наличии в деформируемом теле мелких геометрических особенностей, являющихся концентраторами напряжений.

Перед исследованием совместного использования крупных и мелких сеток для двумерных задач имеет смысл рассмотреть одномерную задачу. Рассмотрим одномерное волновое уравнение на отрезке [0, 1] для струны, жестко закрепленной на обоих концах и покрытой крупной и мелкой сеткой. Задана начальная скорость струны:

[1-0082070;, 0< ;< 0.1.

0,

0,1 < ; < 1.

Покроем мелкой сеткой отрезки [0, 0.1] и [0.9, 1], на остальной части струны введем крупную сетку с шагом Н = 0.005, в 5 раз больше чем на мелкой. Шаг по времени будем определять из условия Куранта - Фридрихса -Леви для крупной сетки, а для устойчивости решения в области мелкой сетки введем стабилизирующий оператор. На рис. 1 приведены профили струны в момент времени ? = 0.75 . Расчеты велись на «эталонной» мелкой сетке с шагом по времени т = 0.0001 при числе Куранта кк = 1 и на смешанной сетке с шагами по времени т = 0.0001 при числе Куранта кк = 1 с применением стабилизирующего оператора на той же сетке с шагом по времени т = 0.0005 при числе Куранта кк = 5 .

Использование схемы со стабилизирующим оператором позволило сократить время подсчета по сравнению со схемой, составленной только из мелких ячеек, почти в 20 раз и для схемы с ячейками крупных и мелких размеров - почти в 4 раза. При этом достигнуто качественное и удовлетворительное количественное совпадение результатов. Осцилляции численных решений на смешанной сетке возникают вследствие явления численной дисперсии на стыке мелкой и крупной сеток и не связаны с применением стабилизирующего оператора. Для борьбы с ними можно применить какую-либо из известных процедур сглаживания [3].

Предположим, мы решаем двумерную волновую задачу теории упругости (плоская де-

Рис. 1

-----11=0.001, кк=1

----И=0.005+Іі=0.001, кк=1

----11=0.001+11=0.005, кк=5

формация), которая описывается системой уравнений Ламе:

д 2щ

д 2ы2

(X + ц)

+

д2м1 + д 2Ы2

дх1д^2

дх2

д 2щ

+ц А^1 АМ2 = р дл2 д 2и2

дл2

(3)

венно больше другого }%2 . Шаг по времени определяется условием Куранта относительно минимального размера пространственной ячейки,

Н2

то есть в данном случае условием т < —^, где

с

Здесь и = (м15 Ы2) - вектор перемещений, А,,ц -константы Ламе, р - плотность, А - оператор Лапласа.

Разностная схема решения системы (1), построенная вариационно-разностным методом [4], имеет вид:

(А. + ц)

ао( )Аіиі + ААи2 А2 А1и1 + А0(1) А22и2

+

Х + 2 ц

с =-------- - скорость продольной волны.

В [1] предлагается дополнить сеточную задачу стабилизирующим оператором А в виде А = Е - аО$2, где Е - тождественный оператор, -Ож - оператор, построенный по трем узлам вдоль короткой стороны ячеек разностной сетки и пропорциональный оператору второй частной д

производной

+ Ц Ад и1 = Ра« и1

и2 и2

(4)

д%2

А?2 и]к = и]к+1 - 2и]к + и]к-1 •

(8)

Здесь ^ и £ - операторы первых, а и А2 - вторых производных соответственно по

Хі и %2, и О0} - операторы усреднения

по тем же координатам. Все они определяются следующими формулами:

(А и)ік — 27“ і+1к — и]-1к ),

л(2)

Параметр а определяется из необходимого условия устойчивости разностной схемы в виде

а = ((й1 / ^2 )2 -1)/4 . Введение данного оператора позволяет повысить временной шаг до условия Куранта относительно максимального размера ячейки, то есть в данном случае схема со стабилизирующим оператором

(2),

2 Л

+ М-)

А Аім1 + ААМ2 ААи1 + А( 1 ^22и2

+

1и)}к —

2Л

(иі+1к иі-1к) +Ц иі = РЩі иі

и2 и2

(9)

(Аіи)/к - 7Г (иі+1к - 2и]к + и]-1 к )

(в22и\к - Л2 ^і'+ік

— —;т (и ,-+ік — 2и .к + и

Л2

к -1к

)

(5)

{р01)и. - — (и._1к + 2мд + и]+1к ),

(о02)М. - — (иjk_ + 2ид + и.к+1 )

Сеточный оператор Лапласа в рассматриваемом случае имеет вид:

£д = + ^22 . (6)

Явная схема типа «крест» (2) будет устойчивой при выполнении условия Куранта - Фрид-рихса - Леви

т < т1п(й1, ^2) / с. (7)

В двумерном случае типичными являются следующие ситуации.

Сетка с вытянутыми ячейками. Пусть для численного решения системы применяется сетка, у которой один из шагов сетки Н сущест-

будет устойчивой при выполнении условия Н

т<.

с

Алгоритмически введение стабилизирующего оператора приводит к необходимости решать для каждого ряда узлов систему линейных алгебраических уравнений с трехдиагональной матрицей, что удобно делать методом прогонки.

Совместное использование крупной и мелкой сеток. Предположим, мы решаем задачу, имеющую некоторую геометрическую особенность (выточку, отверстие, жесткое включение и т.п.) малого размера по отношению к геометрическим параметрам как самой конструкции, так и происходящего в ней волнового процесса. Для обеспечения достаточной точности решения в окрестностях особенности, являющейся концентратором напряжений, необходимо ввести более мелкую сетку. Для того чтобы расчеты можно было вести с временным шагом, соответствующим крупной сетке, введем одномерные стабилизирующие операторы для каждого направления.

1^л2

В схеме (9) оператор А зададим в виде А = (Е - Е - «22), где Е - тождест-

венный оператор, Ля и - операторы, построенные по трем узлам вдоль сторон ячеек мелкой разностной сетки и пропорциональные вторым производным по данным направлениям.

В данном случае = и}+1 к - 1и]к + и]-1 к,

2и]к ~ и]к+1 2и]к + и]к-1 •

и2у=0_ 0 и2у=Ь = 0 иі|х=0 = 0

и7\ = 0 ,

х=0

х <е[1,2і,1,8і] ує[ 0,71, і]

=0,

хє[1,2і,1,8і] ує[0,7і, і]

=0 •

Для расчета пластина покрывалась квадратной разностной сеткой 10x30 ячеек. В местах жесткого закрепления каждая ячейка сетки разбивалась более мелкой сеткой 5x5 ячеек, при этом применение стабилизирующего оператора на мелкой сетке позволило оставить неизменным шаг по времени. На рис. 3а приведено зна-

Коэффициенты

соответственно равны: = ((#1 /к1 )2 -1)/4,

а2 = ((И2 /И2)2 -1)/4 (см. [1]). Здесь ИъН2 -шаги крупной сетки, к15 к2 - шаги мелкой сетки.

Рассмотрим задачу о деформировании упругой бесконечной полосы с симметрично расположенными закреплениями (рис. 2). Полоса имеет ширину 3Ь и длину симметричной части Ь, квадратные закрепления размером 0.6Ьх0.6Ь расположены периодически с шагом 2Ь вдоль оси симметрии полосы.

Выделим расчетную область задачи с учетом наличия периодичности и симметрии и сформулируем для нее начальные и граничные условия:

ди

щ\т=0 = 0 и2|т=0 = 0 "Т7Iт=0 =-^0 , дt

чение давления а =

а11 +а

22

2

в момент вре-

Пластина размерами Ьх3Ь с жестким закреп-

х е [1,21,1,8 Ь] лением в области имеет мгновен-

уе[ 0,7 Ь,Ь ]

но приложенную скорость У0 на линии у = 3 Ь , направленную вдоль оси у. Коэффициент Пуассона для материала V = 0.3 .

мени /с / Ь = 1 на крупной сетке 10x30 ячеек с шагом Н = к2 = 0.1 Ь , на рис. 3б - на сетке 50x150 ячеек с шагом к = к2 = 0.02Ь .

Для получения детальной картины в области закрепления Ь < х < 2Ь, 0 < у < Ь вводилась

более мелкая сетка с шагом к = к2 = 0.02Ь , при этом применение стабилизирующего оператора на мелкой сетке позволило оставить неизменным шаг по времени. Граничные условия задавались так, как показано на рис. 4.

На границах стыка крупной и мелкой сетки в совпадающих узлах брались значения из крупной сетки, в промежуточных узлах значения задавались из линейной интерполяции соседних значений крупной сетки.

На рис. 5 изображены линии уровня давления в области, разбитой мелкой сеткой: а) схема с сочетанием крупной и мелкой сетки, б) схема только с мелкой сеткой.

Рассмотрим аналогичную задачу о деформировании упругой бесконечной полосы, геометрия которой приводилась выше, с условием, что полоса вместо закреплений имеет симметрично расположенные отверстия (см. рис. 2). Метод, описанный выше, приводит к неустойчивости задачи в области границ отверстий. Данная неустойчивость связана с реализацией статических граничных условий на краю отверстия. Проведенные исследования показали, что неус-

и

и

1

2

О

возмущение

' г

1

_ 1.1. 4. и —1.,

Ч-М-!- -гггт ■Т-1-И- 41-!-!

:Г1:1ч

^ X ^ ^ - п: т :||й: .ни 1

-ГГ Г 7 ■1Й: I!; хшШ

Ч-Г7-1

■I ! !! ё!!! ш \ !! 4\\Ы

21_

закрепление

значения в узлах берутся из крупной сетки

вычисленные значения передаются в узлы крупной сетки

щ = О

и2 - О

значения в промежуточных узлах строятся как линейные интерполяции ближайших значений крупной сетки

Рис. 4

тойчивости можно избежать, скорректировав ускорения граничных точек значениями ускорений этих точек, взятых из решения с меньшим числом Куранта, распределив это добавочное ускорение на соседние узлы по закону, зависящему от порядка измельчения сетки и основанному на интерполяции Лагранжа. В данной задаче ускорения граничных точек корректировались значениями ускорений соседних точек, взятых из основной сетки (с крупным шагом по пространственной координате). Добавочное ускорение распределялось на соседние узлы мелкой сетки по закону, зависящему от порядка измельчения сетки (см. рис. 4). В данной задаче ускорение в каждом узле граничной области корректировалось по формуле

= а

1 +

п=5

п

мелк кр

а-,] -а7

кр кр

к=0,тФк ащ,] ак

(10)

где а*1’а^у - значения у границ отверстия,

взятые из крупной сетки 30x10.

Вычисления на крупной сетке проводились во всей области пластины, но значения в области, покрытой мелкой сеткой, участвовали только в корректировке ускорений точек на границе, в окончательном расчете не учитывались.

На рис. 6 изображены картины давления в области, покрытой мелкой сеткой: а) схема с сочетанием крупной 30x10 и мелкой 150x50 сетки, б) «эталонное» решение на сетке 300x100 ячеек.

*

:к

'—ОС

а)

х

б)

Рис. 5

На рис. 7 показаны графики изменения давления в точке «максимума» (180, 30) в безразмерном времени т = № / Ь (на оси отложены шаги по времени) для схемы с сочетанием крупной 30x10 и мелкой 150x50 сетки и «эталонного» решения на сетке 300x100 ячеек.

На основе приведенных численных результатов можно сделать следующие выводы.

Решения, полученные на основе схемы со стабилизирующим оператором, позволили сократить объем вычислений для приведенных выше задач почти в 3 раза.

В задаче с периодическими закреплениями удалось получить качественное и удовлетворительное количественное совпадение численных результатов, полученных на основе данной схемы, и эталонного решения на мелкой сетке.

В задаче с периодическими отверстиями можно говорить лишь об удовлетворительном качественном соответствии решений. Для получения приемлемых результатов необходимо решить проблему статических граничных условий в рамках данного подхода.

Схема со стабилизирующим оператором позволяет адекватно описывать длинноволновые динамические процессы (те, которые хорошо аппроксимируются на крупной сетке) и квази-статические мелкомасштабные процессы в областях мелкой сетки. В описание высокочастотных процессов данная схема может вносить значительные искажения, поэтому рамки применимости данного подхода ограничиваются классом задач, в которых высокочастотная составляющая волновых процессов несущественна.

5 ООО 10.00 15 00 КОС

25.00

50.00

;=. со

а)

2L

Ч ■

mv

г

I "V.

б)

И------------------Г-

Рис. б

Рис. 7

В целом, данный подход позволяет значительно повысить эффективность решения данного класса разномасштабных волновых задач механики сплошных сред в рамках явных схем.

Работа выполнена при поддержке РФФИ (грант № 09-01-97019-р_поволжье_а).

Список литературы

1. Баженов В.Г., Пирогов С.А., Чекмарев Д.Т. Явная схема со стабилизирующим оператором

для решения нестационарных задач динамики конструкций // Изв. РАН. МТТ. 2002. № 5. С. 120-130.

2. Марчук Г.И. Методы расщепления. М.: Наука, 1988. 263 с.

3. Самарский А.А., Попов Ю.П.Разностные схемы газовой динамики. М.: Наука, 1980. 351 с.

4. Баженов В.Г., Чекмарев Д.Т. Решение задач нестационарной динамики пластин и оболочек вариационно-разностным методом. Н. Новгород: Изд-во ННГУ, 2000. 118 с.

IMPROVING THE EFFICIENCY OF THE EXPLICIT CROSS SCHEME FOR THE NUMERICAL SOLUTION OF DIFFERENT SCALE WAVE PROBLEMS OF CONTINUUM MECHANICS

K.A. Kastalskaya, D. T. Chekmarev

We consider a method for increasing the efficiency of an explicit numerical scheme of the cross type by introducing an implicit stabilizing operator with a simple structure. The method is applied in the case of simultaneous use of both small-scale and large-scale meshes. We consider its implementation for plane dynamic problems of elasticity theory. Some examples of numerical calculations are given showing the effectiveness of the method.

Keywords: explicit scheme, regularization, split method, different scale processes.