Текст подпрограммы и версий
aig4r_c.zip  aig4d_c.zip  aig4c_c.zip 
Тексты тестовых примеров
taig4r_c.zip  taig4d_c.zip  taig4c_c.zip 

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

Назначение

Обращение вещественной матрицы общего вида методом Гаусса с выбором ведущего элемента по столбцу.

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

Для заданной вещественной kвaдpaтной матрицы А порядка N выполняется треугольная факторизация  L- 1А = U, где U - верхняя треугольная матрица. Далее вычисляется матрица  U- 1,  и затем вычисляется  А- 1 = U- 1L-1.

Дж. Форсайт, М. Малькольм, К. Моулер. Машинные методы математических вычислений. М.: Мир, 1980.

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

    int aig4r_c (real *a, integer *m, integer *n, integer *nlead,
            real *z__, integer *ierr)

Параметры

a - вещественный двумерный массив размера m на n, в котором задается исходная матрица; на выходе содержит вычисленную обращенную матрицу;
m - первая размерность массива a в вызывающей программе (тип: целый);
n - порядок обращаемой матрицы (тип: целый);
nlead - целый вектор длины n, содержащий на выходе информацию о выполненных в ходе факторизации перестановках (см. замечания по использованию);
z - вещественный рабочий вектор длины n;
ierr - целая переменная, содержащая на выходе информацию о прохождении счета; при этом:
ierr=65 - если m ≤ 0 или n ≤ 0;
ierr=66 - если в процессе работы подпрограммы произошло переполнение (это говорит о том, что некоторые элементы матрицы U или матрицы  A- 1  превосходят по абсолютной величине максимальное представимое на данной машине число);
ierr=-k - если матрица A вырождена (k - номер последней из строк матрицы U, содержащих нулевой диагональный элемент).

Версии

aig4d_c - обращение методом Гаусса с выбором ведущего элемента по столбцу вещественной матрицы A, заданной с удвоенной точностью.
aig4c_c - обращение методом Гаусса с выбором ведущего элемента по столбцу комплексной матрицы A.

Вызываемые подпрограммы

afg5r_c - подпрограмма треугольной факторизации матрицы A.
utafsi_c - подпрограмма выдачи диагностических сообщений.

Замечания по использованию

  1. 

В подпрограмме aig4d_c массивы a и z имеют тип double, для факторизации матрицы A вызывается подпрограмма afg5d_c.

  2. 

В подпрограмме aig4c_c массивы a и z имеют тип complex, для факторизации матрицы A вызывается подпрограмма afg5c_c.

  3. 

На выходе k - й элемент вектора nlead равен номеру строки, переставленной на k - м шаге факторизации с k - й строкой матрицы A. Поскольку факторизация Гаусса требует n - 1 шагов, то nlead (n) = n.

  4.  Если вырабатывается значение ierr ≠ 0, то выдается соответствующее диагностическое сообщение и происходит выход из подпрограммы.

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

int main(void)
{
    /* Initialized data */
    static float a[16] /* was [4][4] */ = { 7.9f,8.5f,4.3f,3.2f,5.6f,-4.8f,
                       4.2f,-1.4f,5.7f,.8f,-3.2f,-8.9f,-7.2f,3.5f,9.3f,3.3f };
    /* System generated locals */
    int i__1, i__2, i__3;

    /* Local variables */
    static int ierr;
    extern int aig4r_c(float *, int *, int *, int *, float *, int *);
    static float b[16] /* was [4][4] */,
               c__[16] /* was [4][4] */;
    static int i__, j, k, m, n, nlead[4];
    static float z__[4];

#define a_ref(a_1,a_2) a[(a_2)*4 + a_1 - 5]
#define b_ref(a_1,a_2) b[(a_2)*4 + a_1 - 5]
#define c___ref(a_1,a_2) c__[(a_2)*4 + a_1 - 5]

    m = 4;
    n = m;
    i__1 = m;
    for (i__ = 1; i__ <= i__1; ++i__) {
        i__2 = n;
        for (j = 1; j <= i__2; ++j) {
            b_ref(i__, j) = a_ref(i__, j);
/* l1: */
        }
/* l2: */
    }
    for (i__ = 1; i__ <= 4; ++i__) {
    printf("\n %16.7e %16.7e %16.7e %16.7e \n",
           a_ref(i__,1), a_ref(i__,2), a_ref(i__,3), a_ref(i__,4));
    }
    aig4r_c(a, &m, &n, nlead, z__, &ierr);

    printf("\n %5i %5i %5i %5i \n",
           nlead[0], nlead[1], nlead[2], nlead[3]);
    for (i__ = 1; i__ <= 4; ++i__) {
    printf("\n %16.7e %16.7e %16.7e %16.7e \n",
           a_ref(i__,1), a_ref(i__,2), a_ref(i__,3), a_ref(i__,4));
    }
    printf("\n %5i \n", ierr);
    i__1 = m;
    for (i__ = 1; i__ <= i__1; ++i__) {
        i__2 = n;
        for (j = 1; j <= i__2; ++j) {
            c___ref(i__, j) = 0.f;
            i__3 = n;
            for (k = 1; k <= i__3; ++k) {
                c___ref(i__, j) = c___ref(i__, j) +
                        b_ref(i__, k) * a_ref(k, j);
/* l3: */
            }
/* l4: */
        }
/* l5: */
    }
    for (i__ = 1; i__ <= 4; ++i__) {
    printf("\n %16.7e %16.7e %16.7e %16.7e \n",
       c___ref(i__,1), c___ref(i__,2), c___ref(i__,3), c___ref(i__,4));
    }
    i__2 = m;
    for (i__ = 1; i__ <= i__2; ++i__) {
        i__1 = n;
        for (j = 1; j <= i__1; ++j) {
            c___ref(i__, j) = 0.f;
            i__3 = n;
            for (k = 1; k <= i__3; ++k) {
                c___ref(i__, j) = c___ref(i__, j) +
                        a_ref(i__, k) * b_ref(k, j);
/* l6: */
            }
/* l7: */
        }
/* l8: */
    }
    for (i__ = 1; i__ <= 4; ++i__) {
    printf("\n %16.7e %16.7e %16.7e %16.7e \n",
       c___ref(i__,1), c___ref(i__,2), c___ref(i__,3), c___ref(i__,4));
    }
    return 0;
} /* main */


Результаты:    nlead = (2, 2, 4, 4)

               |  0.05056    0.05429    0.00629     0.03500 |
      a  =  |  0.05189    0.08460    0.07212    -0.00030 |
               | -0.00841    0.04319    0.02021    -0.12113 |
               | -0.04971    0.02797    0.07900     0.05773 |