Решение систем линейных уравнений

Исходными данными для решения этой задачи являются:

-форма пластины (двумерная область);

-сетка, покрывающая данную область;

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

стины (30°C, 40°C, 50°C, 80°C, 100°C или 150°C).

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

Все величины представляются в системе измерения СИ.

Все величины – вещественного типа.

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

	∂ 2U	+	∂ 2U	= 0 ,	(1.22)
	∂x2		∂y 2

Этому уравнению удовлетворяет распределение температуры во всей внутренней области пластины, а распределение температур на границах пластины задано граничными условиями (рис. 1.11).

Для вычисления значений температуры T во внутренних узлах пластины (в точках 1, 2, 3, 4, 5, 6, 7) достаточно записать уравнения для каждого узла как среднее арифметическое значений T в соседних точках по горизонтали и вертикали:

31

		=	50 + 80 + 50 + T2
	T1
			4
				T1 + T3 + 45 + 90
		=
	T2
			4
		= 40 + 100 + T2 + T4
	T
	3			4
				T3 + T5 + 35 + 100
	T4	=
			4
		=	T4 + T6 + 30 + 100
	T5
			4
				T5 + T7 + 100 + 100
		=
	T6
			4
		= T6 + 150 + 100 + 100
	T
	7			4

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

4T1 − T2 = 180
− T			+ 4T	2		− T	= 135
	1					3
− T2			+ 4T3			− T4 = 140
			+ 4T4			− T5 = 135
− T3								(1.23)
− T			+ 4T			− T	= 130
		4			4		6
− T5 + 4T6						− T7 = 200
			+ 4T7 = 350
− T6

В матричном виде система линейных уравнений (1.23) запишет-

ся

AT = b ,

где

32

4		− 1	0	0	0	0	0
− 1 4			− 1 0		0	0	0
0 − 1 4 − 1 0						0	0
A =	0	0	− 1 4 − 1 0				0		,
	0	0	0	− 1 4 − 1			0
	0	0	0	0	− 1 4		− 1
		0	0	0	0	− 1	4
0

T1

T2

T3 T = T4 ,

T5T6

T7

180

135 140

b= 135 .

130200

350

Это и есть математическая модель задачи, ее решение сводится к решению системы линейных уравнений (1.23).

Задача 2. По длинной квадратной трубе с квадратным отверстием течет горячая жидкость. Наружный размер трубы равен 8 см, а внутренний – 4 см.

Нижняя половина трубы помещена в ледяную ванну и имеет постоянную температуру 0ºC. Верхняя плоскость имеет постоянную температуру 100ºC.

Температура наружной поверхности трубы на участке между ледяной ванной и верхней плоскостью трубы изменяется линейно:

T (0, y) = T (8, y) = 100( y − 4) (от 0O С до 100O С). 4

Жидкость внутри трубы имеет температуру 200ºC.

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

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

33

Рис. 1.12. Сечение трубы

Исходными данными для решения этой задачи являются:

-наружный и внутренний размер трубы,

-температура нижней и верхней границы трубы,

-температура жидкости,

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

Результатом является распределение температуры внутри тру-

бы.

Все величины представлены в системе измерения СИ. Все величины – вещественного типа.

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

	∂ 2T	+	∂ 2T	= 0.	(1.24)
	∂x 2		∂y 2

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

34

-нижняя половина трубы имеет постоянную температуру 0ºC;

-верхняя половина трубы имеет постоянную температуру 100ºC;

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

Составим математическую модель задачи. Для вычисления температурного режима для заданного сечения трубы достаточно разбить его на n интервалов по направлениям x и y (рис. 1.13). Шаг сетки бу-

дет равен h = 8 . n

Координаты узлов	сетки:	x = x0 + ih ,	(i = 0, 1,K, n) ;
y = y0 + jh, ( j = 0, 1,K, n) .	Для простоты обозначим температуру

в нижнем левом узле на границе области как T0,0 , а температуру в точке с координатами (x0 + ih , y0 + jh ) как Ti, j .

Рис. 1.13. Сетка, покрывающая область распространения тепла

35

Если n = 8 , то шаг сетки равен 1, а граничные условия запишут-

ся так:
Ti,0	= 0;	если i = 0,1,K, 8;
		0, если j = 1, 2,K, 4
T0, j = Tn, j = 100				( j − 4), если	j = 5,K, 8 (линейный закон)
		4
Ti,8	= 100;	если i = 1, 2,K, 8			(1.25)
Ti,2	= Ti,6	= 200;		если i = 2, 3,K, 6
T2, j = T6, j = 200;				если j = 2, 3,K, 6

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

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

	T	=	1	(T	+ T	+ T	+ T	)	(1.26)
	i, j	4		i −1, j	i, j −1	i +1, j	i, j +1
для i = 1, 2, K, 7;		j = 1;
i = 1, 2, K, 7;		j = 7;
i = 1;	j = 2, 3, 4, 5, 6;
i = 7;	j = 2, 3, 4, 5, 6.

Если в этих уравнениях учесть граничные условия (1.25), то получим систему из 24-x линейных уравнений (по количеству внутренних узлов), имеющую единственное решение.

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

36

Замечание. Труба в задаче 3 может иметь разное сечение, например, такое как на рис. 1.14. Алгоритм построения математической модели не изменится.

Рис. 1.14. Сечение трубы

37

2.АЛГОРИТМЫ МЕТОДОВ РЕШЕНИЯ СИСТЕМ ЛИНЕЙНЫХ УРАВНЕНИЙ

Способы решения систем линейных уравнений делятся на две группы:

1.прямые методы, представляющие собой конечные алгоритмы для вычисления корней системы (решение систем с помощью обратной матрицы, правило Крамера, метод Гаусса и др.),

2.итерационные методы, позволяющие получить решение системы с заданной точностью путем сходящихся итерационных процессов (метод итерации, метод Зейделя и др.).

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

	a11x1 + a12 x2 + K+ a1n xn = b1
	a21x1 + a22 x2 + K+ a2n xn = b2	(2.1)
	K K K K K K K
	an1x1 + an2 x2 + K+ ann xn = bn

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

так:

Ax = b ,

	a		a	K a		x
		11	12	1n		1
где A =	a21		a22	K a2n		x2
					, x =	K	,
		K	K	K K
			an2
	an1			K ann		xn

b1

b2 b = K .

bn

2.1.МЕТОД ГАУССА

Метод Гаусса (метод исключения) для решения систем линейных уравнений относится к точным методам. Идея метода Гаусса состоит в том, что система (2.1) путем последовательного исключения неизвестных приводится к системе с треугольной матрицей, из которой и определяются значения неизвестных.

38

Процесс исключения неизвестных состоит в следующем:

Пусть a11≠0. Разделим первое уравнение на a11. Затем вычтем из каждого i–го (i≥2) уравнения, полученного после деления, первое, умноженное на ai1 . В результате, после преобразований x1 окажется

исключенным из всех уравнений кроме первого. По той же схеме исключается x2 (разделив второе уравнение на a22≠0), затем x3 и т.д.

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

~	~	~	~
x1 + a12 x2	+ a13 x3	+ K+ a1n xn = b1
	~	~	~
x2 + a23 x3		+ K+ a2n xn = b2
		~	~	(2.2)
	x3 + K+ a3n xn = b2

KKKKKKKK

~

xn = bn

Особенность этой системы – в строках с номером i все коэффициенты aij при j<i равны нулю. Эту систему уравнений треугольного вида решить уже просто. Из последнего уравнения определяется xn,

далее, подставляя его в предпоследнее уравнение, получаем xn-1				и т.д.
Общая формула определения неизвестных имеет вид
~(i −1)	n	~ (i −1)
	− ∑		x j , i = n, n − 1, K, 1
xi = bi		a ij		(2.3)

j =i +1

Приведение системы (2.1) к треугольному виду (2.2) называется прямым ходом метода Гаусса. Процесс исключения k-го неизвестно-

го называется k-м шагом прямого хода. Элементы a	, a	(1)	, K, a	(n−1)
11		22		nn

называются ведущими.

Общие формулы пересчета коэффициентов системы на k-м шаге имеют вид:

39

~ (k −1)

=

akj(k −1)

j = k , K, n

~ (k −1)

=

bk(k −1)

akj

,

bk

,

a

a (k −1)

(k −1)

kk

kk

k

(k −1)

(k −1)

~ (k −1)

aij

= aij

− aik

akj

,

b

k

= b

(k −1)

− a

(k −1)

~

(k −1)

,

ik

b

k

i

i

i = k + 1, K, n ;

j = k , K, n

Определение неизвестных по формулам (2.3) называется обрат-

ным ходом метода Гаусса.

В методе Гаусса происходит деление строк на соответствующие ведущие элементы, поэтому, если на каком-то k-м шаге на главной

диагонали окажется нулевой элемент akk(k −1) = 0 , то среди элементов

aik(k −1) (i=k+1,…, n) следует найти ненулевой и перестановкой строк

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

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

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

На рис. 2.1 представлена блок-схема прямого хода – исключение i-го неизвестного по методу Гаусса. На рис. 2.2 представлена блоксхема обратного хода – определение неизвестных по методу Гаусса

Примечание. В этом алгоритме А – матрица коэффициентов при неизвестных, В – вектор свободных членов. В самом начале алгоритма необходимо сформировать два массива: двумерный А1 = А и одномерный Х = В.

40