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

Назначение

Умножение матрицы, обратной к заданной в компактной форме вещественной ленточной матрице A, на прямоугольную матрицу и оценка числа обусловленности матрицы A.

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

Для заданной в компактной форме ленточной вещественной

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

      RCOND  =  1 / ( || A ||1 - || A-1 ||1 )  ,   где
                                           N
        || A ||1 =      max     {   ∑   | ai j |   }
                         j=1,...N       i=1 

и затем для заданной вещественной прямоугольной матрицы B размера N на NN вычисляется прямоугольная матрица С = А^- 1 В размера N на NN путем решения NN систем линейных алгебраических уравнений А*С (J) = В(J), где С(J) и В(J) суть J - е столбцы матриц C и B,  J = 1,...,NN.

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

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

    int ammhr_c (real *a, integer *ma, integer *n, integer *ml,
        integer *mu, integer *nlead, real *b, integer *mb, integer *nn,
        real *rcond, real *z, integer *ierr)

   Параметры 

   Версии 

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

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

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

int main(void)
{
    /* System generated locals */
    int i__1, i__2, i__3;

    /* Local variables */
    static int ierr;
    static float a[25] /* was [5][5] */,
                 b[50] /* was [5][10] */;
    static int i__, j, k, n, nlead[5];
    static float z__[5];
    extern int ammhr_c(float *, int *, int *, int *, int *, int *,
                       float *, int *, int *, float *, float *, int *);
    static float rcond;
    static int j0, j1, ma, mb, ml, nn, mu;

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

    ma = 5;
    n = 5;
    ml = 1;
    mu = 1;
    mb = 5;
    nn = 10;
    i__1 = ma;
    for (i__ = 1; i__ <= i__1; ++i__) {
        i__2 = n;
        for (j = 1; j <= i__2; ++j) {
            a_ref(i__, j) = 0.f;
            b_ref(i__, j) = 0.f;
            b_ref(i__, j + 5) = 0.f;
/* L1: */
        }
/* L2: */
    }
    i__1 = ma;
    for (i__ = 1; i__ <= i__1; ++i__) {
/* Computing max */
        i__2 = 1, i__3 = i__ - ml;
        j0 = max(i__2,i__3);
/* Computing min */
        i__2 = n, i__3 = i__ + mu;
        j1 = min(i__2,i__3);
        i__2 = j1;
        for (j = j0; j <= i__2; ++j) {
            k = j - i__ + ml + 1;
            a_ref(i__, k) = (float) (i__ * 10 + j);
            b_ref(i__, j) = a_ref(i__, k);
            b_ref(i__, j + 5) = a_ref(i__, k);
/* L3: */
        }
/* L4: */
    }
    for (i__ = 1; i__ <= 5; ++i__) {
         printf("\n  %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));
    }
    for (i__ = 1; i__ <= 5; ++i__) {
         printf("\n  %12.4e %12.4e %12.4e %12.4e %12.4e \n",
                b_ref(i__, 1), b_ref(i__, 2), b_ref(i__, 3),
                b_ref(i__, 4), b_ref(i__, 5));
         printf("\n  %12.4e %12.4e %12.4e %12.4e %12.4e \n",
                b_ref(i__, 6), b_ref(i__, 7), b_ref(i__, 8),
                b_ref(i__, 9), b_ref(i__, 10));
    }
    ammhr_c(a, &ma, &n, &ml, &mu, nlead, b, &mb, &nn, &rcond, z__, &ierr);

    printf("\n  %5i %5i %5i %5i %5i \n",
           nlead[0], nlead[1], nlead[2], nlead[3], nlead[4]);
    for (i__ = 1; i__ <= 5; ++i__) {
         printf("\n  %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));
    }
    printf("\n  %5i \n", ierr);
    printf("\n  %16.7e \n", rcond);
    for (i__ = 1; i__ <= 5; ++i__) {
         printf("\n  %12.4e %12.4e %12.4e %12.4e %12.4e \n",
                b_ref(i__, 1), b_ref(i__, 2), b_ref(i__, 3),
                b_ref(i__, 4), b_ref(i__, 5));
         printf("\n  %12.4e %12.4e %12.4e %12.4e %12.4e \n",
                b_ref(i__, 6), b_ref(i__, 7), b_ref(i__, 8),
                b_ref(i__, 9), b_ref(i__, 10));
    }
    return 0;
} /* main */


Результаты:

                 |   0         21.0     22.0   23.0    0 |
                 | -0.524   32.0     33.0   34.0    0 |
   a_ref  =  | -0.015   43.0     44.0   45.0    0 |
                 |  0.292   54.0     55.0     0       0 |
                 | -0.228    0.569    0        0       0 |

   nlead  =  (2, 3, 4, 5, 5) ,   rcond  =  1.47362E-03 ,

                 | 1.0  0.0  0.0  0.0  0.0  1.0  0.0  0.0  0.0  0.0 |
                 | 0.0  1.0  0.0  0.0  0.0  0.0  1.0  0.0  0.0  0.0 |
   b_ref  =  | 0.0  0.0  1.0  0.0  0.0  0.0  0.0  1.0  0.0  0.0 |
                 | 0.0  0.0  0.0  1.0  0.0  0.0  0.0  0.0  1.0  0.0 |
                 | 1.0  0.0  0.0  0.0  1.0  0.0  0.0  0.0  0.0  1.0 |