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