3.4. Решение систем линейных уравнений методом квадратного корня

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

Метод квадратного корня применяется для ре­ше­­ния сис­темы линейных уравнений, ко­эф­фи­ци­ен­ты ко­торой об­ра­зуют эрмитову сим­мет­ри­чес­кую мат­ри­цу (эр­ми­то­ва мат­­рица совпадает с ком­­плек­с­­но-со­пря­женной тран­с­­по­ни­ро­ваннойA^*=A).

Представим матрицу Ав видеА=S^*D S, гдеS - верх­няя треугольная матрица с по­ло­жи­тель­ны­ми эле­мен­тами на главной диагонали;S^*- транс­по­ни­ро­ван­ная к матрицеS;D- ди­а­го­наль­ная мат­ри­ца, на ди­агонали ко­торой на­хо­дят­ся чис­ла (+1) или (— 1).

Если матрицы SиD найдены, то заданная сис­­те­маAХ= = F может быть решена следующим пу­­тем:

AX = S^* D S X = (S^* D) S X = B Y = F, (1.28)

где S^*D = Весть нижняя треугольная матрица и Y = S X- вспомогательный вектор.

Таким об­ра­зом, ре­шение системы АX = Fрав­но­­силь­но ре­ше­нию двух треугольных систем ВY = Fи SX =Y.

Пусть S = (s_ik) при i > k и s_ik № 0, s_ii > 0; S^* = (); D = =(d_ik), d =1; i № k. Тогда из срав­не­ния мат­риц АиS^*DSпо­лучим

.

Ограничение в сумме получается из учета то­го фак­­­та, что в матрице Sниже главной ди­а­го­на­ли эле­­менты об­ра­ща­ются в нуль. Последнее ра­венс­т­во можно пе­ре­пи­сать не­сколько иначе, при­няв для оп­ределенностиk= min(k,l):

;

,

откуда окончательно получим формулы для вы­чис­­ления элементов матриц SиD:

(1.29)

Единственным условием возможности опре­де­ле­­ния s_ikявляетсяs_kk № 0. Для по­стро­е­ния мат­­риц по­лагаемk = 1 и последовательно вы­чис­ля­­емвсе эле­менты первой строкиsпо формуле (1.29); за­т­ем по­лагаемk= 2- и определяем эле­мен­ты вто­рой стро­­ки и т.д. Когдаk=n, тогда най­де­ны все эле­­мен­тыматрицS иD,а следовательно иS^*. За­тем последовательно вы­пол­ня­­ем вы­чис­ления (1.28):

S^*Z = F; D Y = Z; S X = Y

обыч­ным ходом по формулам

(1.30)

при i = 2, 3, ..., n.

Попутно заметим, что определитель матрицы А можно вычислить из выражения

. (1.31)

Этот метод экономичен, требует n³/3 ариф­ме­­тических действий и при большихn(n> 50) вдвое бы­стрее метода Га­усса. Если за основу про­цедуры при­нять алгоритм П5.6 Дья­конова [1987], то тогда метод квад­ратного кор­ня мо­жет быть ре­а­­ли­зован c помощью сле­дующей про­це­ду­ры:

procedure KVK (M, N : integer;

aa: array of real; bb : array of real;

var C: array of real);

type mas1=array [0..4] of real;

mas=array [0..4] of mas1;

var a : mas; j, k, i : integer; s, c1 : real;

begin

i := 0;

for j := 1 to m do

for k := 1 to n do

begin Inc (i); a[j,k] := aa [i]; end;

for j := 1 to N do

begin

for k := j to N do

begin

s := 0.0;

for i := 1 to M do s := s + a[i,j] * a[i,k];

c[k] := s;

end;

c1 := 0.0;

for i := 1 to M do с1 := c1 + a[i,j]*Bb[i];

for i := j to N do a[i,j] := c[i];

c[J] := c1;

end;

a[1,1] := sqrt (a[1,1]);

for j := 2 to N do a[1,j] := a[j,1] / a[1,1];

for i := 2 to n do

begin

s := 0.0;

for k := 1 to i-1 dо s := s + a[k,i]*a[k,i];

a[i,i] := sqrt (a[i,i] - S);

for j := i+1 to N do

begin

s := 0.0;

for k := 1 to i-1 do

S := S + A[K,I]*A[K,J];

A[I,J] := (A[J,I] - S) / A[I,I];

end;

end;

c[1] := c[1] / a[1,1];

for i := 2 to N do

begin s := 0.0;

for k := 1 to i-1do

s := s + a[k,i] * c[k];

c[i] := (c[i] - s) / a[i,i];

end;

c[n] := c[n] / a[n,n];

for i := n-1 downto 1 do

begin s := 0.0;

for k := i+1 to N do s := s + a[i,k] * c[k];

c[i] := (c[i] - s) / a[i,i];

end;

еnd.

Формальные параметры процедуры.Входные: m, N(типinteger) -mуравнений иnне­из­вест­ных оп­ре­де­ляют раз­мер матрицыА.Для этой про­це­ду­рыобя­за­тельно вы­пол­не­ние m > n, т.е. ко­ли­чествоурав­­не­ний долж­но быть боль­­­­ше ко­личества не­из­вест­­ных (сис­­те­ма пе­ре­о­пре­делена), пре­дель­­ныйраз­­­ре­шен­ный слу­чай m = n; а - мат­рица, со­став­­­лен­­наяиз ко­эф­фи­ци­ен­тов при не­из­вест­­ных;В - мас­сив, со­став­­ленный из стол­бца сво­бод­ных чле­­нов.Вы­ход­ные:С - массив, в ко­­то­ром со­дер­жит­ся решение сис­темы.

Для проверки правильности работы про­це­ду­ры ре­­ша­лась система 4´4 линейных уравнений ме­то­дом квад­рат­ных корней с точностью 0.0001:

0.68X₁ + 0.05X₂ + 0.11X₃ + 0.08X₄ = 2.15 ;

0.05X₁ + 0.13X₂ + 0.27X₃ + 0.80X₄ = 0.44 ;

0.11X₁ + 0.27X₂ + 0.28X₃ + 0.06X₄ = 0.83 ;

0.08X₁ + 0.80X₂ + 0.06X₃ + 0.12X₄ = 1.16 ,

решение которой, полученное методом квад­рат­ных кор­ней, будет следующим

Х ₁ = 2.97; Х₂ = 1.11; Х₃ = 0.74; Х₄ = -0.07.