pzgebrd_ge2gb.c 3.94 KB
Newer Older
1
/**
2 3
 *
 * @file pzgebrd_ge2gb.c
4
 *
Mathieu Faverge's avatar
Mathieu Faverge committed
5 6
 * @copyright 2009-2014 The University of Tennessee and The University of
 *                      Tennessee Research Foundation. All rights reserved.
7 8
 * @copyright 2012-2014 Bordeaux INP, CNRS (LaBRI UMR 5800), Inria,
 *                      Univ. Bordeaux. All rights reserved.
9
 *
Mathieu Faverge's avatar
Mathieu Faverge committed
10
 ***
11
 *
Mathieu Faverge's avatar
Mathieu Faverge committed
12
 * @brief Chameleon zgebrd_ge2gb parallel algorithm
13
 *
Mathieu Faverge's avatar
Mathieu Faverge committed
14
 * @version 1.0.0
15 16 17 18 19
 * @author Hatem Ltaief
 * @author Azzam Haidar
 * @date 2010-11-15
 * @precisions normal z -> s d c
 *
20
 */
21 22
#include "control/common.h"

23
void morse_pzgebrd_ge2gb(MORSE_desc_t *A, MORSE_desc_t *T, MORSE_desc_t *D,
24 25 26 27
                         MORSE_sequence_t *sequence, MORSE_request_t *request)
{
    int k;
    int tempkm, tempkn;
28
    MORSE_desc_t *A1, *A2, *T1, *D1 = NULL;
29

30 31 32 33
    if (A->m >= A->n){
       for (k = 0; k < A->nt; k++) {
           tempkm = k == A->mt-1 ? A->m-k*A->mb : A->mb;
           tempkn = k == A->nt-1 ? A->n-k*A->nb : A->nb;
34

35 36 37 38 39 40 41 42 43 44 45 46 47
           A1 = morse_desc_submatrix(A, k*A->mb,     k*A->nb, A->m-k*A->mb, tempkn);
           A2 = morse_desc_submatrix(A, k*A->mb, (k+1)*A->nb, A->m-k*A->mb, A->n-(k+1)*A->nb);
           T1 = morse_desc_submatrix(T, k*T->mb,     k*T->nb, T->m-k*T->mb, tempkn);
           if ( D != NULL ) {
               D1 = morse_desc_submatrix(D, k*D->mb,     k*D->nb, D->m-k*D->mb, tempkn);
           }

           morse_pzgeqrf( A1, T1, D1,
                          sequence, request);

           morse_pzunmqr( MorseLeft, MorseConjTrans,
                          A1, A2, T1, D1,
                          sequence, request);
48

49 50
           if (k+1 < A->nt){
              tempkn = k+1 == A->nt-1 ? A->n-(k+1)*A->nb : A->nb;
51

52 53 54 55 56 57
              A1 = morse_desc_submatrix(A,     k*A->mb, (k+1)*A->nb, tempkm,           A->n-(k+1)*A->nb);
              A2 = morse_desc_submatrix(A, (k+1)*A->mb, (k+1)*A->nb, A->m-(k+1)*A->mb, A->n-(k+1)*A->nb);
              T1 = morse_desc_submatrix(T,     k*T->mb, (k+1)*T->nb, tempkm,           T->n-(k+1)*T->nb);
              if ( D != NULL ) {
                  D1 = morse_desc_submatrix(D, k*D->mb, (k+1)*D->nb, tempkm,           D->n-(k+1)*D->nb);
              }
58

59 60 61 62 63 64
              morse_pzgelqf( A1, T1, D1,
                             sequence, request);

              morse_pzunmlq( MorseRight, MorseConjTrans,
                             A1, A2, T1, D1,
                             sequence, request);
65 66 67 68
           }
       }
    }
    else{
69 70 71 72 73 74 75 76 77 78 79 80
       for (k = 0; k < A->mt; k++) {
           tempkm = k == A->mt-1 ? A->m-k*A->mb : A->mb;
           tempkn = k == A->nt-1 ? A->n-k*A->nb : A->nb;

           A1 = morse_desc_submatrix(A,     k*A->mb, k*A->nb, tempkm,           A->n-k*A->nb);
           A2 = morse_desc_submatrix(A, (k+1)*A->mb, k*A->nb, A->m-(k+1)*A->mb, A->n-k*A->nb);
           T1 = morse_desc_submatrix(T,     k*T->mb, k*T->nb, tempkm,           T->n-k*T->nb);
           if ( D != NULL ) {
               D1 = morse_desc_submatrix(D, k*D->mb, k*D->nb, tempkm,           D->n-k*D->nb);
           }
           morse_pzgelqf( A1, T1, D1,
                          sequence, request);
81

82 83 84
           morse_pzunmlq( MorseRight, MorseConjTrans,
                          A1, A2, T1, D1,
                          sequence, request);
85

86 87
           if (k+1 < A->mt){
              tempkm = k+1 == A->mt-1 ? A->m-(k+1)*A->mb : A->mb;
88

89 90 91 92 93 94
              A1 = morse_desc_submatrix(A, (k+1)*A->mb,     k*A->nb, A->m-(k+1)*A->mb, tempkn);
              A2 = morse_desc_submatrix(A, (k+1)*A->mb, (k+1)*A->nb, A->m-(k+1)*A->mb, A->n-(k+1)*A->nb);
              T1 = morse_desc_submatrix(T, (k+1)*T->mb,     k*T->nb, T->m-(k+1)*T->mb, tempkn);
              if ( D != NULL ) {
                  D1 = morse_desc_submatrix(D, (k+1)*D->mb, k*D->nb, D->m-(k+1)*D->mb, tempkn);
              }
95

96 97
              morse_pzgeqrf( A1, T1, D1,
                             sequence, request);
98

99 100 101
              morse_pzunmqr( MorseLeft, MorseConjTrans,
                             A1, A2, T1, D1,
                             sequence, request);
102 103 104 105
           }
       }
    }
}