книги / Основы математического моделирования и численные методы

T{xx,x2,x) =f { x x,x2),x =Q\		(2.35)
:Г(0,х>,т) = ц(х2 ,г);	T(d>x2>т) = у(х2 ,т);	(2.36)
r(x l50 ,T) = v ( ^ 5T);	Г(х,,£,т) = х(х„т).	(2.37)

Объем вычислений при решении задач (2.34)-(2.37) с помощью явной разностной схемы, если число сеточных неизвестных по каж­ дой переменной равно N, будет пропорционален количеству узлов сетки N 2 Явная схема условно устойчива. Возникают существен­ ные ограничения на шаг по времени, что, в свою очередь, приводит к неоправданно большому объему вычислений N p+2 (р - количество измерений):

r \< l/( tf+ % ) /2 /a .

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

каждом временном слое необходимо решить систему N p уравне­ ний. Для решения этой системы методом Гаусса необходимо про­ вести N 3p~2 действий. При р> 2 расчет по неявной схеме даже ме­ нее эффективен, чем по явной.

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

экономичными.

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

3.МЕТОД КОНЕЧНЫХ ЭЛЕМЕНТОВ

3.1.Суть метода конечных элементов

Метод конечных элементов (МКЭ) - это численный метод, ко­ торый используется для решения дифференциальных уравнений с частными производными задач прикладной физики: механика твердого тела, механика жидкости, электро- и магнитостатика, элек­ тро- и магнитодинамика, тепломассоперенос и др. [3—7]. В МКЭ ис­ следуемая область делится на конечные элементы, в каждом из ко­ торых задается аппроксимирующая кусочно-гладкая функция, кото­ рая за пределами своего конечного элемента равна нулю. В качестве аппроксимирующего полинома очень часто задается полином пер­ вой степени. Решением задачи МКЭ являются значения функций в узлах конечных элементов, которые определяются из решения сис­ темы линейных алгебраических уравнений.

Суть МКЭ может быть продемонстрирована на примере рас­ пределения электрического потенциала ср(х) по толщине диэлек­ трической стенки (область 0Ô), показанном на рис. 3.1 [5]. Иссле­ дуемая область разбивается на конечные элементы (рис. 3.2). Непре­ рывная функция ср(х), определенная в расчетной области, заменяется непрерывной аппроксимирующей функцией и(х) , пред­ ставляющей собой сумму аппроксимирующих кусочно-гладких функций м(Аг) (см. рис. 3.2):

и(*) = 5 У * )

Для рассматриваемого примера количество конечных элемен­ тов равняется четырем, а количество узлов - пяти (см. рис. 3.2). Бу­

дем считать, что значения функции и(х)		в узлах известны и равны
соответственно значениям	Ux, U2,	U5. В пределах своего ко­
нечного элемента функция	представляет собой полином первой

степени, коэффициенты которого вычисляются по значениям функ-

ции в узловых точках Ui и Uj (рис. 3.3), за пределами текущего ко­

нечного элемента	= 0. Функция	представляет собой кусоч­
но-линейную функцию.

Рис. 3.1. Распределение электриче­	Рис. 3.2. Деление области на элементы,
ского потенциала ф по х	узловые точки и аппроксимирующая
	функция и(х)

В общем случае аппроксимирующая кусочно-гладкая функция и(х) заранее не известна. В МКЭ значения в узлах С/,, U2, U3,

должны быть подобраны таким образом, чтобы обеспечивалось «наилучшее» приближение к искомому распределению ср(дг). Для этого используют метод Галеркина, метод Ритца, метод наименьших

квадратов и другие методы [3-8], которые минимизируют некото­ рую величину, связанную с физической сущностью задачи. В ре­ зультате этого задача сводится к решению системы линейных алгеб­ раических уравнений относительно узловых неизвестных и(х).

МКЭ широко используется также при решении двух- и трехмерных задач.

3.2. Дискретизация области

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

3.3. Типы конечных элементов

Простейший конечный элемент имеет число узлов больше на единицу мерности элемента и называется симплекс-элементом. Ко­ нечный элемент, имеющий большее число узлов, называется ком­ плекс-элементом [3, 5].

На рис. 3.4 приведены симплекс-элементы различной размер­ ности.

I

к

а

б

J

в

Рис. 3.4. Симплекс-элементы: а - одномерный; б-двухмерный (треугольник); в-трехмерный (тетраэдр)

Симплекс-элементам соответствуют полиномы, содержащие константу и линейные члены:

- одномерный (см. рис. 3.4, а):

U = OL1+а 2х ;

-двухмерный (см. рис. 3.4, б):

и = а, + а 2х + а Зу;

-трехмерный (см. рис. 3.4, в):

и =а х +а 2х + а Зу + а 4z.

Нелинейные элементы высокого порядка, или комплекс-эле­ менты должны иметь число узлов, равное числу коэффициентов в аппроксимирующем полиноме этого конечного элемента. Например, одномерный квадратичный элемент должен содержать три узла (рис. 3.5), так как его полином содержит три неизвестных коэффи­ циента:

и =ос, + а 2дс + а 3л:2,

a кубический - четыре узла:

и =а, + а 2х +а 3х2+ а 4х3

1	2	3	1	2	3	4
	а				б

Рис. 3.5. Одномерные нелинейные элементы высокого порядка а - квадратичный; б - кубический

5

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

и = а, + а 2х + а гу +а 4х2 +а 5у 2 +а 6х у .

Рис. 3.6. Квадратичный треугольный элемент

Нелинейный двухмерный элемент с наименьшим числом узлов - это четырех­

		угольник (рис. 3.7):
		и =ос, + a 2x + a 3j ; + a 4xy.
	Рис. 3.7. Нелинейный	3.4. Одномерный симплекс-элемент
	четырехугольный	На рис. 3.8 приведен одномерный сим­
	элемент
		плекс-элемент, который имеет два узла /

и j [3, 5]. Длина конечного элемента равняется L , а значения функ­ ции в узлах - соответственно U, и U j.

Начало системы координат расположено вне элемента. Поли­ номиальная функция и для скалярной величины имеет вид

и = а , + а 2х.	(3.1)
Здесь х€ L.
Коэффициенты полинома а х и а 2	определяются через значе­

ния функции в узлах. В результате имеем систему из двух уравнений с двумя неизвестными:

	f£/,.= щ + а 2Х (\
	[UJ =a.i+a.2X j ,
решение которой дает
а ,	UtXj - UjXf	а ,	UJ -и,	(3.2)
	= -------------— ,		= — -------- .
1	L	2	L

С учетом формулы (3.2) выражение (3.1) примет вид

(и,хг и,хЛ [Uj-иЛ

L	) к L	;
которое может быть переписано в виде
и - Г x t		(3.3)

\У

Вформуле (3.3) выражения в скобках называются функциями формы и обозначаются буквой N:

Х . - х	x - X t
N , = ^ —	и Nj =
	L
На рис. 3.9 приведены графики функции формы, из которых
видно, что в узле i функция	= 1, а в узле j - N, = 0. Функция

формы N j —0 в узле /, Nj = 1 в узле j .

Рис. 3.9. Функции формы одномерного симплекс-элемента

Уравнение (3.3) может быть записано как
u =NiUi + NJUJ .		(3.4)
Матричная форма записи уравнения (3.4) имеет вид
»	= [ * М .	<3-5)
где [У] - матрица-строка	функций формы конечного элемента
(матрица функций формы),	[N] = [У, Nj J ; {С/} -	вектор-столбец
узловых значений функции м, {£/} =
Выражение (3.5) является основополагающим		и применимо

к любым конечным элементам. Число элементов матриц [N] и {£/}

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

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

Двухмерный симплекс-элемент показан на рис. 3.10 [3, 5].

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

<Uj = a t + a 2X j + a 2Yj ;

Uk = ax+ a2X k +a 3Yk,

из решения которой получаем

~ X ‘VJ )u <+(x i ï - x & ) u i +(Х Ъ

Нумерация узлов в МКЭ производится против часовой стрелки,

начиная от некоторого произвольно выбранного узла i .
Координаты узлов конечных элементов обозначим	через
(X k,Yk) , а узловые значения функции и -	через
Ult Uj и Uk (см. рис. ЗЛО).
Интерполяционный полином имеет вид
и = а х +а 2х +а 3у .	(3.6)

После подстановки узловых значений функции и и соответст­ вующих координат узлов получаем систему трех уравнений

Г/, = ccj + а 2Х,. + а 3^ ;

(3-7)

<3-8)

щ = ^ [ ( X J ~ х , ) и , Ц х к - x , ) u 1 + ( x ,- х к)и к] . (3.9)

Определитель системы связан с площадью треугольного конеч­

ного элемента А соотношением
1	X,	Yt
1	X J	Yj
1	Х к	Yk

Выражения (3.7)-(3.9) можно записать в виде