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;
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],sizeof(file[0]),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],sizeof(file[1]),&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++) {
57: MatSetValues(B,1,&i,1,&i,&one,ADD_VALUES);
58: }
59: MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY);
60: MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY);
61: }
63: PetscObjectTypeCompare((PetscObject)A,MATELEMENTAL,&isElemental);
64: if (isElemental) {
65: Ae = A;
66: Be = B;
67: isDense = isAij = isSbaij = PETSC_FALSE;
68: } else { /* Convert AIJ/DENSE/SBAIJ matrices into Elemental matrices */
69: if (size == 1) {
70: PetscObjectTypeCompare((PetscObject)A,MATSEQDENSE,&isDense);
71: PetscObjectTypeCompare((PetscObject)A,MATSEQAIJ,&isAij);
72: PetscObjectTypeCompare((PetscObject)A,MATSEQSBAIJ,&isSbaij);
73: } else {
74: PetscObjectTypeCompare((PetscObject)A,MATMPIDENSE,&isDense);
75: PetscObjectTypeCompare((PetscObject)A,MATMPIAIJ,&isAij);
76: PetscObjectTypeCompare((PetscObject)A,MATMPISBAIJ,&isSbaij);
77: }
79: if (rank == 0) {
80: if (isDense) {
81: printf(" Convert DENSE matrices A and B into Elemental matrix... \n");
82: } else if (isAij) {
83: printf(" Convert AIJ matrices A and B into Elemental matrix... \n");
84: } else if (isSbaij) {
85: printf(" Convert SBAIJ matrices A and B into Elemental matrix... \n");
86: } else SETERRQ(PetscObjectComm((PetscObject)A),PETSC_ERR_SUP,"Not supported yet");
87: }
88: MatConvert(A, MATELEMENTAL, MAT_INITIAL_MATRIX, &Ae);
89: MatConvert(B, MATELEMENTAL, MAT_INITIAL_MATRIX, &Be);
91: /* Test accuracy */
92: MatMultEqual(A,Ae,5,&flg);
94: MatMultEqual(B,Be,5,&flg);
96: }
98: if (!isElemental) {
99: MatDestroy(&Ae);
100: MatDestroy(&Be);
102: /* Test MAT_REUSE_MATRIX which is only supported for inplace conversion */
103: MatConvert(A, MATELEMENTAL, MAT_INPLACE_MATRIX, &A);
104: //MatView(A,PETSC_VIEWER_STDOUT_WORLD);
105: }
107: MatDestroy(&A);
108: MatDestroy(&B);
109: PetscFinalize();
110: return 0;
111: }
113: /*TEST
115: build:
116: requires: elemental
118: test:
119: requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
120: args: -fA ${DATAFILESPATH}/matrices/EigenProblems/Eigdftb/dftb_bin/graphene_xxs_A_aij -fB ${DATAFILESPATH}/matrices/EigenProblems/Eigdftb/dftb_bin/graphene_xxs_B_aij
121: output_file: output/ex174.out
123: test:
124: suffix: 2
125: nsize: 8
126: requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
127: args: -fA ${DATAFILESPATH}/matrices/EigenProblems/Eigdftb/dftb_bin/graphene_xxs_A_aij -fB ${DATAFILESPATH}/matrices/EigenProblems/Eigdftb/dftb_bin/graphene_xxs_B_aij
128: output_file: output/ex174.out
130: test:
131: suffix: 2_dense
132: nsize: 8
133: requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
134: 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
135: output_file: output/ex174_dense.out
137: test:
138: suffix: 2_elemental
139: nsize: 8
140: requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
141: 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
142: output_file: output/ex174_elemental.out
144: test:
145: suffix: 2_sbaij
146: nsize: 8
147: requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
148: 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
149: output_file: output/ex174_sbaij.out
151: test:
152: suffix: complex
153: requires: complex double datafilespath !defined(PETSC_USE_64BIT_INDICES)
154: args: -fA ${DATAFILESPATH}/matrices/nimrod/small_112905
155: output_file: output/ex174.out
157: test:
158: suffix: complex_2
159: nsize: 4
160: requires: complex double datafilespath !defined(PETSC_USE_64BIT_INDICES)
161: args: -fA ${DATAFILESPATH}/matrices/nimrod/small_112905
162: output_file: output/ex174.out
164: test:
165: suffix: dense
166: requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
167: 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
169: test:
170: suffix: elemental
171: requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
172: 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
174: test:
175: suffix: sbaij
176: requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
177: 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
179: TEST*/