1 static char help[] = "Tests MatConvert(), MatLoad() for MATSCALAPACK interface.\n\n";
2 /*
3 Example:
4 mpiexec -n <np> ./ex244 -fA <A_data> -fB <B_data> -orig_mat_type <type> -orig_mat_type <mat_type>
5 */
6
7 #include <petscmat.h>
8
main(int argc,char ** args)9 int main(int argc, char **args)
10 {
11 Mat A, Ae, B, Be;
12 PetscViewer view;
13 char file[2][PETSC_MAX_PATH_LEN];
14 PetscBool flg, flgB, isScaLAPACK, isDense, isAij, isSbaij;
15 PetscScalar one = 1.0;
16 PetscMPIInt rank, size;
17 PetscInt M, N;
18
19 PetscFunctionBeginUser;
20 PetscCall(PetscInitialize(&argc, &args, nullptr, help));
21 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
22 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
23
24 /* Load PETSc matrices */
25 PetscCall(PetscOptionsGetString(NULL, NULL, "-fA", file[0], PETSC_MAX_PATH_LEN, NULL));
26 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, file[0], FILE_MODE_READ, &view));
27 PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
28 PetscCall(MatSetOptionsPrefix(A, "orig_"));
29 PetscCall(MatSetType(A, MATAIJ));
30 PetscCall(MatSetFromOptions(A));
31 PetscCall(MatLoad(A, view));
32 PetscCall(PetscViewerDestroy(&view));
33
34 PetscCall(PetscOptionsGetString(NULL, NULL, "-fB", file[1], PETSC_MAX_PATH_LEN, &flgB));
35 if (flgB) {
36 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, file[1], FILE_MODE_READ, &view));
37 PetscCall(MatCreate(PETSC_COMM_WORLD, &B));
38 PetscCall(MatSetOptionsPrefix(B, "orig_"));
39 PetscCall(MatSetType(B, MATAIJ));
40 PetscCall(MatSetFromOptions(B));
41 PetscCall(MatLoad(B, view));
42 PetscCall(PetscViewerDestroy(&view));
43 } else {
44 /* Create matrix B = I */
45 PetscInt rstart, rend, i;
46 PetscCall(MatGetSize(A, &M, &N));
47 PetscCall(MatGetOwnershipRange(A, &rstart, &rend));
48
49 PetscCall(MatCreate(PETSC_COMM_WORLD, &B));
50 PetscCall(MatSetOptionsPrefix(B, "orig_"));
51 PetscCall(MatSetSizes(B, PETSC_DECIDE, PETSC_DECIDE, M, N));
52 PetscCall(MatSetType(B, MATAIJ));
53 PetscCall(MatSetFromOptions(B));
54 PetscCall(MatSetUp(B));
55 for (i = rstart; i < rend; i++) PetscCall(MatSetValues(B, 1, &i, 1, &i, &one, ADD_VALUES));
56 PetscCall(MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY));
57 PetscCall(MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY));
58 }
59
60 PetscCall(PetscObjectTypeCompare((PetscObject)A, MATSCALAPACK, &isScaLAPACK));
61 if (isScaLAPACK) {
62 Ae = A;
63 Be = B;
64 isDense = isAij = isSbaij = PETSC_FALSE;
65 } else { /* Convert AIJ/DENSE/SBAIJ matrices into ScaLAPACK matrices */
66 if (size == 1) {
67 PetscCall(PetscObjectTypeCompare((PetscObject)A, MATSEQDENSE, &isDense));
68 PetscCall(PetscObjectTypeCompare((PetscObject)A, MATSEQAIJ, &isAij));
69 PetscCall(PetscObjectTypeCompare((PetscObject)A, MATSEQSBAIJ, &isSbaij));
70 } else {
71 PetscCall(PetscObjectTypeCompare((PetscObject)A, MATMPIDENSE, &isDense));
72 PetscCall(PetscObjectTypeCompare((PetscObject)A, MATMPIAIJ, &isAij));
73 PetscCall(PetscObjectTypeCompare((PetscObject)A, MATMPISBAIJ, &isSbaij));
74 }
75
76 if (rank == 0) {
77 if (isDense) {
78 printf(" Convert DENSE matrices A and B into ScaLAPACK matrix... \n");
79 } else if (isAij) {
80 printf(" Convert AIJ matrices A and B into ScaLAPACK matrix... \n");
81 } else if (isSbaij) {
82 printf(" Convert SBAIJ matrices A and B into ScaLAPACK matrix... \n");
83 } else SETERRQ(PetscObjectComm((PetscObject)A), PETSC_ERR_SUP, "Not supported yet");
84 }
85 PetscCall(MatConvert(A, MATSCALAPACK, MAT_INITIAL_MATRIX, &Ae));
86 PetscCall(MatConvert(B, MATSCALAPACK, MAT_INITIAL_MATRIX, &Be));
87
88 /* Test accuracy */
89 PetscCall(MatMultEqual(A, Ae, 5, &flg));
90 PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_ARG_NOTSAMETYPE, "A != A_scalapack.");
91 PetscCall(MatMultEqual(B, Be, 5, &flg));
92 PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_ARG_NOTSAMETYPE, "B != B_scalapack.");
93 }
94
95 if (!isScaLAPACK) {
96 PetscCall(MatDestroy(&Ae));
97 PetscCall(MatDestroy(&Be));
98
99 /* Test MAT_REUSE_MATRIX which is only supported for inplace conversion */
100 PetscCall(MatConvert(A, MATSCALAPACK, MAT_INPLACE_MATRIX, &A));
101 //PetscCall(MatView(A,PETSC_VIEWER_STDOUT_WORLD));
102 }
103
104 PetscCall(MatDestroy(&A));
105 PetscCall(MatDestroy(&B));
106 PetscCall(PetscFinalize());
107 return 0;
108 }
109
110 /*TEST
111
112 build:
113 requires: scalapack
114
115 test:
116 requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
117 args: -fA ${DATAFILESPATH}/matrices/EigenProblems/Eigdftb/dftb_bin/graphene_xxs_A_aij -fB ${DATAFILESPATH}/matrices/EigenProblems/Eigdftb/dftb_bin/graphene_xxs_B_aij
118 output_file: output/ex244.out
119
120 test:
121 suffix: 2
122 nsize: 8
123 requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
124 args: -fA ${DATAFILESPATH}/matrices/EigenProblems/Eigdftb/dftb_bin/graphene_xxs_A_aij -fB ${DATAFILESPATH}/matrices/EigenProblems/Eigdftb/dftb_bin/graphene_xxs_B_aij
125 output_file: output/ex244.out
126
127 test:
128 suffix: 2_dense
129 nsize: 8
130 requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) !defined(PETSCTEST_VALGRIND)
131 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
132 output_file: output/ex244_dense.out
133
134 test:
135 suffix: 2_scalapack
136 nsize: 8
137 requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) !defined(PETSCTEST_VALGRIND)
138 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
139 output_file: output/empty.out
140
141 test:
142 suffix: 2_sbaij
143 nsize: 8
144 requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
145 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
146 output_file: output/ex244_sbaij.out
147
148 test:
149 suffix: complex
150 requires: complex double datafilespath !defined(PETSC_USE_64BIT_INDICES)
151 args: -fA ${DATAFILESPATH}/matrices/nimrod/small_112905
152 output_file: output/ex244.out
153
154 test:
155 suffix: complex_2
156 nsize: 4
157 requires: complex double datafilespath !defined(PETSC_USE_64BIT_INDICES)
158 args: -fA ${DATAFILESPATH}/matrices/nimrod/small_112905
159 output_file: output/ex244.out
160
161 test:
162 suffix: dense
163 requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
164 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
165
166 test:
167 suffix: scalapack
168 requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
169 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
170 output_file: output/empty.out
171
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
176
177 TEST*/
178