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

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

Типичное и наиболее простое уравнение гиперболического типа – это волновое уравнение:

. (3.61)

В связи с этим уравнением будем рассматривать следующие задачи.

1. 
   Задача Коши. В области найти функцию, удовлетворяющую уравнению (3.61), а на прямой y=0 – начальным условиям
   

. (3.62)

2. 
   Смешанная граничная задача. В области найти функцию , которая в этой области удовлетворяет уравнению (3.61), а на границе Г области D при y=0 – начальным условиям (3.62), а при – одному из трех граничных условий:
   

а) условиям первого рода

; (3.63)

б) условиям второго рода

; (3.64)

в) условиям третьего рода

(3.65)

Решение задачи Коши

Выберем прямоугольную сетку, положив

Рассмотрим трехслойный шаблон

Используя взятый шаблон, можно получить разностную схему

, (3.66)

где

Разностная схема, аппроксимирующая задачу (3.61), (3.62) с погрешностью порядка :, (3.70)

где оператор L_h вычисляется так же, как и в (3.66), а

Покажем теперь, как по значениям вычислить значения . В силу (3.66) и (3.70) имеем

(3.71)

Разностное уравнение в схеме (3.70) перепишем в виде

(3.72)

При n=1 по формуле (3.72) вычислим значения ; значения известны в силу (3.71). Затем по (3.72) при n=2 вычисляем значения через уже известные значения и так далее.

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

Пусть уравнение

(3.73)задано в области

.

Будем считать, что к уравнению (3.73) присоединены начальные условия

(3.74) и граничные условия третьего рода

(3.75)

Выберем равномерную прямоугольную сетку, положив

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

и множество Г_h граничных узлов

На множестве уравнение (3.73) и начальные условия (3.74) аппроксимируются разностной схемой вида (3.70). Для замены граничных условий на прямых x= и x= воспользуемся формулами типа (3.68):

(3.76)

где _1n=₁(y), _2n=₂(y_n) и так далее. В заключение заметим, что достаточным для устойчивости по начальным данным разностной схемы является условие

.

18. 
    Общие понятия метода конечных элементов.
    

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

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

При построении дискретной модели поступают следующим образом.

1. 
   В рассматриваемой области фиксируется конечное число точек. Эти точки называются узловыми точками или узлами.
   
2. 
   Значения непрерывной величины в каждой узловой точке считается переменной, которая должна быть определена.
   
3. 
   Область определения непрерывной величины разбивается на конечное число подобластей, называемых элементами. Эти элементы имеют общие узловые точки и в совокупности аппроксимируют форму области.
   
4. 
   Непрерывная величина аппроксимируется на каждом элементе многочленом, который определяется с помощью узловых значений этой величины. Для каждого элемента определяется свой многочлен, но многочлены подбираются таким образом, чтобы сохранилась непрерывность величины вдоль границ элементов. Полином, связанный с каждым элементом, называют функцией элемента.
   
5. 
   Объединение конечных элементов в ансамбль. В этом ансамбле узловые значения неизвестной функции должны быть отрегулированы таким образом, чтобы обеспечить наилучшее приближение к истинному непрерывному распределению. Этот этап приводит к алгебраической системе линейных уравнений относительно узловых значений. Эта система и является моделью искомой непрерывной функции.
   
6. 
   Решение полученной системы, то есть нахождение узловых значений.
   
7. 
   Нахождение значения искомой величины в любой точке области по узловым значениям и функциям элементов.
   

Основная концепция МКЭ может быть наглядно проиллюстрирована на одномерном примере заданного распределения температуры в стержне (см. рис. 9).

Рис. 9. Применение МКЭ к заданному распределению
температуры в одномерном стержне

Разбиение области на элементы может быть проведено двумя различными способами. Можно, например, ограничить каждый элемент двумя соседними узловыми точками, образовав 4 элемента, или разбить область на два элемента, каждый из которых содержит 3 узла (рис. 10).

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

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

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

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

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

Рис. 11а. Элементы в форме треугольника (3 узла)

и четырехугольника (4 узла)

Рис. 11б. Элемент в форме треугольника (6 узлов)
Кроме того, избыточное число узлов позволяет рассматривать элементы с криволинейными границами.

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

19. 
    Дискретизация области и нумерация узлов.
    

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

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

Нумерация узлов – следующая процедура этапа выделения конечных элементов. Порядок нумерации имеет в данном случае большое значение, так как влияет на эффективность метода.

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

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

При нумерации узлов предпочтителен способ, обеспечивающий минимальную разность между номерами узлов в каждом отдельном элементе. Если наибольшую по всей области разность между номерами узлов для отдельного элемента обозначить через R, а число степеней свободы – через Q, то ширина полосы

.

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

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


Рис. 12. Два различных способа нумерации узлов
При первом способе R=14, при втором R=6.

Ширина полосы для этих способов при одной степени свободы в узле получается равной соответственно 15 и 7, а при двух степенях свободы – 30 и 14. Рациональная нумерация в случае б) сокращает объем оперативной памяти примерно в два раза по отношению к случаю а).

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

20. 
    Линейные интерполяционные полиномы.
    

Классификация конечных элементов может быть проведена в соответствии с порядком многочленов − функций этих элементов. При этом рассматриваются три следующие группы элементов:

• 
  симплекс-элементы,
  
• 
  комплекс-элементы,
  
• 
  мультиплекс-элементы.
  

Симплекс-элементам соответствуют многочлены первой степени. Комплекс-элементам – многочлены более высокого порядка.

В симплекс-элементе число узлов равно размерности пространства + 1. В комплекс-элементе число узлов больше этой величины.

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

21. 
    Одномерный симплекс-элемент.
    

Одномерный симплекс-элемент − это прямолинейный отрезок длины L с двумя узлами, по одному на концах отрезка (рис. 13). Узлы обозначаются индексами i и j , узловые значения − Ф_i , Ф_j .

Рис. 13. Одномерный симплекс-элемент

Функция элемента  имеет вид

. (3.77)

Коэффициенты и легко определяются.

При и формула (3.77) дает

.

Аналогично

.

Решая два последних уравнения относительно , получим

, (3.78)

. (3.79)

Подставляя найденные значения ₁ и ₂ в формулу (1), получаем для  выражение

, (3.80)

которое можно переписать в виде

. (3.80^а)

Линейные функции от х в формуле (3.80^а) называются функциями формы или интерполяционными функциями. Будем обозначать их через N_i и N_j :

, . (3.81)

Здесь индексы i и j у N обозначают узел, к которому относится функция формы. Теперь соотношение (3.80) можно записать в матричном виде:

, (3.82)

где − матричная строка и − вектор-столбец.

Как видно из формулы (3.81), функция N_i равна 1 в узле с номером i и равна 0 в узле с номером j. Аналогично функция N_j равна 0 в i-ом узле и 1 в j-ом узле.

Эти значения характерны для функций формы. Они равны 1 в одном определенном узле и 0 во всех других узлах.

22. 
    Двумерный симплекс-элемент.
    

Двумерный симплекс-элемент изображен на рис. 14.


Рис. 14. Двумерный симплекс-элемент

Будем нумеровать узлы против часовой стрелки. Интерполяционный полином имеет вид . (3.83)

Условия в узлах i, j, k приводят к системе уравнений

(3.84)

Решение этой системы дает

,

.

Здесь А − площадь треугольника i j k, которая связана с определителем системы (3.84) следующим образом

. (3.85)

Подставляя значения ₁, ₂, ₃ в формулу (3.83), можно преобразовать выражение для  к виду, подобному (3.82):, (3.86)

где, (3.86^а)

, (3.86^б)

, (3.86^в)

Легко показать, что значение N_i в i-ом узле равно 1, N_i=0 во втором и третьем узлах, а также во всех точках прямой, проведенной через эти узлы.

23. 
    Местная система координат.
    

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

Рассмотрим треугольный элемент, в котором скалярная величина представлена в виде:, (3.86)

где

, (3.86^а)

, (3.86^б)

, (3.86^в)

а функции формы определяются формулами (3.86^а, 3.86^б, 3.86^в).

Поместим начало локальной системы в центре элемента (рис. 15):

t

y

x

Y

j

i

k

s

X
Рис. 15. Местная система координат

Запишем формулы преобразования координат:

(3.87)

Здесь и – координаты центра треугольника

(3.88)

Функция формы N_i после подстановки (3.87) примет вид

. (3.89)

Учитывая теперь (3.86^а) и соотношение (3.88), можно получить, что

.

Таким образом, функция формы в местной системе координат принимает вид

. (3.90^а)

Аналогично получаем

, (3.90^б)

. (3.90^в)

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

, (3.91)

где R и R* − соответственно старая и новая области интегрирования, |J| − модуль якобиана преобразования системы координат, который равен отношению площадей в двух системах координат .Так как обе системы прямоугольные и масштабы измерения в них совпадают, то |J| = 1.

Кроме того, заметим, что формы элементов R и R* сохраняются.

Таким образом, соотношение (3.91) принимает вид

. (3.92)

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