Actual source code: ex244.cxx


  2: static char help[] = "Tests MatConvert(), MatLoad() for MATSCALAPACK interface.\n\n";
  3: /*
  4:  Example:
  5:    mpiexec -n <np> ./ex244 -fA <A_data> -fB <B_data> -orig_mat_type <type> -orig_mat_type <mat_type>
  6: */

  8: #include <petscmat.h>

 10: int main(int argc, char **args)
 11: {
 12:   Mat         A, Ae, B, Be;
 13:   PetscViewer view;
 14:   char        file[2][PETSC_MAX_PATH_LEN];
 15:   PetscBool   flg, flgB, isScaLAPACK, isDense, isAij, isSbaij;
 16:   PetscScalar one = 1.0;
 17:   PetscMPIInt rank, size;
 18:   PetscInt    M, N;

 21:   PetscInitialize(&argc, &args, (char *)0, help);
 22:   MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
 23:   MPI_Comm_size(PETSC_COMM_WORLD, &size);

 25:   /* Load PETSc matrices */
 26:   PetscOptionsGetString(NULL, NULL, "-fA", file[0], PETSC_MAX_PATH_LEN, NULL);
 27:   PetscViewerBinaryOpen(PETSC_COMM_WORLD, file[0], FILE_MODE_READ, &view);
 28:   MatCreate(PETSC_COMM_WORLD, &A);
 29:   MatSetOptionsPrefix(A, "orig_");
 30:   MatSetType(A, MATAIJ);
 31:   MatSetFromOptions(A);
 32:   MatLoad(A, view);
 33:   PetscViewerDestroy(&view);

 35:   PetscOptionsGetString(NULL, NULL, "-fB", file[1], PETSC_MAX_PATH_LEN, &flgB);
 36:   if (flgB) {
 37:     PetscViewerBinaryOpen(PETSC_COMM_WORLD, file[1], FILE_MODE_READ, &view);
 38:     MatCreate(PETSC_COMM_WORLD, &B);
 39:     MatSetOptionsPrefix(B, "orig_");
 40:     MatSetType(B, MATAIJ);
 41:     MatSetFromOptions(B);
 42:     MatLoad(B, view);
 43:     PetscViewerDestroy(&view);
 44:   } else {
 45:     /* Create matrix B = I */
 46:     PetscInt rstart, rend, i;
 47:     MatGetSize(A, &M, &N);
 48:     MatGetOwnershipRange(A, &rstart, &rend);

 50:     MatCreate(PETSC_COMM_WORLD, &B);
 51:     MatSetOptionsPrefix(B, "orig_");
 52:     MatSetSizes(B, PETSC_DECIDE, PETSC_DECIDE, M, N);
 53:     MatSetType(B, MATAIJ);
 54:     MatSetFromOptions(B);
 55:     MatSetUp(B);
 56:     for (i = rstart; i < rend; i++) MatSetValues(B, 1, &i, 1, &i, &one, ADD_VALUES);
 57:     MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
 58:     MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
 59:   }

 61:   PetscObjectTypeCompare((PetscObject)A, MATSCALAPACK, &isScaLAPACK);
 62:   if (isScaLAPACK) {
 63:     Ae      = A;
 64:     Be      = B;
 65:     isDense = isAij = isSbaij = PETSC_FALSE;
 66:   } else { /* Convert AIJ/DENSE/SBAIJ matrices into ScaLAPACK matrices */
 67:     if (size == 1) {
 68:       PetscObjectTypeCompare((PetscObject)A, MATSEQDENSE, &isDense);
 69:       PetscObjectTypeCompare((PetscObject)A, MATSEQAIJ, &isAij);
 70:       PetscObjectTypeCompare((PetscObject)A, MATSEQSBAIJ, &isSbaij);
 71:     } else {
 72:       PetscObjectTypeCompare((PetscObject)A, MATMPIDENSE, &isDense);
 73:       PetscObjectTypeCompare((PetscObject)A, MATMPIAIJ, &isAij);
 74:       PetscObjectTypeCompare((PetscObject)A, MATMPISBAIJ, &isSbaij);
 75:     }

 77:     if (rank == 0) {
 78:       if (isDense) {
 79:         printf(" Convert DENSE matrices A and B into ScaLAPACK matrix... \n");
 80:       } else if (isAij) {
 81:         printf(" Convert AIJ matrices A and B into ScaLAPACK matrix... \n");
 82:       } else if (isSbaij) {
 83:         printf(" Convert SBAIJ matrices A and B into ScaLAPACK matrix... \n");
 84:       } else SETERRQ(PetscObjectComm((PetscObject)A), PETSC_ERR_SUP, "Not supported yet");
 85:     }
 86:     MatConvert(A, MATSCALAPACK, MAT_INITIAL_MATRIX, &Ae);
 87:     MatConvert(B, MATSCALAPACK, MAT_INITIAL_MATRIX, &Be);

 89:     /* Test accuracy */
 90:     MatMultEqual(A, Ae, 5, &flg);
 92:     MatMultEqual(B, Be, 5, &flg);
 94:   }

 96:   if (!isScaLAPACK) {
 97:     MatDestroy(&Ae);
 98:     MatDestroy(&Be);

100:     /* Test MAT_REUSE_MATRIX which is only supported for inplace conversion */
101:     MatConvert(A, MATSCALAPACK, MAT_INPLACE_MATRIX, &A);
102:     //MatView(A,PETSC_VIEWER_STDOUT_WORLD);
103:   }

105:   MatDestroy(&A);
106:   MatDestroy(&B);
107:   PetscFinalize();
108:   return 0;
109: }

111: /*TEST

113:    build:
114:       requires: scalapack

116:    test:
117:       requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
118:       args: -fA ${DATAFILESPATH}/matrices/EigenProblems/Eigdftb/dftb_bin/graphene_xxs_A_aij -fB ${DATAFILESPATH}/matrices/EigenProblems/Eigdftb/dftb_bin/graphene_xxs_B_aij
119:       output_file: output/ex244.out

121:    test:
122:       suffix: 2
123:       nsize: 8
124:       requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
125:       args: -fA ${DATAFILESPATH}/matrices/EigenProblems/Eigdftb/dftb_bin/graphene_xxs_A_aij -fB ${DATAFILESPATH}/matrices/EigenProblems/Eigdftb/dftb_bin/graphene_xxs_B_aij
126:       output_file: output/ex244.out

128:    test:
129:       suffix: 2_dense
130:       nsize: 8
131:       requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
132:       args: -fA ${DATAFILESPATH}/matrices/EigenProblems/Eigdftb/dftb_bin/graphene_xxs_A_aij -fB ${DATAFILESPATH}/matrices/EigenProblems/Eigdftb/dftb_bin/graphene_xxs_B_aij -orig_mat_type dense
133:       output_file: output/ex244_dense.out

135:    test:
136:       suffix: 2_scalapack
137:       nsize: 8
138:       requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
139:       args: -fA ${DATAFILESPATH}/matrices/EigenProblems/Eigdftb/dftb_bin/graphene_xxs_A_aij -fB ${DATAFILESPATH}/matrices/EigenProblems/Eigdftb/dftb_bin/graphene_xxs_B_aij -orig_mat_type scalapack
140:       output_file: output/ex244_scalapack.out

142:    test:
143:       suffix: 2_sbaij
144:       nsize: 8
145:       requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
146:       args: -fA ${DATAFILESPATH}/matrices/EigenProblems/Eigdftb/dftb_bin/graphene_xxs_A -fB ${DATAFILESPATH}/matrices/EigenProblems/Eigdftb/dftb_bin/graphene_xxs_B -orig_mat_type sbaij
147:       output_file: output/ex244_sbaij.out

149:    test:
150:       suffix: complex
151:       requires: complex double datafilespath !defined(PETSC_USE_64BIT_INDICES)
152:       args: -fA ${DATAFILESPATH}/matrices/nimrod/small_112905
153:       output_file: output/ex244.out

155:    test:
156:       suffix: complex_2
157:       nsize: 4
158:       requires: complex double datafilespath !defined(PETSC_USE_64BIT_INDICES)
159:       args: -fA ${DATAFILESPATH}/matrices/nimrod/small_112905
160:       output_file: output/ex244.out

162:    test:
163:       suffix: dense
164:       requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
165:       args: -fA ${DATAFILESPATH}/matrices/EigenProblems/Eigdftb/dftb_bin/graphene_xxs_A_aij -fB ${DATAFILESPATH}/matrices/EigenProblems/Eigdftb/dftb_bin/graphene_xxs_B_aij -orig_mat_type dense

167:    test:
168:       suffix: scalapack
169:       requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
170:       args: -fA ${DATAFILESPATH}/matrices/EigenProblems/Eigdftb/dftb_bin/graphene_xxs_A_aij -fB ${DATAFILESPATH}/matrices/EigenProblems/Eigdftb/dftb_bin/graphene_xxs_B_aij -orig_mat_type scalapack

172:    test:
173:       suffix: sbaij
174:       requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
175:       args: -fA ${DATAFILESPATH}/matrices/EigenProblems/Eigdftb/dftb_bin/graphene_xxs_A -fB ${DATAFILESPATH}/matrices/EigenProblems/Eigdftb/dftb_bin/graphene_xxs_B -orig_mat_type sbaij

177: TEST*/