3 Глава Параллельные вычисления при численном решении уравнений математической физики

Глава 1. Параллельные вычисления при численном решении уравнений математической физики.

1. Пример из механики

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

m₁u¨_n
= −k(v_n
^{+ w}n−1
— 2u_n
) − f ⁽¹⁾(t),

n
m₂v¨_n
= −k(w_n
+ u_n
— 2v_n
) − f ⁽²⁾(t),
(1)

n
m₃w¨_n
= −k(u
n+1
+ v_n
— 2w_n
) − f ⁽³⁾(t),

n

n

n

n
представляющей собой результат применения второго закона Нью- тона, где u_n, v_n, w_n суть вертикальные смещения частиц, m₁, m₂, m₃ их массы, k является общей жесткостью пружин, f ⁽¹⁾, f ⁽²⁾, f ⁽³⁾

± ± ±∞
суть внешние силы и n = 0, 1, 2, ..., . Происхождение этой
системы уравнений можно проиллюстрировать с помощью аппрок- симации конечными разностями хорошо известного уравнения ко- лебаний струны (см. [6], [7])

ρu¨ = Tu_xx + F (t, x),
T

mu¨_n =
_h (u_n+1 + u_n−1 − 2u_n) + f_n(t), m = ρh, f_n(t) = hF (t, x_n),
^un+1^−un ₋ ^un^−un−1

u_xx(x_n) ∼
^{h h} , x_n = hn.
h

−
Мы предполагаем, что временная зависимость внешних сил име- ет вид гармонической функции exp( iωt) с частотой ω. Пользуясь трансляционной инвариантностью предположим, что

n
f ^(j)(t) = f_j exp(−iωt + iKLn),
где j = 1, 2, 3, K является «квазиимпульсомk - свободным парамет- ром (−π/L < K < π/L), L - период массива. Мы ищем решение

Эту систему можно решить методом Крамера

2. Пример применения метода Гаусса

Для решения системы (4) мы также можем использовать метод Гаусса - метод последовательных исключений неизвестных. Рас- смотрим пример:

x₁ + x₂ + 3x₄ = 4, 2x₁ + x₂ − x₃ + x₄ = 1,

_
3x₁ − x₂ − x₃ + 2x₄ = −3,
−x₁ + 2x₂ + 3x₃ − x₄ = 4,
Прямые исключения дают




1 1 0 3

_
2 1 −1 1
3 −1 −1 2


−1 2 3 −1

4

1



−
3 _{ ,}

4

1 1 0 3

— − −



_


0 1 1 5

— − −
0 4 1 7
0 3 3 2

4

7



−

−
15 ^{ ,}

8

1 1 0 3

_ − − −
0 1 1 5
0 0 3 13
0 0 0 −13

4

−
7
13
−13
 .

Теперь рассмотрим обратные замены:
−13x₄ = −13, ⇒ x₄ = 1, 3x₃ + 13x₄ = 13 ⇒ x₃ = 0,
−x₂ − x₃ − 5x₄ = −7, ⇒ x₂ = 2,

3. Метод исключения Гаусса в общем случае

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

a₁₁x₁ + a₁₂x₂ + ... + a_1nx_n = f₁, a₂₁x₁ + a₂₂x₂ + ... + a_2nx_n = f₂,

_
.................................................,
a_n1x₁ + a_n2x₂ + ... + a_nnx_n = f_n,

который в матричной форме выглядит так: Ax = f , или

(6)

      ₍₇₎
^_a₁₁ a₁₂ ... a_1n^{ }x₁^





a_n1 a_n2 ... a_nn

x_n

_f1 

f_n

a₂₁ a₂₂ ... a_2n x₂ f₂

... ... ... ... ^{ } ... <sup>

</sup>...^
_ = .

ǁ ǁ
Здесь A = a_ij - матрица с индексами i, j = 1, .., n, x =

/
(x₁, x₂, ..., x_n)^T и f = (f₁, f₂, ..., f_n)^T суть векторы. Вектор x неизве- стен. Предположим, что det A = 0.

≈
Если n имеет большое значение, то методом Крамера потребу- ется много времени для выполнения вычислений, приблизительное количество операций имеет порядок n!. Напротив, метод Гаусса го- раздо быстрее, так как ему нужно приблизительно n³/3 операций. Для примера, если n = 10, тогда 10! = 3628800 в сравнении с 10³/3 333. Таким образом метод исключений Гаусса более эф- фективен при численном анализе.
После n шагов по прямому устранению неизвестных мы полу- чаем систему с верхней треугольной матрицей U



	^u12	^u13	...
	1	^u23	...
	0	1	...
...	...	...	...
0	0	0	...

₁

Сделав следующий шаг, получим

^0 + a

k+1,n

k+1

(k)
k+1,k+1

x_k+1 + ... + a^(k)

x_n = f ^(k) ,

(9)
x_k + u_k,k+1x_k+1 + ... + u_knx_n = y_k,

nn

n
где
(k)

^0 + a

...........................................................,
n,k+1
^xk+1
+ ... + a^(k)x_{n
= f} (k)_.


u_kj =
(k 1)

, y =

_a −
kj

_a −
(k 1) ^k
kk
_f (k−1)

_a −

(k 1)
^k , (10)
kk

a^(k) = a^(k−1) − a^(k−1)u_kj, i, j = k + 1, ..., n, (11)
^{ij ij} ik

i

i

ik
_f (k) _{= f} (k−1) _{− a}(k−1)_{yk, (12)}
и a⁽⁰⁾ = a_1j, j = 1, 2, ..., n, f ⁽⁰⁾ = f₁.
1j 1
Эти формулы представляют собой первую часть алгоритма. Они описывают последовательное исключение неизвестных x_i. Неиз- вестные могут быть легко найдены из полученной системы с по-

≤ ≤ −
для 1 i n 1, если i = n, мы имеем x_n = y_n. Эти формулы обеспечивают вторую часть алгоритма. Используя этот алгоритм, мы можем построить компьютерный код.

4. Трехдиагональные системы

Рассмотрим систему 4 × 4

      <sub>(14)

</sub> = .

a₁ b₂ c₂ 0 x₂ f₂
b₁ c₁ 0 0  x₁

₀ 0 a₃ b₄

_x4

_f1




f₄

0 a₂ b₃ c₃^{ }x₃^
_f3

Это и есть пример трехдиагональной системы. Решение таких си- стем приводит к конкретным простым результатам по гауссовско- му исключению неизвестных. Применим метод Гаусса для решения такой системы. Прямое исключение на каждом шаге приводит к системе следующего общего вида:


1 c₁/d₁ 0 0
 _x1 _y1

_x2 y₂

0 1 c₂/d₂ 0

^{  } = ^{ } , (15)

0 0 1 c₃/d₃^{ }x₃^
_y3

где
0 0 0 1 x₄ y₄

На первом этапе: d₁ = b₁, y₁ = f₁/d₁.

— −

— −

— −
На втором этапе: d₂ = b₂ a₁c₁/d₁, y₂ = (f₁ y₁a₁)/d₂. На третьем этапе: d₃ = b₃ a₂c₂/d₂, y₃ = (f₃ y₂a₂)/d₃. На четвертом этапе: d₄ = b₄ a₃c₃/d₃, y₄ = (f₄ y₃a₃)/d₄.


Наконец, процедура обратной замены дает ответ
x₄ = y₄,

_
x₃ = y₃ − x₄c₃/d₃, x₂ = y₂ − x₃c₂/d₂, , x₁ = y₁ − x₂c₁/d₁.

×
Ясно, что произойдет в общем случае для n n трехдиагональ- ной системы.
На первом этапе: d₁ = b₁, y₁ = f₁/d₁

5.