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

Назначение

Вычисление всех собственных значений вещественной матрицы общего вида

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

На первом этапе подпрограмма приводит матрицу общего вида A к верхней форме Хессенберга:
A =  Q*H*QT, где Q - ортогональная матрица, а H - матрица в виде верхней формы Хессенберга.
На втором этапе матрица H приводится к форме Шура:
H =  S*T*ST, где S - матрица, состоящая из векторов Шура матрицы H, а T - матрица в виде формы Шура. Собственные значения матрицы A являются элементами диагонали матрицы T.

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

         CALL  PDGEEV1 ( N, A, IA, JA, DESCA, WR, WI, WORK, LWORK, INFO)

   Параметры
N - порядок распределенной подматрицы sub (A); N ≥ 0 (глобальный входной параметр, тип целый);
A - массив двойной точности, распределенный по процессам блочно - циклическим образом (см.), глобальная размерность которого (N, N), а локальная размерность ( LLD_A, LOCc ( JA + N - 1)),
на входе - это матрица A общего вида;
на выходе - содержание матрицы A не сохраняется (не определено) (локальный входной параметр, рабочее пространство);
IA - глобальный номер строки матрицы A, который указывает на начало подматрицы, которая должна быть обработана; в настоящей реализации всегда IA = 1 (глобальный входной параметр, тип целый);
JA - глобальный номер столбца матрицы A, который указывает на начало подматрицы, которая должна быть обработана; в настоящей реализации всегда JA = 1 (глобальный входной параметр, тип целый);
DESCA - дескриптор распределенной матрицы A (одномерный массив длины DLEN);
если тип дескриптора двумерный    ( DTYPE_A = 1),      DLEN ≥ 9;
DESCA содержит информацию о размещении A в памяти; полное описание DESCA см. в разделе документации "Дескрипторы глобальных массивов" (глобальный и локальный входной параметр, тип целый);
WR - одномерный массив двойной точности длины N; если  INFO = 0, содержит вещественные части собственных значений матрицы A (глобальный выходной параметр);
WI - одномерный массив двойной точности длины N; если  INFO = 0, содержит мнимые части собственных значений матрицы A (глобальный выходной параметр);
WORK - одномерный рабочий массив двойной точности длины LWORK;
на выходе, в элементе WORK (1) возвращается необходимая длина рабочего пространства
(локальное рабочее пространство, локальный выходной параметр);
LWORK - задаваемая длина рабочего пространства WORK;
Минимальная требуемая величина LWORK может быть вычислена самой подпрограммой, если к ней обратиться со значением LWORK = -1; при этом подпрограмма не производит никаких других вычислений, а требуемое значение LWORK возвращается в элементе массива WORK (1);
(локальный и глобальный входной параметр, тип целый);
INFO - целая переменная, диагностирующая результат работы подпрограммы (глобальный выходной параметр)
= 0 - успешное завершение работы;
< 0 - если i - ый фактический параметр подпрограммы является массивом и его j - ый элемент имеет недопустимое значение, тогда INFO = - ( i * 100 + j ), если i - ый фактический параметр является скаляром и имеет недопустимое значение, тогда INFO = - i;
> 0 - то не все собственные значения смогли быть вычислены после выполнения 30 * N итераций;
если INFO = i, где 1 ≤ i ≤ N, то элементы векторов WR, WI с индексами от i+1 до N содержат успешно вычисленные собственные значения; 
   Версии
PZGEEV1 -   вычисление всех собственных значений матрицы общего вида для случая комплексных данных двойной точности,
PSGEEV1 -   вычисление всех собственных значений матрицы общего вида для случаев вещественных данных одинарной точности, 
   Вызываемые подпрограммы

Здесь указаны только базовые подпрограммы (2 - ого уровня), которые вызываются из целевой подпрограммы (1 - ого уровня).

PDGEHRD - PZGEHRD   приведение матрицы общего вида к верхней форме Хессенберга для вещественных и комплексных данных двойной точности
PDLAHQR - PZLAHQR   вычисление всех собственных значений матрицы общего вида, приведенной к форме Хессенберга (PDGEHRD,PZGEHRD) для вещественных и комплексных данных двойной точности
PDLAMCH - вычисление машинных параметров для арифметики с плавающей запятой 
   Замечания по использованию
  1.  В настоящей реализации подпрограмм PDGEEV1 и PZGEEV1 нельзя задавать размер блока (NB), на которые разбивается матрица A при распределении ее частей по паралельным процессам, меньше 7. NB является частью дескриптора DESCA. См. раздел документации "Дескрипторы глобальных массивов"
  2.  В подпрограмме PZGEEV1 параметры A, W и WORK имеют тип COMPLEX*16
Первый оператор подпрограммы PZGEEV1 имеет вид:
SUBROUTINE PZGEEV1( N, A, IA, JA, DESCA, W, WORK, LWORK, INFO ),
где комплексный параметр W содержит вычисленные собственные значения
  3.  В подпрограмме PSGEEV1 параметры A, WR, WI и WORK имеют тип REAL
  4.  Используются подпрограммы BLACS_GRIDINFO ( из пакета BLACS),
DCOPY, DSCAL ( из пакета BLAS),
LSAME, INDXG2P, NUMROC, CHK1MAT, PDLANSY, PXERBLA ( из библиотеки ScaLAPACK_TOOLS)

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

Вычисление всех собственных значений матрицы A, которая имеет вид: 

    |   1   0   0   0   0   0   0   0   0   0    0    0    0    0    0    0    0    0    0    1   |
    |   0   1   0   0   0   0   0   0   0   0    0    0    0    0    0    0    0    0    0    2   |
    |   0   0   1   0   0   0   0   0   0   0    0    0    0    0    0    0    0    0    0    3   |
    |   0   0   0   1   0   0   0   0   0   0    0    0    0    0    0    0    0    0    0    4   |
    |   0   0   0   0   1   0   0   0   0   0    0    0    0    0    0    0    0    0    0    5   |
    |   0   0   0   0   0   1   0   0   0   0    0    0    0    0    0    0    0    0    0    6   |
    |   0   0   0   0   0   0   1   0   0   0    0    0    0    0    0    0    0    0    0    7   |
    |   0   0   0   0   0   0   0   1   0   0    0    0    0    0    0    0    0    0    0    8   |
    |   0   0   0   0   0   0   0   0   1   0    0    0    0    0    0    0    0    0    0    9   |
    |   0   0   0   0   0   0   0   0   0   1    0    0    0    0    0    0    0    0    0  10   |
    |   0   0   0   0   0   0   0   0   0   0    1    0    0    0    0    0    0    0    0  11   |
    |   0   0   0   0   0   0   0   0   0   0    0    1    0    0    0    0    0    0    0  12   |
    |   0   0   0   0   0   0   0   0   0   0    0    0    1    0    0    0    0    0    0  13   |
    |   0   0   0   0   0   0   0   0   0   0    0    0    0    1    0    0    0    0    0  14   |
    |   0   0   0   0   0   0   0   0   0   0    0    0    0    0    1    0    0    0    0  15   |
    |   0   0   0   0   0   0   0   0   0   0    0    0    0    0    0    1    0    0    0  16   |
    |   0   0   0   0   0   0   0   0   0   0    0    0    0    0    0    0    1    0    0  17   |
    |   0   0   0   0   0   0   0   0   0   0    0    0    0    0    0    0    0    1    0  18   |
    |   0   0   0   0   0   0   0   0   0   0    0    0    0    0    0    0    0    0    1  19   | 
    |   1   2   3   4   5   6   7   8   9  10  11  12  13  14  15  16  17  18  19  20  |

Порядок матрицы N = 20.
Пусть матрица разбивается на блоки с размером блоков NB = 7.

Осуществляется пропуск программы на 4 процессах, которые образуют решетку 2 на 2.

Фрагмент фортранного текста вызывающей программы.
(полный текст теста можно получить в tdgeev1.zip, а принятые обозначения - посмотреть в разделе документации "Обозначения и упрощения в примерах по использованию подпрограмм комплекса").

Матрица A инициализируется посредством обращения к подпрограмме PDELSET из библиотеки ScaLAPACK_TOOLS. 

      PROGRAM TDGEEV1
      include 'mpif.h'
      INTEGER                         LWORK, MAXN, LDA, MAXPROCS, NOUT
      DOUBLE PRECISION    ZERO, MONE
      PARAMETER                  ( LDA = 20, LWORK = 231, MAXN = 100,
    $                                         MAXPROCS = 16, NOUT =6, ZERO = 0.0D+0,
    $                                         MONE = -1.0D+0 )
      INTEGER                        CTXT, I, IAM, INFO, MYCOL, MYROW, N, NB,
    $                                        NPCOL, NPROCS, NPROW, IA, JA
*
      INTEGER                        DESCA( 9 )
      DOUBLE PRECISION   A( LDA, LDA ), W( MAXN ),
    $                                        WORK( LWORK ), WORK1( 7 )
*
      EXTERNAL                    BLACS_EXIT, BLACS_GET, BLACS_GRIDEXIT,
    $                                        BLACS_GRIDINFO, BLACS_GRIDINIT, BLACS_PINFO,
    $                                        BLACS_SETUP, DESCINIT, PDMATINIT, PDLAPRNT,
    $                                        PDGEHRD, PDLAHQR
*
      N = 20
      NB = 7
      NPROW = 2
      NPCOL = 2
      IA = 1
      JA = 1
*
      CALL  BLACS_PINFO( IAM, NPROCS )
      IF( ( NPROCS .LT. 1 ) ) THEN
         CALL  BLACS_SETUP( IAM, NPROW*NPCOL )
      END IF
*
      CALL  BLACS_GET( -1, 0, CTXT )
      CALL  BLACS_GRIDINIT( CTXT, 'R', NPROW, NPCOL )
      CALL  BLACS_GRIDINFO( CTXT, NPROW, NPCOL, MYROW, MYCOL )

      IF( MYROW .EQ. -1 ) GO TO 20

      CALL  DESCINIT( DESCA, N, N, NB, NB, 0, 0, CTXT, LDA, INFO )

* Построение матрицы A

      CALL  PDMATINIT( N, A, IA, JA, DESCA, INFO )

* Вычисление собственных значений

      CALL  PDGEEV1( N, A, IA, JA, DESCA, WR, WI, WORK, LWORK, INFO )

      CALL  BLACS_GRIDEXIT( CTXT )
 20 CONTINUE

      CALL  BLACS_EXIT( 0 )
      STOP
      END
*
      SUBROUTINE  PDMATINIT( N, A, IA, JA, DESCA, INFO )
*
* PDMATINIT генерирует и распределяет матрицу A по решетке процессов
*
      DOUBLE PRECISION   ONE
      PARAMETER                 ( ONE = 1.0D+0 )
*
      INTEGER                        IA, INFO, JA, N
      INTEGER                        DESCA( * )
      DOUBLE PRECISION   A( * )
      INTEGER                        I, J
      EXTERNAL                    PDELSET
      INTRINSIC                     DBLE
*
      INFO = 0
*
      IF( IA .NE. 1 ) THEN
         INFO = -3
      ELSE IF( JA .NE. 1 ) THEN
         INFO = -4
      END IF
*
      DO 2 J = 1, N
      DO 1 I = 1, N
      CALL  PDELSET( A, I, J, DESCA, 0.0D0)
   1 CONTINUE
   2 CONTINUE
      DO 4 J = 1, N
      DO 3 I = 1, N
         IF( I .EQ. J ) THEN
            CALL  PDELSET( A, I, J, DESCA, 1.0D0)
         END IF
   3 CONTINUE
   4 CONTINUE
      J = N
      DO 5 I = 1, N
      CALL  PDELSET( A, I, J, DESCA, DBLE(I) )
   5 CONTINUE
      I = N         
      DO 6 J = 1, N
      CALL  PDELSET( A, I, J, DESCA, DBLE(J) )
   6 CONTINUE
*
      RETURN
      END

Результаты:

 Значение  INFO  =   0

 Собственные значения:

  W(1) = -40.0989D+00
  W(2) =  61.0989D+00
  W(3) =   1.D+00
  W(4) =   1.D+00
  W(5) =   1.D+00
  W(6) =   1.D+00
  W(7) =   1.D+00
  W(8) =   1.D+00
  W(9) =   1.D+00
  W(10) =  1.D+00
  W(11) =  1.D+00
  W(12) =  1.D+00
  W(13) =  1.D+00
  W(14) =  1.D+00
  W(15) =  1.D+00
  W(16) =  1.D+00
  W(17) =  1.D+00
  W(18) =  1.D+00
  W(19) =  1.D+00
  W(20) =  1.D+00