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

^Э[лг]г Э[ЛГ]Л

Jpc[jV f[W ]rû& -^S+ p i		dr	rdr [T}-	j*[Af]r qyrdr •
A	CTÜ л		dr
	+ '»«[ЛгГ №	} | „ ,	- r (a [ A f /„ =0.		(3.54)
Введем обозначения:
-			\	т
	[М] - матрица демпфирования, [Л/] = рс j[W]			[N]rdr ;
-			d[JV]r э[лгр
	[£] - матрица коэффициентов, [F] = A.j				rdr +
			dr		dr
+ r»a[Aff[JV];
				n>	T
-	{F} - вектор-столбец свободных членов, {F} = J[N] qvrdr ■

С учетом введенных обозначений уравнение (3.54) запишется

как

(3.55)

Производную по времени аппроксимируем по неявной схеме:

Э{Г}	{г’ } - { г - 1}
Эт	(3.56)
	К

Тогда с учетом уравнения (3.56) уравнение (3.55) запишется как

— [М]({г” } - {г""1}) +[х]{г"} ={ f "}

или

1

^

1

[К] + — [М] {Г”} = — [M ]{ r’“1}+ {F m}.

h „ J 1 J AL 1 J i J

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

R j - r	r - R j
[Af]=[jV, Nj] =	L(e)
L{e)

где L{e) —длина конечного элемента.

Для одномерного симплекс-элемента матрица градиентов за­ пишется как

dN,	dNj
M - j f a < 1 “ dr	dr

Поскольку

1

jie) [ - 1 i]

г = N,R + N/R.			и	fL aL,bdx = . a 'M				NZ(e),
'	'	J J				(a +b +1 )!
, f э[лт]г э[лгр
Д	—			rdr = | [i?]r [5] rdr =
		- э Г у			№
		î	- i'
		- î	î	J ( W + t f j R j ) d r =
И				Ri
1	'1	- Г	^ ( 4		+ Яу)	" 1		- Г
к	-1	1		2		" Де) -1		1
J[W ]r		[W]rrfr=			NJfr	NjNj-r	dr =
AW			Z.M NtN jr			N /r
					Л^ЛГЯ,. + JV(wy2 Æy
= 1			Ду					d r =
					N,Nj2+ N /Rj
/<*>
		L(e)	ЗД+Д,		д , + л /
		1 2 [ л , + д у ^ + з д у

	’ у , 2Д +N {NjRj "
J W " * - =	I NtNjRi+N/Rj
I < * >	L<»)
	0	0
	глф Г ] г [ЛГ] = г*а	1
	0

01

W A [ 4 = « V O

1	1
to J * +	____

6

4. МЕТОДЫ РЕШЕНИЯ СИСТЕМ АЛГЕБРАИЧЕСКИХ

УРАВНЕНИЙ

4.1. Метод переменных направлений

Метод переменных направлений - один из лучших методов для решения краевых задач в двумерной постановке. Для записи уравне­ ния (2.34) в разностном виде формируется сетка с узлами в следую­ щих координатах:

хХп= п \ (» = 0 ,1,2 ,...,#), х2т = mhi (« = 0 ,1,2 ,...,АО, ту =УЛО‘= 0 ,1,2 ,...).

Затем на сетке выбирается шаблон, который содержит полуцелый слой т = т + Г |/ 2 , и на нем составляется следующая разностная схема:

K m	- :Unm ) ! Л = ( Л 1«ит + Л 2Мл т + “ « т ) 1 Ъ	( 4 - О
(^л/и	^пт) / Л —(Л]И„т + Л2 йти +сопт) / 2.	(4.2)

Из уравнений (4.1), (4.2) видно, что переход отJ-го временного слоя к (/+ 1 )-му временному слою происходит за два шага по вре­ мени. Сперва при помощи уравнения (4.1) определяются промежу­ точные значения искомой функции йпт. Уравнение (4.1) содержит

три неизвестных мл_, т, ипт, ип+1т ; оставшиеся значения и опреде­ лены на исходном слое. Таким образом, уравнение (4.1) неявно по

координате х, и явно - по х2. При любом фиксированном т урав­ нение (4.1) может быть решено методом одномерной прогонки по направлению хх. Затем схема (4.2), включающая неизвестные ип тА,

ипт, м„ т+1, неявна по направлению х2 и явна по направлению х1. Следовательно, решение системы (4.2) можно получить одномерной прогонкой по направлению х2. Тогда переход с у-го временного слоя на (/+ 1 )-й осуществляется посредством одномерного метода прогонки сперва в продольном направлении, а потом в поперечном. Этим фактом и обусловливается название данного метода.

Каким образом происходит аппроксимация граничных усло­ вий (2.36), (2.37)? Сперва определяются разностные граничные ус­

ловия на целом слое:
«*> =v(*,„.T). ым = х (*,„,т),	(1£»£ЛГ-1).	(4.3)
Затем для прогонки по направлению	х2 требуется определить
значение йпт для п = 0 и п = N. Однако задавать щт=Ц		и ïiNm= ÿ

нецелесообразно, так как значения « не вполне соответствуют мо­ менту т , в этом случае возникает погрешность. По этой причине при аппроксимации граничных условий на полуцелом слое опреде­ ляют промежуточное значение йпт путем вычитания уравнения (4 .2 ) из формулы (4.1):

«лт ~ («/ют «пт) ! 2 —^ 2 («пм ~ «пт)Л! 4 .	(4.4)
Разность в узлах п = 0 и п = N записывается следующим об­
разом:
) / 2 - Л 2 (Ц« -Р » )Л / 4;
«м» = (Y« + Ym) / 2 - Д2 (У. - ут)л / 4,	(4.5)
( 1 < т < М - 1 ) .

4.2. Метод установления

Приведенные ранее экономичные разностные схемы были сформулированы для нестационарных уравнений. Тем не менее в некоторых ситуациях их можно использовать и для решения ста­ ционарных задач. Математическая физика говорит, что возможно рассматривать решение стационарного уравнения как пришедшее к установившемуся состоянию решение соответствующего ему не­ стационарного уравнения при соблюдении условия / —><». В этом случае сам процесс выхода на установившееся состояние сущест­ венной роли не играет и зачастую лишен какого-либо физического смысла, а величина т становится параметром разностной схемы и не является параметром времени. Признаком окончания счета яв­ ляется условие

||^лт ^пт 1 ~ ® •

4.3.Методы решения сеточных уравнений

Система линейных алгебраических уравнений

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

Рассмотрим систему из п уравнений с п неизвестными. В урав­ нении каждый член содержит только одно неизвестное, каждое не­ известное представлено в уравнении в первой степени. Такую сис­ тему принято называть линейной.

Если искомых неизвестных всего два, то каждое уравнение графически изображается линией; три неизвестных - плоскостью;

5х + 7у = 12;

(4.8)

l x + 10y = 17.

Прямые, описываемые системой (4.8), близки друг к другу. Ре­ шением данной системы будет х = 1 , у = 1 , однако при значениях не­ известных х = 2,415, у = 0 получим

(5х +1у = 12,75; [7л + 10у = 16,905.

Эта точка (х = 2,415, у = 0) не будет принадлежать ни одной из описываемых прямых. В таком случае определить численное реше­ ние затруднительно, а его точность будет сомнительной. Система, состоящая из трех и более уравнений, может оказаться плохо обу­ словленной, даже если ни одна из плоскостей не располагается поч­ ти параллельно к другой плоскости.

Условно принято выделять два метода решения систем линей­ ных уравнений: прямые (конечные) и итерационные (бесконечные). Прямые методы могут давать весьма точное решение (с погрешно­ стью округления), если оно существует, посредством конечного числа алгебраических действий. Итерационные методы требуют бесконечного количества алгебраических действий, приводящих в итоге к точному решению. Таким образом, при применении ите­ рационных методов возникает ошибка ограничения, отсутствующая при использовании прямых методов. Однако это не означает, что прямые методы всегда являются наиболее точными. Ситуация опре­ деляется исключительно величинами ошибок округления и ограни­ чения.

4.4.Прямые методы

4.4.1.Метод Гаусса (метод исключения)

Для примера рассмотрим систему, состоящую из трех уравне­ ний с тремя неизвестными:

° \ 1*1 + а12*2 + а13*3 = h '■>	(4.9)
а2Ххх+а22х2 +а2ЪХз=Ь2\	(4.10)
а31хх+а32х2 +а33х3 =Ьз.	(4.11)

Определим множитель т2 как т2 = а2ХI ап , домножим уравне­ ние (4.9) на т2 и вычтем его из уравнения (4.10):

( t f 2 j ~ т 2а \ \ ) х \ ( Л 22 —

"*"(^23 ~ т 2а \ ъ ) х ъ

Первая скобка станет равна нулю, другие переопределим че­ рез а', тогда

ü22x2“Н&23Х3 — *

(4*12)

Заменим уравнение (4.10) на уравнение (4.12). Введем множи­ тель щ = а 3 1 1ап , домножим уравнение (4.9) на щ и вычтем из уравнения (4.11):

К" ™зЩ\)хх +(ап ~ тза\г) х2 + (язз -»*за1з)*з = Ьг -щ Ь х.

Заменим уравнение (4.11 ) на полученное выражение, получится следующая система:

axxxx+ax2x2+ax3x3=bx;		(4.13)
а22Х2	а23Х3= ^ >	(4.14)
&22Х2	*33*3 = *3 •	(4.15)

Новая сформулированная система уравнений (4.13)—(4,15) эк­ вивалентна системе уравнений (4.9)-(4.11), но отличается от нее тем, что в двух последних уравнениях отсутствует член с неизвест­ ным хх. Тогда можно определить неизвестные х2 и % , а хх вычислится в результате подстановки найденных значений в уравне­ ние (4.13). Исключим х2 из уравнения (4.15).

(«32 - »*3«22 )■*2 + («33 - "Ь«23 )*3 = К ~ ” Ч Ь 2 •

Получится новая эквивалентная система

x l + a [2x 2 + a l3x 3 = b l ;	(4.16)
«22*2 + «23*3 =*2 ;	(4.17)

//I//

«33*3 =*3

(4.18)

Такая система называется треугольной. Процесс определения неизвестных существенно упрощается. Сперва вычисляется хг из уравнения (4.18), его значение подставляется в уравнение (4.17) и определяется х2. Затем из уравнения (4.16) по ранее найденным х2

и х3 находится последнее неизвестное х, :

_	Ъ3	_	Ь 2 — & 2 3 Х 3 .		_ ^1 ~ а П Х 2 ~«13*3
•*3	// >	л 2	/	»	л \
	«33			«22	«П

Если «зз = 0, система уравнений является вырожденной. Пример:

x + y +z = 4;

«2х +Зу + z - 9;

х —у —z =- 2 ;

/и2 - 2 ; ( 2 - 2)х + (3 - 2)у +z = 4; => y - z = 1;

x + y + z = 4;

*y	- z	= l;
—2 y —2z —- 6 ;
= -2 ; (-2 + 2)y + (-2 -		2)z = - 6 + 2; => -Az = -4;
x + y + z = 4;
• y - z = l;	=> z= l;y = 2 ; x= 1 .
-4 z = -4;

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

Рассматриваемый метод решения можно обобщить в систему из л уравнений с л неизвестными. Обозначим неизвестные хх,х2,...,х„.

Тогда уравнения запишутся в виде

'аихх+ап х2+... +ах„хп =Ьх; a21*l +а22Х2 + — +а2пХп = ^2>

anixi +a„2x2 + -+ a nnx„=bn.

Введем (и - 1) множителей mi (/ = 1...л):

Вычтем из каждого /-го уравнения первое уравнение, домно женное на /и,. Обозначим

ац =atj -т ,а{]\ b ^ b j - m h ,

где /' = 2 ,... n ,j = 1 ,... л.

Модифицированная система будет выглядеть как

a l l * l + a 12*2 + ‘ " + а \пХп ~ W *

О + Chlx 2 + — + а2пХп = К ;

О + ЛП2 ^ 2 + •••+атхп —Ъп.

Точно так же можно исключить из (и - 2) уравнений неизвест­ ное xi, потом из (л - 3) уравнений хз и т.д. На некотором к-м этапе в процессе исключения неизвестного хк множители будут выглядеть как