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