Текст подпрограммы и версий
ask1r_p.zip , ask1e_p.zip , ask1c_p.zip 		Тексты тестовых примеров
task1r_p.zip , task1e_p.zip , task1c_p.zip

Подпрограмма:  ASK1R (модуль ASK1R_p)

Назначение

Решение вещественной системы линейных алгебраических уравненений AX = b с клеточно - теплицевой матрицей A специального вида.

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

Матрица А состоит из Р на Р клеток размера М на М. В силу клеточно - теплицевости А полностью определяется своими первой клеточной строкой и первым клеточным столбцом.

Обозначим SK, RK, 1 ≤ К ≤ Р соответственно К - ю клетку первой клеточной строки и К - ю клетку первого клеточного столбца матрицы А. В данном случае под специальным видом матрицы А понимается выполнение соотношений 

     RK = L SK N ,   1 ≤ K ≤ P ,

где L и N - некоторые фиксированные симметричные матрицы перестановок.

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

С.Н.Воеводина, Решение системы уравнений с клеточно-теплицевыми матрицами. Сб. "Вычислительные методы и программирование", вып.24, Изд-во МГУ, 1975.

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

procedure ASK1R(var A1 :Array of Real; var B :Array of Real;
                var Y :Array of Real; var B3 :Array of Real;
                var L :Array of Integer; var N :Array of Integer;
                var C :Array of Real; var Q :Array of Integer;
                M :Integer; P :Integer);


   Параметры
A1 - трехмерный М на М на Р массив, в котором задается первая клеточная строка матрицы системы (тип: вещественный);
B - двумерный М на Р массив, в котором задается вектор правой части исходной системы (тип: вещественный);
Y - двумерный М на Р массив, в котором запоминается вектор вычисленного решения (тип: вещественный);
B3 - трехмерный М на М на Р массив, используемый подпрограммой как рабочий (тип: вещественный);
L, N - векторы длины М, в которых задаются две матрицы перестановок (тип: целый);
C - трехмерный М на М на 3 массив, используемый подпрограммой как рабочий (тип: вещественный);
Q - вектор длины М, используемый как рабочий (тип: целый);
M - заданный порядок клетки (тип: целый);
P - заданное число клеток в клеточной строке или столбце (тип: целый). 
   Версии
ASK1E - решение системы линейных алгебраических уравнений А*x = b с клеточно - теплицевой матрицей специального вида для вещественных А, b, заданных с расширенной (Extended) точностью.
ASK1C - решение системы линейных алгебраических уравнений А*x = b с клеточно - теплицевой матрицей специального вида для комплексных А, b. 
   Вызываемые подпрограммы:  нет 
   Замечания по использованию
  1. 

В массиве В3 по окончании работы подпрограммы запоминается последний клеточный столбец матрицы, являющейся обратной для матрицы исходной системы.

  2. 

Вектор L определяет симметричную матрицу перестановок следующим образом: K - я компонента вектора L есть номер строки, которая переставляется с К - ой. Аналогично, К - я компонента вектора N есть номер столбца, переставляемого с К - м столбцом.

  3. 

В подпрограмме АSК1E массивы А1, В, Y, В3, С имеют тип Extended.

  4.  В подпрограмме АSК1С массивы А1, В, Y, В3, С имеют тип Complex.

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

1.
Unit TASK1R_p;
interface
uses
SysUtils, Math, { Delphi }
Lstruct, Lfunc, UtRes_p, ASK1R_p;

function TASK1R: String;

implementation

function TASK1R: String;
var
M,P,I,J,K,_i :Integer;
Y :Array [0..8] of Real;
B3 :Array [0..26] of Real;
C :Array [0..26] of Real;
Q :Array [0..2] of Integer;
const
A1 :Array [0..26] of Real = ( 1.0,1.0,0.0,0.0,1.0,1.0,1.0,0.0,1.0,0.0,0.0,0.0,
0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,
1.0,0.0,0.0 );
B :Array [0..8] of Real = ( 3.0,2.0,2.0,2.0,2.0,2.0,3.0,2.0,2.0 );
L :Array [0..2] of Integer = ( 1,3,2 );
N :Array [0..2] of Integer = ( 3,2,1 );
begin
Result := '';  { результат функции }
M := 3;
P := 3;
ASK1R(A1,B,Y,B3,L,N,C,Q,M,P);
Result := Result + #$0D#$0A;
for _i:=0 to 8 do
 begin
  Result := Result + Format('%20.16f ',[Y[_i]]);
  if ( ((_i+1) mod 4)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + #$0D#$0A;
for _i:=0 to 8 do
 begin
  Result := Result + Format('%20.16f ',[B[_i]]);
  if ( ((_i+1) mod 4)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + #$0D#$0A;
for _i:=0 to 26 do
 begin
  Result := Result + Format('%20.16f ',[A1[_i]]);
  if ( ((_i+1) mod 4)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + #$0D#$0A;
for _i:=0 to 26 do
 begin
  Result := Result + Format('%20.16f ',[B3[_i]]);
  if ( ((_i+1) mod 4)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
UtRes('TASK1R',Result);  { вывод результатов в файл TASK1R.res }
exit;
end;

end.


Результат:

       Y  =   (1., 1., 1., 1., 1., 1., 1., 1., 1.)

2. 
Unit task1e_p;
interface
uses
SysUtils, Math, { Delphi }
Lstruct, Lfunc, UtRes_p, ASK1E_p;

function task1e: String;

implementation

function task1e: String;
var
M,P,I,J,K,_i :Integer;
Y :Array [0..5] of Extended;
B3 :Array [0..11] of Extended;
C :Array [0..11] of Extended;
Q :Array [0..1] of Integer;
const
A1 :Array [0..11] of Extended = ( 2.0e0,1.0e0,3.0e0,2.0e0,4.0e0,2.0e0,6.0e0,4.0e0,
-2.0e0,-1.0e0,-3.0e0,-2.0e0 );
B :Array [0..5] of Extended = ( 10.0e0,6.0e0,25.0e0,15.0e0,10.0e0,6.0e0 );
L :Array [0..1] of Integer = ( 1,2 );
N :Array [0..1] of Integer = ( 1,2 );
begin
Result := '';  { результат функции }
M := 2;
P := 3;
ASK1E(A1,B,Y,B3,L,N,C,Q,M,P);
Result := Result + #$0D#$0A;
for _i:=0 to 5 do
 begin
  Result := Result + Format('%20.16f ',[Y[_i]]);
  if ( ((_i+1) mod 4)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + #$0D#$0A;
for _i:=0 to 5 do
 begin
  Result := Result + Format('%20.16f ',[B[_i]]);
  if ( ((_i+1) mod 4)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + #$0D#$0A;
for _i:=0 to 11 do
 begin
  Result := Result + Format('%20.16f ',[A1[_i]]);
  if ( ((_i+1) mod 4)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + #$0D#$0A;
for _i:=0 to 11 do
 begin
  Result := Result + Format('%20.16f ',[B3[_i]]);
  if ( ((_i+1) mod 4)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
UtRes('task1e',Result);  { вывод результатов в файл task1e.res }
exit;
end;

end.

Результат:

       Y  =   (1., 1., 1., 1., 1., 1.)

3.
Unit task1c_p;
interface
uses
SysUtils, Math, { Delphi }
Lstruct, Lfunc, UtRes_p, ASK1C_p;

function task1c: String;

implementation

function task1c: String;
var
M,P,I,J,K,_i :Integer;
Y :Array [0..5] of Complex;
B3 :Array [0..11] of Complex;
C :Array [0..11] of Complex;
Q :Array [0..1] of Integer;
const
A1 :Array [0..11] of Complex = ( ( re:2.0; im:0.0 ),( re:1.0; im:0.0 ),( re:3.0; 
im:0.0 ),( re:2.0; im:0.0 ),( re:4.0; im:0.0 ),( 
re:2.0; im:0.0 ),( re:6.0; im:0.0 ),( re:4.0; 
im:0.0 ),( re:-2.0; im:0.0 ),( re:-1.0; im:0.0 ),
( re:-3.0; im:0.0 ),( re:-2.0; im: 0.0 ) );
B :Array [0..5] of Complex = ( ( re:10.0; im:0.0 ),( re:6.0; im:0.0 ),( 
re:25.0; im:0.0 ),( re:15.0; im:0.0 ),( re:10.0; 
im:0.0 ),( re:6.0; im: 0.0 ) );
L :Array [0..1] of Integer = ( 1,2 );
N :Array [0..1] of Integer = ( 1,2 );
begin
Result := '';  { результат функции }
P := 3;
M := 2;
ASK1C(A1,B,Y,B3,L,N,C,Q,M,P);
Result := Result + #$0D#$0A;
for _i:=0 to 5 do
 begin
  Result := Result + Format('%20.16f %20.16f ',[Y[_i].re,Y[_i].im]);
  if ( ((_i+1) mod 2)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + #$0D#$0A;
for _i:=0 to 5 do
 begin
  Result := Result + Format('%20.16f %20.16f ',[B[_i].re,B[_i].im]);
  if ( ((_i+1) mod 2)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + #$0D#$0A;
for _i:=0 to 11 do
 begin
  Result := Result + Format('%20.16f %20.16f ',[A1[_i].re,A1[_i].im]);
  if ( ((_i+1) mod 2)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + #$0D#$0A;
for _i:=0 to 11 do
 begin
  Result := Result + Format('%20.16f %20.16f ',[B3[_i].re,B3[_i].im]);
  if ( ((_i+1) mod 2)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
UtRes('task1c',Result);  { вывод результатов в файл task1c.res }
exit;
end;

end.

Результат:

       Y  =  ( (1., 0.),  (1., 0.),  (1., 0.),  (1., 0.),  (1., 0.),  (1., 0.) )