Текст подпрограммы и версий
asg6r_c.zip , asg6c_c.zip , asg6d_c.zip 		Тексты тестовых примеров
tasg6r_c.zip , tasg6c_c.zip , tasg6d_c.zip

Подпрограмма:  asg6r_c

Назначение

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

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

Находится решение системы A*X = B, где А - невырожденая квадратная матрица порядка N, В - заданный вектор длины N. Для решения используется нормализованное приведение матрицы системы А к верхней треугольной форме R с помощью последовательности преобразований отражения 

 (1)     QN-1 ... Q2*Q1*A*S = R ,

где Qi - соответствующие матрицы отражения, S - результирующая матрица перестановок, R - верхняя треугольная матрица.

Из полученной треугольной системы 

     R*Y = QN-1 ... Q2*Q1*B

находится решение Y, по которому затем определяется искомое решение Х = S*Y.

Дж.Х.Уилкинсон. Алгебраическая проблема собственных значений. Изд."Наука", М., 1970.

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

    int asg6r_c (real *a, real *b, real *x, real *t, integer *is,
            integer *n, integer *l)

   Параметры
a - вещественный двумерный массив размера n на n, в котором задается исходная матрица; в результате работы подпрограммы в массиве a на соответствующих местах запоминаются наддиагональные элементы вычисленной верхней треугольной матрицы R, в остальной части массива a в последовательных столбцах запоминаются векторы, порождающие соответствующие матрицы отражения;
b - вещественный вектор длины n, в котором задается правая часть системы;
x - вещественный вектор длины n, в котором запоминается полученное решение системы;
t - вещественный вектор длины n, в котором запоминаются диагональные элементы верхней треугольной матрицы R;
is - целый вектор длины n, в котором запоминаются перестановки столбцов при проведении нормализованного процесса (1); в is (k) запоминается номер столбца, переставленного с k - ым столбцом на k - ом шаге преобразования;
n - порядок исходной матрицы (тип: целый);
l - задает режим работы подпрограммы (тип: целый); при этом:
l = 1 - если система с данной матрицей A решается впервые;
l ≠ 1 - если повторно решается система с матрицей A, но другой правой частью; в этом случае факторизация (1) матрицы A не выполняется. 
   Версии
asg6c_c - решение комплексной системы линейных алгебраических уравнений методом отражений.
asg6d_c - решение системы линейных алгебраических уравнений, заданной с удвоенной точностью, методом отражений. 
   Вызываемые подпрограммы:  нет 
   Замечания по использованию
  1. 

В подпрограмме asg6c_c параметры a, b, x, t имеют тип complex.

  2. 

В подпрограмме asg6d_c параметры a, b, t, x имеют тип double.

  3.  При повторном решении системы с той же матрицей A, но с другой правой частью, информация, полученная ранее в массивах a, t, is не должна портиться.

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

int main(void)
{
    /* Initialized data */
    static float b[6] = { 1.f,2.f,3.f,4.f,5.f,1.f };

    /* System generated locals */
    int i__1, i__2;

    /* Local variables */
    extern int asg6r_c(float *, float *, float *, float *, int *,
                       int *, int *);
    static float a[36] /* was [6][6] */;
    static int i__, j, k, l, n;
    static float t[6], x[6];
    static int kh, is[6];

#define a_ref(a_1,a_2) a[(a_2)*6 + a_1 - 7]

    n = 6;
    l = 1;
    k = n + 1;
    i__1 = n;
    for (j = 1; j <= i__1; ++j) {
        i__2 = n;
        for (i__ = j; i__ <= i__2; ++i__) {
            kh = k - i__;
            a_ref(i__, j) = (float) kh;
/* l21: */
            a_ref(j, i__) = (float) kh;
        }
/* l20: */
    }
    for (i__ = 1; i__ <= 6; ++i__) {
         printf("\n  %12.4e %12.4e %12.4e %12.4e %12.4e %12.4e \n",
                a_ref(i__, 1), a_ref(i__, 2), a_ref(i__, 3),
                a_ref(i__, 4), a_ref(i__, 5), a_ref(i__, 6));
    }
    for (i__ = 0; i__ <= 4; i__+= 2) {
         printf("\n  %12.4e %12.4e \n", b[i__], b[i__+1]);
    }
    asg6r_c(a, b, x, t, is, &n, &l);

    for (i__ = 0; i__ <= 34; i__+= 2) {
         printf("\n  %16.7e %16.7e \n", a[i__], a[i__+1]);
    }
    for (i__ = 0; i__ <= 4; i__+= 2) {
         printf("\n  %16.7e %16.7e \n", t[i__], t[i__+1]);
    }
    for (i__ = 0; i__ <= 4; i__+= 2) {
         printf("\n  %16.7e %16.7e \n", b[i__], b[i__+1]);
    }
    for (i__ = 0; i__ <= 4; i__+= 2) {
         printf("\n  %16.7e %16.7e \n", x[i__], x[i__+1]);
    }
    printf("\n  %5i %5i %5i %5i %5i %5i \n",
           is[0], is[1], is[2], is[3], is[4], is[5]);
    return 0;
} /* main */


Результат:

      x  =  (-1., 0., 0., 0., 5., -3.)