1 2 #include <petsc/private/matimpl.h> /*I "petscmat.h" I*/ 3 4 static PetscErrorCode MatMultEqual_Private(Mat A, Mat B, PetscInt n, PetscBool *flg, PetscInt t, PetscBool add) 5 { 6 Vec Ax = NULL, Bx = NULL, s1 = NULL, s2 = NULL, Ay = NULL, By = NULL; 7 PetscRandom rctx; 8 PetscReal r1, r2, tol = PETSC_SQRT_MACHINE_EPSILON; 9 PetscInt am, an, bm, bn, k; 10 PetscScalar none = -1.0; 11 #if defined(PETSC_USE_INFO) 12 const char *sops[] = {"MatMult", "MatMultAdd", "MatMultTranspose", "MatMultTransposeAdd", "MatMultHermitianTranspose", "MatMultHermitianTransposeAdd"}; 13 const char *sop; 14 #endif 15 16 PetscFunctionBegin; 17 PetscValidHeaderSpecific(A, MAT_CLASSID, 1); 18 PetscValidHeaderSpecific(B, MAT_CLASSID, 2); 19 PetscCheckSameComm(A, 1, B, 2); 20 PetscValidLogicalCollectiveInt(A, n, 3); 21 PetscValidBoolPointer(flg, 4); 22 PetscValidLogicalCollectiveInt(A, t, 5); 23 PetscValidLogicalCollectiveBool(A, add, 6); 24 PetscCall(MatGetLocalSize(A, &am, &an)); 25 PetscCall(MatGetLocalSize(B, &bm, &bn)); 26 PetscCheck(am == bm && an == bn, PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Mat A,Mat B: local dim %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT, am, bm, an, bn); 27 #if defined(PETSC_USE_INFO) 28 sop = sops[(add ? 1 : 0) + 2 * t]; /* t = 0 => no transpose, t = 1 => transpose, t = 2 => Hermitian transpose */ 29 #endif 30 PetscCall(PetscRandomCreate(PetscObjectComm((PetscObject)A), &rctx)); 31 PetscCall(PetscRandomSetFromOptions(rctx)); 32 if (t) { 33 PetscCall(MatCreateVecs(A, &s1, &Ax)); 34 PetscCall(MatCreateVecs(B, &s2, &Bx)); 35 } else { 36 PetscCall(MatCreateVecs(A, &Ax, &s1)); 37 PetscCall(MatCreateVecs(B, &Bx, &s2)); 38 } 39 if (add) { 40 PetscCall(VecDuplicate(s1, &Ay)); 41 PetscCall(VecDuplicate(s2, &By)); 42 } 43 44 *flg = PETSC_TRUE; 45 for (k = 0; k < n; k++) { 46 PetscCall(VecSetRandom(Ax, rctx)); 47 PetscCall(VecCopy(Ax, Bx)); 48 if (add) { 49 PetscCall(VecSetRandom(Ay, rctx)); 50 PetscCall(VecCopy(Ay, By)); 51 } 52 if (t == 1) { 53 if (add) { 54 PetscCall(MatMultTransposeAdd(A, Ax, Ay, s1)); 55 PetscCall(MatMultTransposeAdd(B, Bx, By, s2)); 56 } else { 57 PetscCall(MatMultTranspose(A, Ax, s1)); 58 PetscCall(MatMultTranspose(B, Bx, s2)); 59 } 60 } else if (t == 2) { 61 if (add) { 62 PetscCall(MatMultHermitianTransposeAdd(A, Ax, Ay, s1)); 63 PetscCall(MatMultHermitianTransposeAdd(B, Bx, By, s2)); 64 } else { 65 PetscCall(MatMultHermitianTranspose(A, Ax, s1)); 66 PetscCall(MatMultHermitianTranspose(B, Bx, s2)); 67 } 68 } else { 69 if (add) { 70 PetscCall(MatMultAdd(A, Ax, Ay, s1)); 71 PetscCall(MatMultAdd(B, Bx, By, s2)); 72 } else { 73 PetscCall(MatMult(A, Ax, s1)); 74 PetscCall(MatMult(B, Bx, s2)); 75 } 76 } 77 PetscCall(VecNorm(s2, NORM_INFINITY, &r2)); 78 if (r2 < tol) { 79 PetscCall(VecNorm(s1, NORM_INFINITY, &r1)); 80 } else { 81 PetscCall(VecAXPY(s2, none, s1)); 82 PetscCall(VecNorm(s2, NORM_INFINITY, &r1)); 83 r1 /= r2; 84 } 85 if (r1 > tol) { 86 *flg = PETSC_FALSE; 87 PetscCall(PetscInfo(A, "Error: %" PetscInt_FMT "-th %s() %g\n", k, sop, (double)r1)); 88 break; 89 } 90 } 91 PetscCall(PetscRandomDestroy(&rctx)); 92 PetscCall(VecDestroy(&Ax)); 93 PetscCall(VecDestroy(&Bx)); 94 PetscCall(VecDestroy(&Ay)); 95 PetscCall(VecDestroy(&By)); 96 PetscCall(VecDestroy(&s1)); 97 PetscCall(VecDestroy(&s2)); 98 PetscFunctionReturn(PETSC_SUCCESS); 99 } 100 101 static PetscErrorCode MatMatMultEqual_Private(Mat A, Mat B, Mat C, PetscInt n, PetscBool *flg, PetscBool At, PetscBool Bt) 102 { 103 Vec Ax, Bx, Cx, s1, s2, s3; 104 PetscRandom rctx; 105 PetscReal r1, r2, tol = PETSC_SQRT_MACHINE_EPSILON; 106 PetscInt am, an, bm, bn, cm, cn, k; 107 PetscScalar none = -1.0; 108 #if defined(PETSC_USE_INFO) 109 const char *sops[] = {"MatMatMult", "MatTransposeMatMult", "MatMatTransposeMult", "MatTransposeMatTransposeMult"}; 110 const char *sop; 111 #endif 112 113 PetscFunctionBegin; 114 PetscValidHeaderSpecific(A, MAT_CLASSID, 1); 115 PetscValidHeaderSpecific(B, MAT_CLASSID, 2); 116 PetscCheckSameComm(A, 1, B, 2); 117 PetscValidHeaderSpecific(C, MAT_CLASSID, 3); 118 PetscCheckSameComm(A, 1, C, 3); 119 PetscValidLogicalCollectiveInt(A, n, 4); 120 PetscValidBoolPointer(flg, 5); 121 PetscValidLogicalCollectiveBool(A, At, 6); 122 PetscValidLogicalCollectiveBool(B, Bt, 7); 123 PetscCall(MatGetLocalSize(A, &am, &an)); 124 PetscCall(MatGetLocalSize(B, &bm, &bn)); 125 PetscCall(MatGetLocalSize(C, &cm, &cn)); 126 if (At) { 127 PetscInt tt = an; 128 an = am; 129 am = tt; 130 }; 131 if (Bt) { 132 PetscInt tt = bn; 133 bn = bm; 134 bm = tt; 135 }; 136 PetscCheck(an == bm && am == cm && bn == cn, PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Mat A, B, C local dim %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT, am, an, bm, bn, cm, cn); 137 138 #if defined(PETSC_USE_INFO) 139 sop = sops[(At ? 1 : 0) + 2 * (Bt ? 1 : 0)]; 140 #endif 141 PetscCall(PetscRandomCreate(PetscObjectComm((PetscObject)C), &rctx)); 142 PetscCall(PetscRandomSetFromOptions(rctx)); 143 if (Bt) { 144 PetscCall(MatCreateVecs(B, &s1, &Bx)); 145 } else { 146 PetscCall(MatCreateVecs(B, &Bx, &s1)); 147 } 148 if (At) { 149 PetscCall(MatCreateVecs(A, &s2, &Ax)); 150 } else { 151 PetscCall(MatCreateVecs(A, &Ax, &s2)); 152 } 153 PetscCall(MatCreateVecs(C, &Cx, &s3)); 154 155 *flg = PETSC_TRUE; 156 for (k = 0; k < n; k++) { 157 PetscCall(VecSetRandom(Bx, rctx)); 158 if (Bt) { 159 PetscCall(MatMultTranspose(B, Bx, s1)); 160 } else { 161 PetscCall(MatMult(B, Bx, s1)); 162 } 163 PetscCall(VecCopy(s1, Ax)); 164 if (At) { 165 PetscCall(MatMultTranspose(A, Ax, s2)); 166 } else { 167 PetscCall(MatMult(A, Ax, s2)); 168 } 169 PetscCall(VecCopy(Bx, Cx)); 170 PetscCall(MatMult(C, Cx, s3)); 171 172 PetscCall(VecNorm(s2, NORM_INFINITY, &r2)); 173 if (r2 < tol) { 174 PetscCall(VecNorm(s3, NORM_INFINITY, &r1)); 175 } else { 176 PetscCall(VecAXPY(s2, none, s3)); 177 PetscCall(VecNorm(s2, NORM_INFINITY, &r1)); 178 r1 /= r2; 179 } 180 if (r1 > tol) { 181 *flg = PETSC_FALSE; 182 PetscCall(PetscInfo(A, "Error: %" PetscInt_FMT "-th %s %g\n", k, sop, (double)r1)); 183 break; 184 } 185 } 186 PetscCall(PetscRandomDestroy(&rctx)); 187 PetscCall(VecDestroy(&Ax)); 188 PetscCall(VecDestroy(&Bx)); 189 PetscCall(VecDestroy(&Cx)); 190 PetscCall(VecDestroy(&s1)); 191 PetscCall(VecDestroy(&s2)); 192 PetscCall(VecDestroy(&s3)); 193 PetscFunctionReturn(PETSC_SUCCESS); 194 } 195 196 /*@ 197 MatMultEqual - Compares matrix-vector products of two matrices. 198 199 Collective 200 201 Input Parameters: 202 + A - the first matrix 203 . B - the second matrix 204 - n - number of random vectors to be tested 205 206 Output Parameter: 207 . flg - PETSC_TRUE if the products are equal; PETSC_FALSE otherwise. 208 209 Level: intermediate 210 211 @*/ 212 PetscErrorCode MatMultEqual(Mat A, Mat B, PetscInt n, PetscBool *flg) 213 { 214 PetscFunctionBegin; 215 PetscCall(MatMultEqual_Private(A, B, n, flg, 0, PETSC_FALSE)); 216 PetscFunctionReturn(PETSC_SUCCESS); 217 } 218 219 /*@ 220 MatMultAddEqual - Compares matrix-vector products of two matrices. 221 222 Collective 223 224 Input Parameters: 225 + A - the first matrix 226 . B - the second matrix 227 - n - number of random vectors to be tested 228 229 Output Parameter: 230 . flg - PETSC_TRUE if the products are equal; PETSC_FALSE otherwise. 231 232 Level: intermediate 233 234 @*/ 235 PetscErrorCode MatMultAddEqual(Mat A, Mat B, PetscInt n, PetscBool *flg) 236 { 237 PetscFunctionBegin; 238 PetscCall(MatMultEqual_Private(A, B, n, flg, 0, PETSC_TRUE)); 239 PetscFunctionReturn(PETSC_SUCCESS); 240 } 241 242 /*@ 243 MatMultTransposeEqual - Compares matrix-vector products of two matrices. 244 245 Collective 246 247 Input Parameters: 248 + A - the first matrix 249 . B - the second matrix 250 - n - number of random vectors to be tested 251 252 Output Parameter: 253 . flg - PETSC_TRUE if the products are equal; PETSC_FALSE otherwise. 254 255 Level: intermediate 256 257 @*/ 258 PetscErrorCode MatMultTransposeEqual(Mat A, Mat B, PetscInt n, PetscBool *flg) 259 { 260 PetscFunctionBegin; 261 PetscCall(MatMultEqual_Private(A, B, n, flg, 1, PETSC_FALSE)); 262 PetscFunctionReturn(PETSC_SUCCESS); 263 } 264 265 /*@ 266 MatMultTransposeAddEqual - Compares matrix-vector products of two matrices. 267 268 Collective 269 270 Input Parameters: 271 + A - the first matrix 272 . B - the second matrix 273 - n - number of random vectors to be tested 274 275 Output Parameter: 276 . flg - PETSC_TRUE if the products are equal; PETSC_FALSE otherwise. 277 278 Level: intermediate 279 280 @*/ 281 PetscErrorCode MatMultTransposeAddEqual(Mat A, Mat B, PetscInt n, PetscBool *flg) 282 { 283 PetscFunctionBegin; 284 PetscCall(MatMultEqual_Private(A, B, n, flg, 1, PETSC_TRUE)); 285 PetscFunctionReturn(PETSC_SUCCESS); 286 } 287 288 /*@ 289 MatMultHermitianTransposeEqual - Compares matrix-vector products of two matrices. 290 291 Collective 292 293 Input Parameters: 294 + A - the first matrix 295 . B - the second matrix 296 - n - number of random vectors to be tested 297 298 Output Parameter: 299 . flg - PETSC_TRUE if the products are equal; PETSC_FALSE otherwise. 300 301 Level: intermediate 302 303 @*/ 304 PetscErrorCode MatMultHermitianTransposeEqual(Mat A, Mat B, PetscInt n, PetscBool *flg) 305 { 306 PetscFunctionBegin; 307 PetscCall(MatMultEqual_Private(A, B, n, flg, 2, PETSC_FALSE)); 308 PetscFunctionReturn(PETSC_SUCCESS); 309 } 310 311 /*@ 312 MatMultHermitianTransposeAddEqual - Compares matrix-vector products of two matrices. 313 314 Collective 315 316 Input Parameters: 317 + A - the first matrix 318 . B - the second matrix 319 - n - number of random vectors to be tested 320 321 Output Parameter: 322 . flg - PETSC_TRUE if the products are equal; PETSC_FALSE otherwise. 323 324 Level: intermediate 325 326 @*/ 327 PetscErrorCode MatMultHermitianTransposeAddEqual(Mat A, Mat B, PetscInt n, PetscBool *flg) 328 { 329 PetscFunctionBegin; 330 PetscCall(MatMultEqual_Private(A, B, n, flg, 2, PETSC_TRUE)); 331 PetscFunctionReturn(PETSC_SUCCESS); 332 } 333 334 /*@ 335 MatMatMultEqual - Test A*B*x = C*x for n random vector x 336 337 Collective 338 339 Input Parameters: 340 + A - the first matrix 341 . B - the second matrix 342 . C - the third matrix 343 - n - number of random vectors to be tested 344 345 Output Parameter: 346 . flg - PETSC_TRUE if the products are equal; PETSC_FALSE otherwise. 347 348 Level: intermediate 349 350 @*/ 351 PetscErrorCode MatMatMultEqual(Mat A, Mat B, Mat C, PetscInt n, PetscBool *flg) 352 { 353 PetscFunctionBegin; 354 PetscCall(MatMatMultEqual_Private(A, B, C, n, flg, PETSC_FALSE, PETSC_FALSE)); 355 PetscFunctionReturn(PETSC_SUCCESS); 356 } 357 358 /*@ 359 MatTransposeMatMultEqual - Test A^T*B*x = C*x for n random vector x 360 361 Collective 362 363 Input Parameters: 364 + A - the first matrix 365 . B - the second matrix 366 . C - the third matrix 367 - n - number of random vectors to be tested 368 369 Output Parameter: 370 . flg - PETSC_TRUE if the products are equal; PETSC_FALSE otherwise. 371 372 Level: intermediate 373 374 @*/ 375 PetscErrorCode MatTransposeMatMultEqual(Mat A, Mat B, Mat C, PetscInt n, PetscBool *flg) 376 { 377 PetscFunctionBegin; 378 PetscCall(MatMatMultEqual_Private(A, B, C, n, flg, PETSC_TRUE, PETSC_FALSE)); 379 PetscFunctionReturn(PETSC_SUCCESS); 380 } 381 382 /*@ 383 MatMatTransposeMultEqual - Test A*B^T*x = C*x for n random vector x 384 385 Collective 386 387 Input Parameters: 388 + A - the first matrix 389 . B - the second matrix 390 . C - the third matrix 391 - n - number of random vectors to be tested 392 393 Output Parameter: 394 . flg - PETSC_TRUE if the products are equal; PETSC_FALSE otherwise. 395 396 Level: intermediate 397 398 @*/ 399 PetscErrorCode MatMatTransposeMultEqual(Mat A, Mat B, Mat C, PetscInt n, PetscBool *flg) 400 { 401 PetscFunctionBegin; 402 PetscCall(MatMatMultEqual_Private(A, B, C, n, flg, PETSC_FALSE, PETSC_TRUE)); 403 PetscFunctionReturn(PETSC_SUCCESS); 404 } 405 406 static PetscErrorCode MatProjMultEqual_Private(Mat A, Mat B, Mat C, PetscInt n, PetscBool rart, PetscBool *flg) 407 { 408 Vec x, v1, v2, v3, v4, Cx, Bx; 409 PetscReal norm_abs, norm_rel, tol = PETSC_SQRT_MACHINE_EPSILON; 410 PetscInt i, am, an, bm, bn, cm, cn; 411 PetscRandom rdm; 412 PetscScalar none = -1.0; 413 414 PetscFunctionBegin; 415 PetscCall(MatGetLocalSize(A, &am, &an)); 416 PetscCall(MatGetLocalSize(B, &bm, &bn)); 417 if (rart) { 418 PetscInt t = bm; 419 bm = bn; 420 bn = t; 421 } 422 PetscCall(MatGetLocalSize(C, &cm, &cn)); 423 PetscCheck(an == bm && bn == cm && bn == cn, PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Mat A, B, C local dim %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT, am, an, bm, bn, cm, cn); 424 425 /* Create left vector of A: v2 */ 426 PetscCall(MatCreateVecs(A, &Bx, &v2)); 427 428 /* Create right vectors of B: x, v3, v4 */ 429 if (rart) { 430 PetscCall(MatCreateVecs(B, &v1, &x)); 431 } else { 432 PetscCall(MatCreateVecs(B, &x, &v1)); 433 } 434 PetscCall(VecDuplicate(x, &v3)); 435 436 PetscCall(MatCreateVecs(C, &Cx, &v4)); 437 PetscCall(PetscRandomCreate(PETSC_COMM_WORLD, &rdm)); 438 PetscCall(PetscRandomSetFromOptions(rdm)); 439 440 *flg = PETSC_TRUE; 441 for (i = 0; i < n; i++) { 442 PetscCall(VecSetRandom(x, rdm)); 443 PetscCall(VecCopy(x, Cx)); 444 PetscCall(MatMult(C, Cx, v4)); /* v4 = C*x */ 445 if (rart) { 446 PetscCall(MatMultTranspose(B, x, v1)); 447 } else { 448 PetscCall(MatMult(B, x, v1)); 449 } 450 PetscCall(VecCopy(v1, Bx)); 451 PetscCall(MatMult(A, Bx, v2)); /* v2 = A*B*x */ 452 PetscCall(VecCopy(v2, v1)); 453 if (rart) { 454 PetscCall(MatMult(B, v1, v3)); /* v3 = R*A*R^t*x */ 455 } else { 456 PetscCall(MatMultTranspose(B, v1, v3)); /* v3 = Bt*A*B*x */ 457 } 458 PetscCall(VecNorm(v4, NORM_2, &norm_abs)); 459 PetscCall(VecAXPY(v4, none, v3)); 460 PetscCall(VecNorm(v4, NORM_2, &norm_rel)); 461 462 if (norm_abs > tol) norm_rel /= norm_abs; 463 if (norm_rel > tol) { 464 *flg = PETSC_FALSE; 465 PetscCall(PetscInfo(A, "Error: %" PetscInt_FMT "-th Mat%sMult() %g\n", i, rart ? "RARt" : "PtAP", (double)norm_rel)); 466 break; 467 } 468 } 469 470 PetscCall(PetscRandomDestroy(&rdm)); 471 PetscCall(VecDestroy(&x)); 472 PetscCall(VecDestroy(&Bx)); 473 PetscCall(VecDestroy(&Cx)); 474 PetscCall(VecDestroy(&v1)); 475 PetscCall(VecDestroy(&v2)); 476 PetscCall(VecDestroy(&v3)); 477 PetscCall(VecDestroy(&v4)); 478 PetscFunctionReturn(PETSC_SUCCESS); 479 } 480 481 /*@ 482 MatPtAPMultEqual - Compares matrix-vector products of C = Bt*A*B 483 484 Collective 485 486 Input Parameters: 487 + A - the first matrix 488 . B - the second matrix 489 . C - the third matrix 490 - n - number of random vectors to be tested 491 492 Output Parameter: 493 . flg - PETSC_TRUE if the products are equal; PETSC_FALSE otherwise. 494 495 Level: intermediate 496 497 @*/ 498 PetscErrorCode MatPtAPMultEqual(Mat A, Mat B, Mat C, PetscInt n, PetscBool *flg) 499 { 500 PetscFunctionBegin; 501 PetscCall(MatProjMultEqual_Private(A, B, C, n, PETSC_FALSE, flg)); 502 PetscFunctionReturn(PETSC_SUCCESS); 503 } 504 505 /*@ 506 MatRARtMultEqual - Compares matrix-vector products of C = B*A*B^t 507 508 Collective 509 510 Input Parameters: 511 + A - the first matrix 512 . B - the second matrix 513 . C - the third matrix 514 - n - number of random vectors to be tested 515 516 Output Parameter: 517 . flg - PETSC_TRUE if the products are equal; PETSC_FALSE otherwise. 518 519 Level: intermediate 520 521 @*/ 522 PetscErrorCode MatRARtMultEqual(Mat A, Mat B, Mat C, PetscInt n, PetscBool *flg) 523 { 524 PetscFunctionBegin; 525 PetscCall(MatProjMultEqual_Private(A, B, C, n, PETSC_TRUE, flg)); 526 PetscFunctionReturn(PETSC_SUCCESS); 527 } 528 529 /*@ 530 MatIsLinear - Check if a shell matrix A is a linear operator. 531 532 Collective 533 534 Input Parameters: 535 + A - the shell matrix 536 - n - number of random vectors to be tested 537 538 Output Parameter: 539 . flg - PETSC_TRUE if the shell matrix is linear; PETSC_FALSE otherwise. 540 541 Level: intermediate 542 @*/ 543 PetscErrorCode MatIsLinear(Mat A, PetscInt n, PetscBool *flg) 544 { 545 Vec x, y, s1, s2; 546 PetscRandom rctx; 547 PetscScalar a; 548 PetscInt k; 549 PetscReal norm, normA; 550 MPI_Comm comm; 551 PetscMPIInt rank; 552 553 PetscFunctionBegin; 554 PetscValidHeaderSpecific(A, MAT_CLASSID, 1); 555 PetscCall(PetscObjectGetComm((PetscObject)A, &comm)); 556 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 557 558 PetscCall(PetscRandomCreate(comm, &rctx)); 559 PetscCall(PetscRandomSetFromOptions(rctx)); 560 PetscCall(MatCreateVecs(A, &x, &s1)); 561 PetscCall(VecDuplicate(x, &y)); 562 PetscCall(VecDuplicate(s1, &s2)); 563 564 *flg = PETSC_TRUE; 565 for (k = 0; k < n; k++) { 566 PetscCall(VecSetRandom(x, rctx)); 567 PetscCall(VecSetRandom(y, rctx)); 568 if (rank == 0) PetscCall(PetscRandomGetValue(rctx, &a)); 569 PetscCallMPI(MPI_Bcast(&a, 1, MPIU_SCALAR, 0, comm)); 570 571 /* s2 = a*A*x + A*y */ 572 PetscCall(MatMult(A, y, s2)); /* s2 = A*y */ 573 PetscCall(MatMult(A, x, s1)); /* s1 = A*x */ 574 PetscCall(VecAXPY(s2, a, s1)); /* s2 = a s1 + s2 */ 575 576 /* s1 = A * (a x + y) */ 577 PetscCall(VecAXPY(y, a, x)); /* y = a x + y */ 578 PetscCall(MatMult(A, y, s1)); 579 PetscCall(VecNorm(s1, NORM_INFINITY, &normA)); 580 581 PetscCall(VecAXPY(s2, -1.0, s1)); /* s2 = - s1 + s2 */ 582 PetscCall(VecNorm(s2, NORM_INFINITY, &norm)); 583 if (norm / normA > 100. * PETSC_MACHINE_EPSILON) { 584 *flg = PETSC_FALSE; 585 PetscCall(PetscInfo(A, "Error: %" PetscInt_FMT "-th |A*(ax+y) - (a*A*x+A*y)|/|A(ax+y)| %g > tol %g\n", k, (double)(norm / normA), (double)(100. * PETSC_MACHINE_EPSILON))); 586 break; 587 } 588 } 589 PetscCall(PetscRandomDestroy(&rctx)); 590 PetscCall(VecDestroy(&x)); 591 PetscCall(VecDestroy(&y)); 592 PetscCall(VecDestroy(&s1)); 593 PetscCall(VecDestroy(&s2)); 594 PetscFunctionReturn(PETSC_SUCCESS); 595 } 596