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

ОЦЕНКА ДИНАМИЧЕСКОЙ ПОГРЕШНОСТИ РАЗНОСТНЫХ СХЕМ ЧИСЛЕННОГО МОДЕЛИРОВАНИЯ В РЕШЕНИЯХ ЗАДАЧАХ МЕХАНИКИ ТВЕРДОГО ТЕЛА

THE DIFFERENCE SCHEMES DYNAMIC ERROR ESTIMATION OF NUMERICAL MODELLING IN A SOLID BODY MECHANICS

TASKS SOLUTIONS

В.Н.Савостьянов, В.В.Немчинов, М.С.Хлыстунов, Ж.Г.Могилюк

V.N.Savostjanov, V.V.Nemchinov, M.S.Hlystunov, J.G.Mogilyuk

МГСУ

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

In article are presented a dynamic errors sources researches results at rods oscillations numerical modelling by a finite difference method. The error dependence from a wave period and a grid step are presented.

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

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

Важный этап развития метода конечных элементов связан с решением задач космических исследований в 1950-х годах, хотя идея МКЭ была разработана советскими учёными в 1936 году, но из-за неразвитости вычислительной техники в довоенный период метод не получил дальнейшего развития.

Существенный толчок в своём развитии МКЭ получил в 1963 году после того, как было доказано, что его можно рассматривать как один из вариантов распространённого в строительной механике метода Рэлея-Ритца, который путём минимизации потенциальной энергии сводит задачу к системе линейных уравнений равновесия.

Область применения МКЭ значительно расширилась, когда было установлено (в 1968 году), что уравнения, определяющие элементы в задачах, могут быть легко получены с помощью вариантов метода взвешенных невязок, таких как метод Галёркина или метод наименьших квадратов. Это сыграло важную роль в теоретическом обосновании МКЭ, так как позволило применять его при решении многих типов дифференциальных уравнений.

В последние десятилетия МКЭ занял ведущее положение и получил широкое применение. Он реализован в больших универсальных компьютерных пакетах программ, которые имеют широкое применение.

Первые программные комплексы, в которых реализован метод конечных элементов, были разработаны в 60-х годах. К ним относятся STRUDL-II, SAP-IV, NONSAP, ASKA, NASTRAN, SESAM-69 и другие.

Начиная с конца 70-х годов в СССР появилось несколько десятков программных комплексов, в которых был реализован МКЭ, в том числе: МИРАЖ, МОРЕ, КАСКАД-2, ПРОЧНОСТЬ-75, МКЭ/20, МАРС, ПАРСЕК, ЛИРА, СПРИНТ, FEA. В последние десятилетия наибольшее распространение из них получили: ABAQUS, ADINA, ASKA/DYNAN, ANSYS, MARC, SCAD, MSC/NASTRAN, NEW NASTRAN, Plaxis v9, Plaxis 3D Foundation, EUFEMI, COSMOS, HERCULE, MODULEF, SAP-7, LS-DYNA, Multiphysics и другие.

МКЭ является одним из вариационных методов и часто трактуется как метод Ритца.

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

Таким образом, алгоритм решения МКЭ такой же, как и в вариационных разностных методах, в которых для получения разностной системы уравнений применяется один из вариационных принципов или методов взвешенных невязок.

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

Отличие МКЭ от конечно - разностного метода Эйлера заключается лишь в том, что в выбранном конечном элементе (участке) функция не обязательно заменяется линейной зависимостью, а может быть принята в виде полинома более высокой степени, что повышает порядок аппроксимации и позволяет при той же разбивке на элементы получать более высокую точность (или при той же точности выбрать в области более крупные элементы).

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

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

Учитывая, что МКЭ сформировался как продукт развития и на теоретической базе конечно-разностного метода Эйлера, при всех своих неоспоримых достоинствах, он обладает и определенными ограничениями, в том числе при решении динамических задач механики твердого тела, когда при оценке текущего напряженно-

1S5

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

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

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

В качестве примера, рассмотрим распространение волны в полубесконечном стержне на рис.1, с выбранной базовой длиною Ь, и со скоростью продольной волны ср, которое описывается одномерным безразмерным волновым уравнением, имеющим вид:

д 2 и д2 и

^Ъг' (1)

где и - смещение, 1=срГ/Ь, и=иУЬ, х=х7Ь

Рис.1. Схема распространения волны в полубесконечном стержне

Данное уравнение имеет общее решение

и = /± х) ,

где У - заданные смещения на левом конце стержня.

Любую гладкую функцию, например f, с хорошей точностью можно представить в виде конечного ряда Фурье

N

У^) = X [Л С08О^) + вк )] (3)

к=1

Внутри стержня бегущая волна (слева направо) будет иметь

N

УЦ ± х) =ХН ^ - х)] + Вк 8т[^к^ - х)]}[ ] (4)

к=1

Для исследования частотных характеристик динамических разностных схем применяют комплексную Фурье-моду в'1 ]кс) , где □ - круговая частота, а к - волновое число. Если подставить Фурье-моду в уравнение (1), то она в силу (2) не изменится, но разностное уравнение (аппроксимация дифференциального) меняет как амплитуду (диссипация), так и скорость распространения данной моды (дисперсия). В силу зависимости решения от времени, ошибка так же будет изменятся во времени.

Аналитический метод даёт, в основном, условия устойчивости, но более наглядный, так как это прямой метод получения численного решения, которое основано на

задании синусоиды на левом конце стержня и можно провести оценку изменения амплитуды заданной частоты и скорости её распространения. Рассмотрим разностную схему:

Ч = Uxx + 0,51^х* е

, где и - скорость, u - смещение, - параметры системы.

ut =(1 ~Z)v + Z<€

-тл - i (ш ± kx)

Если подставить комплексную Фурье-моду e в разностную схему, то по-

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

¿a,M = 1 + 2ху2 sin2 ^kh / 2) ± i2rsin(kh / 2^1 -0,25(4 + 2xf r2 sin2 (kh / 2) Откуда следует

2r sin (kh / 2)^1 - 0,25 (4 + 2^)2 r2 sin2 (kh / 2) " 1 -(£ + 2^)r2 sin2 (kh/2)

e"lyM = 1 + 2 [2 - + 2x)\ r2 sin2 (kh /2)

Выражения показывают, что при r<1, ^+2>2 разностная схема устойчива. Зададим шаг по дискретной сетке h=1, тогда минимальная синусоида должна опираться на 5 точек (см. рис.2).

Если задавать Фурье - моду в виде

sin(®t) = sin(2^vt) = sin ^ j,

где □ - круговая частота, v - частота, T - период.

Рис. 2. Один период бесконечных колебаний с максимально возможной частотой согласно теореме Котельникова для Т=5 и при Ъ=Д1=1.

Нагружая подобными импульсами левый конец стержня и рассматривая, например, расчётную область ¿=1000 единиц, проведём расчёт по времени до ¿=1250, т.е. с учётом распространения волны за условную границу стержня, используя граничные условия не отражения, то есть стабилизации волнового процесса во времени, получим решение убедительно иллюстрирующее на рис.3 виртуальную (реально не существующую) «диссипацию» по амплитуде Аа динамических смещений и «дисперсию» по скорости распространения данной моды конечно-разностной модели динамического процесса.

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

Ao

0 25 50 75 100

Рис. 3. Виртуальная затухающая волновая синусоида (выделена пунктиром) на

фоне реального волнового процесса

t

%

Рис.4. Диаграммы роста динамической погрешности при уменьшении периода (повышении частоты) и увеличении времени распространения сигнала.

Кривые диаграммы на рис. 4 подтверждают, что конечно-разностными методами достаточно точно моделируются только длинноволновые приближения.

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

Литература:

1. Зенкевич О. Метод конечных элементов в технике. -М.: МИР, 1975

2. Ильин В.П., Карпов В.В., Масленников A.M. Численные методы решения задач строительной механики. - Минск.: Вышэйшая школа, 1990

3. Теличенко В.И., Хлыстунов М.С. Теория модального анализа микросейсмических процессов и моделирование геодеформационных процессов в основаниях. - Москва-Иваново: Вестник отделения строительных наук РААСН, том 1, 2010, стр.178

Literature:

1. Zenkevich O. Metod konechnyh elementov v tehnike. -M.: MIR, 1975

2. Il'in V.P., Karpov V.V., Maslennikov A.M. Chislennye metody resheniya zadach stroitel'noi mehaniki. - Minsk.: Vysheishaya shkola, 1990

3. Telichenko V.I., Hlystunov M.S. Teoriya modal'nogo analiza mikroseismicheskih processov i modelirovanie geodeformacionnyh processov v osnovaniyah. - Moskva-Ivanovo: Vestnik ot-deleniya stroitel'nyh nauk RAASN, tom 1, 2010, str.178

Ключевые слова: строительная механика, механика твердого тела, напряженно-деформированное состояние, динамические нагрузки, волновые процессы, численное моделирование, метод конечных элементов, конечно-разностный метод, динамическая погрешность, оценка

Keywords: building mechanics, a solid body mechanics, intense-deformed state, dynamic loadings, wave processes, numerical modelling, finite element method, finite difference method, dynamic error, estimation

129337, Москва, Ярославское ш.26, тел.769-73-87, [email protected]

Рецензент: д.т.н., проф. Николаев В.П., зам. научного руководителя ОАО «НИИЭнергетических сооружений Росгидро», зам. директора НТЦ Сооружений, конструкций и материалов