Текст подпрограммы и версий
ash4r_p.zip , ash4e_p.zip 		Тексты тестовых примеров
tash4r_p.zip , tash4e_p.zip

Подпрограмма:  ASH4R (модуль ASH4R_p)

Назначение

Решение разреженной линейной системы с симметричной положительно определенной матрицей, заданной своим треугольным разложением в формате RR (U) O.

Математическое описание

Описание формата RR (U) U приведено в описании подпрограммы AM21R . Формат RR (U) О совпадает с форматом RR (U) U за тем исключением, что в нем элементы каждой строки матрицы упорядочены, т.е. расположены в порядке возрастания номеров столбцов.

Пусть дана система линейных алгебраических уравнений АХ = b, где А - симметричная положительно определенная разреженная матрица. Пусть, кроме того, матрица А предварительно приведена подпрограммой AFH7R к виду A = UTDU, где U - верхняя треугольная матрица с единичными диагональными элементами, представленная в формате RR (U) О, а D - диагональная матрица. Иными словами, предполагается, что известно упорядоченное строчное треугольное разложение матрицы А системы в формате RR (U) О. Требуется найти решение системы Ах = UTDUх = b, т.е. выполнить прямой ход и обратную подстановку.

Полагая  z = DUx,  w = Ux, получаем следующие три линейные системы: 

               UTz  =  b
               Dw  =  z
               Ux  =  w

Решения этих систем находим последовательно: вначале  z, затем  w и, наконец,  x. Поскольку матрицы U и UT треугольные, а матрица D диагональная, то векторы  z,  w и  x вычисляются при помощи несложного алгоритма. Рассмотрим этот алгоритм на примере следующей системы: 

   |  1   0   0  |     |  a   0   0  |     |  1    d   e  |     |  x1  |          |  b1  |
   |  d   1   0  |     |  0   b   0  |     |  0    1   f  |     |  x2  |    =    |  b2  | 
   |  e   f    1  |     |  0   0   c  |     |  0    0   1  |     |  x3  |          |  b3  |

В этом примере система UTz = b имеет вид: 

   |  1   0   0  |     |  z1  |          |  b1  |
   |  d   1   0  |     |  z2  |    =    |  b2  | 
   |  e   f    1  |     |  z3  |          |  b3  | 

Последовательной подстановкой получаем

               z1  =  b1
               z2  =  b2 - d z1
               z3  =  b3 - e z1 - f z2

Вторая система  Dw = z имеет вид: 

   |  a   0   0  |     |  w1  |          |  z1  |
   |  0   b   0  |     |  w2  |    =    |  z2  | 
   |  0   0   c  |     |  w3  |          |  z3  | 

Отсюда имеем:

               w1  =  z1 / a
               w2  =  z2 / b
               w3  =  z3 / c

Так как в памяти машины треугольное разложение хранится таким образом, что вместо D известна сразу же обратная матрица D- 1 (см. подпрограмму AFH7R ), то вместо делений при определении  wi ( i = 1, 2, 3) следует использовать более быстрые операции умножения. Наконец, решение последней системы 

   |  1   d   e  |     |  x1  |          |  w1  |
   |  0   1   f   |     |  x2  |    =    |  w2  | 
   |  0   0   1  |     |  x3  |          |  w3  | 

вычисляется обратной подстановкой 

               x3  =  w3
               x2  =  w2 - f x3
               x1  =  w1 - d x2 - e x3

С.Писсанецки. Технология разреженных матриц. - М.: Мир, 1998

Использование 

procedure ASH4R(var IU :Array of Integer; var JU :Array of Integer;
                var UN :Array of Real; var DI :Array of Real;
                N :Integer; var B :Array of Real;
                var X :Array of Real);


   Параметры
     IU, JU -
      UN   		заданные портрет и ненулевые элементы верхней треугольной матрицы U с единичной диагональю в формате RR (U) O;
DI - вещественный одномерный массив длины N, содержащий элементы матрицы, обратной к диагональной матрице D;
N - заданный порядок системы (тип: целый);
B - вещественный одномерный массив длины N, содержащий компоненты вектора правой части системы;
X - вещественный одномерный массив длины N, содержащий вычисленные компоненты вектора решения системы 
   Версии
ASH4E - решение разреженной линейной системы с симметричной положительно определенной матрицей, заданной своим треугольным разложением в формате RR (U) O, в режиме расширенной (Extended) точности; при этом параметры UN, DI, B и X должны иметь тип Extended. 
   Вызываемые подпрограммы:  нет 
   Замечания по использованию:  нет

Пример использования 

Unit TASH4R_p;
interface
uses
SysUtils, Math, { Delphi }
Lstruct, Lfunc, UtRes_p, ASH4R_p;

function TASH4R: String;

implementation

function TASH4R: String;
var
N,_i :Integer;
X :Array [0..4] of Real;
const
IU :Array [0..5] of Integer = ( 1,2,3,4,5,5 );
JU :Array [0..3] of Integer = ( 5,5,5,5 );
UN :Array [0..3] of Real = ( 0.125,0.8,0.6666667,2.0 );
DI :Array [0..4] of Real = ( 0.0625,1.6,0.3333333,2.0,60.0 );
B :Array [0..4] of Real = ( -4.0,-4.0,7.0,3.0,7.0 );
begin
Result := '';  { результат функции }

{      ТЕСТ ДЛЯ ПОДПРОГРАММЫ ASH4R }

N := 5;
ASH4R(IU,JU,UN,DI,N,B,X);
Result := Result + Format('%s',[' X=']);
Result := Result + #$0D#$0A;
for _i:=0 to 4 do
 begin
  Result := Result + Format('%20.16f ',[X[_i]]);
  if ( ((_i+1) mod 4)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
UtRes('TASH4R',Result);  { вывод результатов в файл TASH4R.res }
exit;
end;

end.


Результаты:  X = (- 0.499996, - 7.99998, 1.00002, 2.00006, 1.99997)