1 2 #include "src/mat/matimpl.h" /*I "petscmat.h" I*/ 3 4 #undef __FUNCT__ 5 #define __FUNCT__ "MatMultEqual" 6 /*@ 7 MatMultEqual - Compares matrix-vector products of two matrices. 8 9 Collective on Mat 10 11 Input Parameters: 12 + A - the first matrix 13 - B - the second matrix 14 - n - number of random vectors to be tested 15 16 Output Parameter: 17 . flg - PETSC_TRUE if the products are equal; PETSC_FALSE otherwise. 18 19 Level: intermediate 20 21 Concepts: matrices^equality between 22 @*/ 23 PetscErrorCode MatMultEqual(Mat A,Mat B,PetscInt n,PetscTruth *flg) 24 { 25 PetscErrorCode ierr; 26 Vec x,s1,s2; 27 PetscRandom rctx; 28 PetscReal r1,r2,tol=1.e-10; 29 PetscInt am,an,bm,bn,k; 30 31 PetscFunctionBegin; 32 ierr = MatGetLocalSize(A,&am,&an);CHKERRQ(ierr); 33 ierr = MatGetLocalSize(B,&bm,&bn);CHKERRQ(ierr); 34 if (am != bm || an != bn) SETERRQ4(PETSC_ERR_ARG_SIZ,"Mat A,Mat B: local dim %D %D %D %D",am,bm,an,bn); 35 PetscCheckSameComm(A,1,B,2); 36 ierr = PetscRandomCreate(A->comm,RANDOM_DEFAULT,&rctx);CHKERRQ(ierr); 37 ierr = VecCreate(A->comm,&x);CHKERRQ(ierr); 38 ierr = VecSetSizes(x,an,PETSC_DECIDE);CHKERRQ(ierr); 39 ierr = VecSetFromOptions(x);CHKERRQ(ierr); 40 41 ierr = VecCreate(A->comm,&s1);CHKERRQ(ierr); 42 ierr = VecSetSizes(s1,am,PETSC_DECIDE);CHKERRQ(ierr); 43 ierr = VecSetFromOptions(s1);CHKERRQ(ierr); 44 ierr = VecDuplicate(s1,&s2);CHKERRQ(ierr); 45 46 *flg = PETSC_TRUE; 47 for (k=0; k<n; k++) { 48 ierr = VecSetRandom(rctx,x);CHKERRQ(ierr); 49 ierr = MatMult(A,x,s1);CHKERRQ(ierr); 50 ierr = MatMult(B,x,s2);CHKERRQ(ierr); 51 ierr = VecNorm(s1,NORM_1,&r1);CHKERRQ(ierr); 52 ierr = VecNorm(s2,NORM_1,&r2);CHKERRQ(ierr); 53 r1 -= r2; 54 if (r1<-tol || r1>tol) { 55 *flg = PETSC_FALSE; 56 ierr = PetscPrintf(PETSC_COMM_SELF,"Error: %d-th MatMult() %g\n",k,r1); 57 break; 58 } 59 } 60 ierr = PetscRandomDestroy(rctx);CHKERRQ(ierr); 61 ierr = VecDestroy(x);CHKERRQ(ierr); 62 ierr = VecDestroy(s1);CHKERRQ(ierr); 63 ierr = VecDestroy(s2);CHKERRQ(ierr); 64 PetscFunctionReturn(0); 65 } 66 67 #undef __FUNCT__ 68 #define __FUNCT__ "MatMultAddEqual" 69 /*@ 70 MatMultAddEqual - Compares matrix-vector products of two matrices. 71 72 Collective on Mat 73 74 Input Parameters: 75 + A - the first matrix 76 - B - the second matrix 77 - n - number of random vectors to be tested 78 79 Output Parameter: 80 . flg - PETSC_TRUE if the products are equal; PETSC_FALSE otherwise. 81 82 Level: intermediate 83 84 Concepts: matrices^equality between 85 @*/ 86 PetscErrorCode MatMultAddEqual(Mat A,Mat B,PetscInt n,PetscTruth *flg) 87 { 88 PetscErrorCode ierr; 89 Vec x,y,s1,s2; 90 PetscRandom rctx; 91 PetscReal r1,r2,tol=1.e-10; 92 PetscInt am,an,bm,bn,k; 93 94 PetscFunctionBegin; 95 ierr = MatGetLocalSize(A,&am,&an);CHKERRQ(ierr); 96 ierr = MatGetLocalSize(B,&bm,&bn);CHKERRQ(ierr); 97 if (am != bm || an != bn) SETERRQ4(PETSC_ERR_ARG_SIZ,"Mat A,Mat B: local dim %D %D %D %D",am,bm,an,bn); 98 PetscCheckSameComm(A,1,B,2); 99 ierr = PetscRandomCreate(A->comm,RANDOM_DEFAULT,&rctx);CHKERRQ(ierr); 100 ierr = VecCreate(A->comm,&x);CHKERRQ(ierr); 101 ierr = VecSetSizes(x,an,PETSC_DECIDE);CHKERRQ(ierr); 102 ierr = VecSetFromOptions(x);CHKERRQ(ierr); 103 104 ierr = VecCreate(A->comm,&s1);CHKERRQ(ierr); 105 ierr = VecSetSizes(s1,am,PETSC_DECIDE);CHKERRQ(ierr); 106 ierr = VecSetFromOptions(s1);CHKERRQ(ierr); 107 ierr = VecDuplicate(s1,&s2);CHKERRQ(ierr); 108 ierr = VecDuplicate(s1,&y);CHKERRQ(ierr); 109 110 *flg = PETSC_TRUE; 111 for (k=0; k<n; k++) { 112 ierr = VecSetRandom(rctx,x);CHKERRQ(ierr); 113 ierr = VecSetRandom(rctx,y);CHKERRQ(ierr); 114 ierr = MatMultAdd(A,x,y,s1);CHKERRQ(ierr); 115 ierr = MatMultAdd(B,x,y,s2);CHKERRQ(ierr); 116 ierr = VecNorm(s1,NORM_1,&r1);CHKERRQ(ierr); 117 ierr = VecNorm(s2,NORM_1,&r2);CHKERRQ(ierr); 118 r1 -= r2; 119 if (r1<-tol || r1>tol) { 120 *flg = PETSC_FALSE; 121 ierr = PetscPrintf(PETSC_COMM_SELF,"Error: %d-th MatMultAdd() %g\n",k,r1); 122 break; 123 } 124 } 125 ierr = PetscRandomDestroy(rctx);CHKERRQ(ierr); 126 ierr = VecDestroy(x);CHKERRQ(ierr); 127 ierr = VecDestroy(y);CHKERRQ(ierr); 128 ierr = VecDestroy(s1);CHKERRQ(ierr); 129 ierr = VecDestroy(s2);CHKERRQ(ierr); 130 PetscFunctionReturn(0); 131 } 132 133 #undef __FUNCT__ 134 #define __FUNCT__ "MatMultTransposeEqual" 135 /*@ 136 MatMultTransposeEqual - Compares matrix-vector products of two matrices. 137 138 Collective on Mat 139 140 Input Parameters: 141 + A - the first matrix 142 - B - the second matrix 143 - n - number of random vectors to be tested 144 145 Output Parameter: 146 . flg - PETSC_TRUE if the products are equal; PETSC_FALSE otherwise. 147 148 Level: intermediate 149 150 Concepts: matrices^equality between 151 @*/ 152 PetscErrorCode MatMultTransposeEqual(Mat A,Mat B,PetscInt n,PetscTruth *flg) 153 { 154 PetscErrorCode ierr; 155 Vec x,s1,s2; 156 PetscRandom rctx; 157 PetscReal r1,r2,tol=1.e-10; 158 PetscInt am,an,bm,bn,k; 159 160 PetscFunctionBegin; 161 ierr = MatGetLocalSize(A,&am,&an);CHKERRQ(ierr); 162 ierr = MatGetLocalSize(B,&bm,&bn);CHKERRQ(ierr); 163 if (am != bm || an != bn) SETERRQ4(PETSC_ERR_ARG_SIZ,"Mat A,Mat B: local dim %D %D %D %D",am,bm,an,bn); 164 PetscCheckSameComm(A,1,B,2); 165 ierr = PetscRandomCreate(A->comm,RANDOM_DEFAULT,&rctx);CHKERRQ(ierr); 166 ierr = VecCreate(A->comm,&x);CHKERRQ(ierr); 167 ierr = VecSetSizes(x,am,PETSC_DECIDE);CHKERRQ(ierr); 168 ierr = VecSetFromOptions(x);CHKERRQ(ierr); 169 170 ierr = VecCreate(A->comm,&s1);CHKERRQ(ierr); 171 ierr = VecSetSizes(s1,an,PETSC_DECIDE);CHKERRQ(ierr); 172 ierr = VecSetFromOptions(s1);CHKERRQ(ierr); 173 ierr = VecDuplicate(s1,&s2);CHKERRQ(ierr); 174 175 *flg = PETSC_TRUE; 176 for (k=0; k<n; k++) { 177 ierr = VecSetRandom(rctx,x);CHKERRQ(ierr); 178 ierr = MatMultTranspose(A,x,s1);CHKERRQ(ierr); 179 ierr = MatMultTranspose(B,x,s2);CHKERRQ(ierr); 180 ierr = VecNorm(s1,NORM_1,&r1);CHKERRQ(ierr); 181 ierr = VecNorm(s2,NORM_1,&r2);CHKERRQ(ierr); 182 r1 -= r2; 183 if (r1<-tol || r1>tol) { 184 *flg = PETSC_FALSE; 185 ierr = PetscPrintf(PETSC_COMM_SELF,"Error: %d-th MatMultTranspose() %g\n",k,r1); 186 break; 187 } 188 } 189 ierr = PetscRandomDestroy(rctx);CHKERRQ(ierr); 190 ierr = VecDestroy(x);CHKERRQ(ierr); 191 ierr = VecDestroy(s1);CHKERRQ(ierr); 192 ierr = VecDestroy(s2);CHKERRQ(ierr); 193 PetscFunctionReturn(0); 194 } 195 196 #undef __FUNCT__ 197 #define __FUNCT__ "MatMultTransposeAddEqual" 198 /*@ 199 MatMultTransposeAddEqual - Compares matrix-vector products of two matrices. 200 201 Collective on Mat 202 203 Input Parameters: 204 + A - the first matrix 205 - B - the second matrix 206 - n - number of random vectors to be tested 207 208 Output Parameter: 209 . flg - PETSC_TRUE if the products are equal; PETSC_FALSE otherwise. 210 211 Level: intermediate 212 213 Concepts: matrices^equality between 214 @*/ 215 PetscErrorCode MatMultTransposeAddEqual(Mat A,Mat B,PetscInt n,PetscTruth *flg) 216 { 217 PetscErrorCode ierr; 218 Vec x,y,s1,s2; 219 PetscRandom rctx; 220 PetscReal r1,r2,tol=1.e-10; 221 PetscInt am,an,bm,bn,k; 222 223 PetscFunctionBegin; 224 ierr = MatGetLocalSize(A,&am,&an);CHKERRQ(ierr); 225 ierr = MatGetLocalSize(B,&bm,&bn);CHKERRQ(ierr); 226 if (am != bm || an != bn) SETERRQ4(PETSC_ERR_ARG_SIZ,"Mat A,Mat B: local dim %D %D %D %D",am,bm,an,bn); 227 PetscCheckSameComm(A,1,B,2); 228 ierr = PetscRandomCreate(A->comm,RANDOM_DEFAULT,&rctx);CHKERRQ(ierr); 229 ierr = VecCreate(A->comm,&x);CHKERRQ(ierr); 230 ierr = VecSetSizes(x,am,PETSC_DECIDE);CHKERRQ(ierr); 231 ierr = VecSetFromOptions(x);CHKERRQ(ierr); 232 233 ierr = VecCreate(A->comm,&s1);CHKERRQ(ierr); 234 ierr = VecSetSizes(s1,an,PETSC_DECIDE);CHKERRQ(ierr); 235 ierr = VecSetFromOptions(s1);CHKERRQ(ierr); 236 ierr = VecDuplicate(s1,&s2);CHKERRQ(ierr); 237 ierr = VecDuplicate(s1,&y);CHKERRQ(ierr); 238 239 *flg = PETSC_TRUE; 240 for (k=0; k<n; k++) { 241 ierr = VecSetRandom(rctx,x);CHKERRQ(ierr); 242 ierr = VecSetRandom(rctx,y);CHKERRQ(ierr); 243 ierr = MatMultTransposeAdd(A,x,y,s1);CHKERRQ(ierr); 244 ierr = MatMultTransposeAdd(B,x,y,s2);CHKERRQ(ierr); 245 ierr = VecNorm(s1,NORM_1,&r1);CHKERRQ(ierr); 246 ierr = VecNorm(s2,NORM_1,&r2);CHKERRQ(ierr); 247 r1 -= r2; 248 if (r1<-tol || r1>tol) { 249 *flg = PETSC_FALSE; 250 ierr = PetscPrintf(PETSC_COMM_SELF,"Error: %d-th MatMultTransposeAdd() %g\n",k,r1); 251 break; 252 } 253 } 254 ierr = PetscRandomDestroy(rctx);CHKERRQ(ierr); 255 ierr = VecDestroy(x);CHKERRQ(ierr); 256 ierr = VecDestroy(y);CHKERRQ(ierr); 257 ierr = VecDestroy(s1);CHKERRQ(ierr); 258 ierr = VecDestroy(s2);CHKERRQ(ierr); 259 PetscFunctionReturn(0); 260 } 261