Текст подпрограммы и версий
afh4c_c.zip , afh4p_c.zip
Тексты тестовых примеров
tafh4c_c.zip , tafh4p_c.zip

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

Назначение

Приведение комплексной эрмитовой матрицы к симметрической трехдиагональной матрице унитарными преобразованиями подобия.

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

Подпрограмма afh4c_c реализует алгоритм приведения эрмитовой матрицы А унитарными преобразованиями к симметрической трехдиагональной форме F = V*Q*АQV, где Q - произведение унитарных эрмитовых матриц, а V - унитарная диагональная матрица.

Мartin R.S., Reinsch С. and Wilkinson J.Н., Нousholder's Тridiagonalization of a Symmetric Мatrix, Numeriche Mатнемаtic, 11, 181-195, 1968.

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

    int afh4c_c (integer *n, real *ar, real *ai, real *d,
            real *e, real *rab)

Параметры

n - порядок исходной матрицы (тип: целый);
ar, ai - вещественные двумерные массивы размера n*n, содержащие соответственно действительную и мнимую части исходной комплексной эрмитовой матрицы; в результате работы подпрограммы в строгом нижнем треугольнике массива ai в последовательных строках запоминаются векторы, порождающие матрицу преобразования Q;
d - вещественный вектор длины n, содержащий диагональные элементы вещественной симметрической трехдиагональной матрицы;
e - вещественный вектор длины n, содержащий в последних n - 1 компонентах поддиагональные элементы симметрической трехдиагональной матрицы; первый элемент e приравнивается к нулю;
rab - вещественный вектор длины 2*n, содержащий диагональную матрицу преобразования V.

Версии

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

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

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

 

Исходные матрицы ar, ai можно задавать нижними треугольниками.

Подпрограмма afh4c_c сохраняет полный верхний треугольник массива ar и строгий верхний треугольник массива ai, в остальных компонентах этих массивов помещается информация об унитарных преобразованиях.

В подпрограмме afh4p_c массивы ar, ai, e, d и rab имееют тип double.

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

int main(void)
{
    /* Initialized data */
    static float ar[16] /* was [4][4] */ = { 1.f,.42f,.54f,.66f,.42f,1.f,.32f,
                                 .44f,.54f,.32f,1.f,.22f,.66f,.44f,.22f,1.f };
    static float ai[16] /* was [4][4] */ = { 0.f,0.f,0.f,0.f,0.f,0.f,0.f,0.f,
                                            0.f,0.f,0.f,0.f,0.f,0.f,0.f,0.f };
    /* Local variables */
    extern int afh4c_c(int *, float *, float *, float *, float *, float *);
    static float d__[4], e[4];
    static int i__, n;
    static float rab[8];

#define ai_ref(a_1,a_2) ai[(a_2)*4 + a_1 - 5]
#define ar_ref(a_1,a_2) ar[(a_2)*4 + a_1 - 5]

    for (i__ = 1; i__ <= 4; ++i__) {
        printf("\n  %14.7e %14.7e %14.7e %14.7e \n",
             ar_ref(i__, 1), ar_ref(i__, 2), ar_ref(i__, 3), ar_ref(i__, 4));
    }
    n = 4;
    afh4c_c(&n, ar, ai, d__, e, rab);

    for (i__ = 1; i__ <= 4; ++i__) {
        printf("\n  %14.7e %14.7e %14.7e %14.7e \n",
             ar_ref(i__, 1), ar_ref(i__, 2), ar_ref(i__, 3), ar_ref(i__, 4));
    }
    for (i__ = 1; i__ <= 4; ++i__) {
        printf("\n  %14.7e %14.7e %14.7e %14.7e \n",
             ai_ref(i__, 1), ai_ref(i__, 2), ai_ref(i__, 3), ai_ref(i__, 4));
    }
    printf("\n  %14.7e %14.7e %14.7e %14.7e \n",
             d__[0], d__[1], d__[2], d__[3]);
    printf("\n  %14.7e %14.7e %14.7e %14.7e \n",
             e[0], e[1], e[2], e[3]);
    return 0;
} /* main */


Результаты:

                      | 1.0000   0.4200   0.5400   0.6600 |
      ar_ref  =   | 0.1758   1.0000   0.3200   0.4400 |
                      | 0.4454   0.5819   1.0000   0.2200 |
                      | 0.6600   0.4400   1.0432   1.0000 |

                      | 0.0   0.0         0.0             0.0       |
      ai_ref  =   | 0.0   0.1243   0.0             0.0       |
                      | 0.0   0.0         0.5181       0.0       |
                      | 0.0   0.0         0.0             0.9267 |
  
                     | 0.667054 |
      d__  =     | 0.650088 |
                     | 1.682857 |
                     | 1.000000 |

                   | 0.000000 |
      e    =     | 0.087921 |
                   | 0.461386 |
                   | 0.823365 |