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