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