|
Текст подпрограммы и версий ash4r_c.zip , ash4d_c.zip |
Тексты тестовых примеров tash4r_c.zip , tash4d_c.zip |
Решение разреженной линейной системы с симметричной положительно определенной матрицей, заданной своим треугольным разложением в формате RR (U) O.
Описание формата RR (U) U приведено в описании подпрограммы am21r_c . Формат RR (U) О совпадает с форматом RR (U) U за тем исключением, что в нем элементы каждой строки матрицы упорядочены, т.е. расположены в порядке возрастания номеров столбцов.
Пусть дана система линейных алгебраических уравнений Аx = b, где А - симметричная положительно определенная разреженная матрица. Пусть, кроме того, матрица А предварительно приведена подпрограммой afh7r_c к виду 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_c ), то вместо делений при определении 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
int ash4r_c (integer *iu, integer *ju, real *un, real *di,
integer *n, real *b, real *x)
Параметры
|
iu, ju - un | заданные портрет и ненулевые элементы верхней треугольной матрицы U с единичной диагональю в формате RR (U) O; |
| di - | вещественный одномерный массив длины n, содержащий элементы матрицы, обратной к диагональной матрице D; |
| n - | заданный порядок системы (тип: целый); |
| b - | вещественный одномерный массив длины n, содержащий компоненты вектора правой части системы; |
| x - | вещественный одномерный массив длины n, содержащий вычисленные компоненты вектора решения системы |
Версии
| ash4d_c - | решение разреженной линейной системы с симметричной положительно определенной матрицей, заданной своим треугольным разложением в формате RR (U) O, в режиме удвоенной точности; при этом параметры un, di, b и x должны иметь тип double. |
Вызываемые подпрограммы: нет
Замечания по использованию: нет
int main(void)
{
/* Initialized data */
static int iu[6] = { 1,2,3,4,5,5 };
static int ju[4] = { 5,5,5,5 };
static float un[4] = { .125f,.8f,.6666667f,2.f };
static float di[5] = { .0625f,1.6f,.3333333f,2.f,60.f };
static float b[5] = { -4.f,-4.f,7.f,3.f,7.f };
/* Local variables */
extern int ash4r_c(int *, int *, float *, float *, int *, float *,
float *);
static int n, i__;
static float x[5];
n = 5;
ash4r_c(iu, ju, un, di, &n, b, x);
for (i__ = 1; i__ <= 5; ++i__) {
printf("\n %15.6e \n", x[i__-1]);
}
return 0;
} /* main */
Результаты: x = (- 0.499996, - 7.99998, 1.00002, 2.00006, 1.99997)