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