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