1 2 static char help[] = "Tests various routines in MatKAIJ format.\n"; 3 4 #include <petscmat.h> 5 #define IMAX 15 6 7 int main(int argc,char **args) 8 { 9 Mat A,B,TA; 10 PetscScalar *S,*T; 11 PetscViewer fd; 12 char file[PETSC_MAX_PATH_LEN]; 13 PetscInt m,n,M,N,p=1,q=1,i,j; 14 PetscMPIInt rank,size; 15 PetscBool flg; 16 17 CHKERRQ(PetscInitialize(&argc,&args,(char*)0,help)); 18 CHKERRMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 19 CHKERRMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size)); 20 21 /* Load AIJ matrix A */ 22 CHKERRQ(PetscOptionsGetString(NULL,NULL,"-f",file,sizeof(file),&flg)); 23 PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_USER,"Must indicate binary file with the -f option"); 24 CHKERRQ(PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&fd)); 25 CHKERRQ(MatCreate(PETSC_COMM_WORLD,&A)); 26 CHKERRQ(MatLoad(A,fd)); 27 CHKERRQ(PetscViewerDestroy(&fd)); 28 29 /* Get dof, then create S and T */ 30 CHKERRQ(PetscOptionsGetInt(NULL,NULL,"-p",&p,NULL)); 31 CHKERRQ(PetscOptionsGetInt(NULL,NULL,"-q",&q,NULL)); 32 CHKERRQ(PetscMalloc2(p*q,&S,p*q,&T)); 33 for (i=0; i<p*q; i++) S[i] = 0; 34 35 for (i=0; i<p; i++) { 36 for (j=0; j<q; j++) { 37 /* Set some random non-zero values */ 38 S[i+p*j] = ((PetscReal) ((i+1)*(j+1))) / ((PetscReal) (p+q)); 39 T[i+p*j] = ((PetscReal) ((p-i)+j)) / ((PetscReal) (p*q)); 40 } 41 } 42 43 /* Test KAIJ when both S & T are not NULL */ 44 45 /* Create KAIJ matrix TA */ 46 CHKERRQ(MatCreateKAIJ(A,p,q,S,T,&TA)); 47 CHKERRQ(MatGetLocalSize(A,&m,&n)); 48 CHKERRQ(MatGetSize(A,&M,&N)); 49 50 CHKERRQ(MatConvert(TA,MATAIJ,MAT_INITIAL_MATRIX,&B)); 51 52 /* Test MatKAIJGetScaledIdentity() */ 53 CHKERRQ(MatKAIJGetScaledIdentity(TA,&flg)); 54 PetscCheck(!flg,PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Error in Test 1: MatKAIJGetScaledIdentity()"); 55 /* Test MatMult() */ 56 CHKERRQ(MatMultEqual(TA,B,10,&flg)); 57 PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_CONV_FAILED,"Error in Test 1: MatMult() for KAIJ matrix"); 58 /* Test MatMultAdd() */ 59 CHKERRQ(MatMultAddEqual(TA,B,10,&flg)); 60 PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_CONV_FAILED,"Error in Test 1: MatMultAdd() for KAIJ matrix"); 61 62 CHKERRQ(MatDestroy(&TA)); 63 CHKERRQ(MatDestroy(&B)); 64 65 /* Test KAIJ when S is NULL */ 66 67 /* Create KAIJ matrix TA */ 68 CHKERRQ(MatCreateKAIJ(A,p,q,NULL,T,&TA)); 69 CHKERRQ(MatGetLocalSize(A,&m,&n)); 70 CHKERRQ(MatGetSize(A,&M,&N)); 71 72 CHKERRQ(MatConvert(TA,MATAIJ,MAT_INITIAL_MATRIX,&B)); 73 74 /* Test MatKAIJGetScaledIdentity() */ 75 CHKERRQ(MatKAIJGetScaledIdentity(TA,&flg)); 76 PetscCheck(!flg,PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Error in Test 2: MatKAIJGetScaledIdentity()"); 77 /* Test MatMult() */ 78 CHKERRQ(MatMultEqual(TA,B,10,&flg)); 79 PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_CONV_FAILED,"Error in Test 2: MatMult() for KAIJ matrix"); 80 /* Test MatMultAdd() */ 81 CHKERRQ(MatMultAddEqual(TA,B,10,&flg)); 82 PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_CONV_FAILED,"Error in Test 2: MatMultAdd() for KAIJ matrix"); 83 84 CHKERRQ(MatDestroy(&TA)); 85 CHKERRQ(MatDestroy(&B)); 86 87 /* Test KAIJ when T is NULL */ 88 89 /* Create KAIJ matrix TA */ 90 CHKERRQ(MatCreateKAIJ(A,p,q,S,NULL,&TA)); 91 CHKERRQ(MatGetLocalSize(A,&m,&n)); 92 CHKERRQ(MatGetSize(A,&M,&N)); 93 94 CHKERRQ(MatConvert(TA,MATAIJ,MAT_INITIAL_MATRIX,&B)); 95 96 /* Test MatKAIJGetScaledIdentity() */ 97 CHKERRQ(MatKAIJGetScaledIdentity(TA,&flg)); 98 PetscCheck(!flg,PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Error in Test 3: MatKAIJGetScaledIdentity()"); 99 /* Test MatMult() */ 100 CHKERRQ(MatMultEqual(TA,B,10,&flg)); 101 PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_CONV_FAILED,"Error in Test 3: MatMult() for KAIJ matrix"); 102 /* Test MatMultAdd() */ 103 CHKERRQ(MatMultAddEqual(TA,B,10,&flg)); 104 PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_CONV_FAILED,"Error in Test 3: MatMultAdd() for KAIJ matrix"); 105 106 CHKERRQ(MatDestroy(&TA)); 107 CHKERRQ(MatDestroy(&B)); 108 109 /* Test KAIJ when T is is an identity matrix */ 110 111 if (p == q) { 112 for (i=0; i<p; i++) { 113 for (j=0; j<q; j++) { 114 if (i==j) T[i+j*p] = 1.0; 115 else T[i+j*p] = 0.0; 116 } 117 } 118 119 /* Create KAIJ matrix TA */ 120 CHKERRQ(MatCreateKAIJ(A,p,q,S,T,&TA)); 121 CHKERRQ(MatGetLocalSize(A,&m,&n)); 122 CHKERRQ(MatGetSize(A,&M,&N)); 123 124 CHKERRQ(MatConvert(TA,MATAIJ,MAT_INITIAL_MATRIX,&B)); 125 126 /* Test MatKAIJGetScaledIdentity() */ 127 CHKERRQ(MatKAIJGetScaledIdentity(TA,&flg)); 128 PetscCheck(!flg,PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Error in Test 4: MatKAIJGetScaledIdentity()"); 129 /* Test MatMult() */ 130 CHKERRQ(MatMultEqual(TA,B,10,&flg)); 131 PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_CONV_FAILED,"Error in Test 4: MatMult() for KAIJ matrix"); 132 /* Test MatMultAdd() */ 133 CHKERRQ(MatMultAddEqual(TA,B,10,&flg)); 134 PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_CONV_FAILED,"Error in Test 4: MatMultAdd() for KAIJ matrix"); 135 136 CHKERRQ(MatDestroy(&TA)); 137 CHKERRQ(MatDestroy(&B)); 138 139 CHKERRQ(MatCreateKAIJ(A,p,q,NULL,T,&TA)); 140 /* Test MatKAIJGetScaledIdentity() */ 141 CHKERRQ(MatKAIJGetScaledIdentity(TA,&flg)); 142 PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Error in Test 5: MatKAIJGetScaledIdentity()"); 143 CHKERRQ(MatDestroy(&TA)); 144 145 for (i=0; i<p; i++) { 146 for (j=0; j<q; j++) { 147 if (i==j) S[i+j*p] = T[i+j*p] = 2.0; 148 else S[i+j*p] = T[i+j*p] = 0.0; 149 } 150 } 151 CHKERRQ(MatCreateKAIJ(A,p,q,S,T,&TA)); 152 /* Test MatKAIJGetScaledIdentity() */ 153 CHKERRQ(MatKAIJGetScaledIdentity(TA,&flg)); 154 PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Error in Test 6: MatKAIJGetScaledIdentity()"); 155 CHKERRQ(MatDestroy(&TA)); 156 } 157 158 /* Done with all tests */ 159 160 CHKERRQ(PetscFree2(S,T)); 161 CHKERRQ(MatDestroy(&A)); 162 CHKERRQ(PetscFinalize()); 163 return 0; 164 } 165 166 /*TEST 167 168 build: 169 requires: !complex 170 171 test: 172 requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 173 output_file: output/ex176.out 174 nsize: {{1 2 3 4}} 175 args: -f ${DATAFILESPATH}/matrices/small -p {{2 3 7}} -q {{3 7}} -viewer_binary_skip_info 176 177 TEST*/ 178