1 #include <../src/mat/impls/shell/shell.h> /*I "petscmat.h" I*/ 2 3 static PetscErrorCode MatMult_Transpose(Mat N, Vec x, Vec y) 4 { 5 Mat A; 6 7 PetscFunctionBegin; 8 PetscCall(MatShellGetContext(N, &A)); 9 PetscCall(MatMultTranspose(A, x, y)); 10 PetscFunctionReturn(PETSC_SUCCESS); 11 } 12 13 static PetscErrorCode MatMultTranspose_Transpose(Mat N, Vec x, Vec y) 14 { 15 Mat A; 16 17 PetscFunctionBegin; 18 PetscCall(MatShellGetContext(N, &A)); 19 PetscCall(MatMult(A, x, y)); 20 PetscFunctionReturn(PETSC_SUCCESS); 21 } 22 23 static PetscErrorCode MatSolve_Transpose_LU(Mat N, Vec b, Vec x) 24 { 25 Mat A; 26 27 PetscFunctionBegin; 28 PetscCall(MatShellGetContext(N, &A)); 29 PetscCall(MatSolveTranspose(A, b, x)); 30 PetscFunctionReturn(PETSC_SUCCESS); 31 } 32 33 static PetscErrorCode MatSolveAdd_Transpose_LU(Mat N, Vec b, Vec y, Vec x) 34 { 35 Mat A; 36 37 PetscFunctionBegin; 38 PetscCall(MatShellGetContext(N, &A)); 39 PetscCall(MatSolveTransposeAdd(A, b, y, x)); 40 PetscFunctionReturn(PETSC_SUCCESS); 41 } 42 43 static PetscErrorCode MatSolveTranspose_Transpose_LU(Mat N, Vec b, Vec x) 44 { 45 Mat A; 46 47 PetscFunctionBegin; 48 PetscCall(MatShellGetContext(N, &A)); 49 PetscCall(MatSolve(A, b, x)); 50 PetscFunctionReturn(PETSC_SUCCESS); 51 } 52 53 static PetscErrorCode MatSolveTransposeAdd_Transpose_LU(Mat N, Vec b, Vec y, Vec x) 54 { 55 Mat A; 56 57 PetscFunctionBegin; 58 PetscCall(MatShellGetContext(N, &A)); 59 PetscCall(MatSolveAdd(A, b, y, x)); 60 PetscFunctionReturn(PETSC_SUCCESS); 61 } 62 63 static PetscErrorCode MatMatSolve_Transpose_LU(Mat N, Mat B, Mat X) 64 { 65 Mat A; 66 67 PetscFunctionBegin; 68 PetscCall(MatShellGetContext(N, &A)); 69 PetscCall(MatMatSolveTranspose(A, B, X)); 70 PetscFunctionReturn(PETSC_SUCCESS); 71 } 72 73 static PetscErrorCode MatMatSolveTranspose_Transpose_LU(Mat N, Mat B, Mat X) 74 { 75 Mat A; 76 77 PetscFunctionBegin; 78 PetscCall(MatShellGetContext(N, &A)); 79 PetscCall(MatMatSolve(A, B, X)); 80 PetscFunctionReturn(PETSC_SUCCESS); 81 } 82 83 static PetscErrorCode MatLUFactor_Transpose(Mat N, IS row, IS col, const MatFactorInfo *minfo) 84 { 85 Mat A; 86 87 PetscFunctionBegin; 88 PetscCall(MatShellGetContext(N, &A)); 89 PetscCall(MatLUFactor(A, col, row, minfo)); 90 PetscCall(MatShellSetOperation(N, MATOP_SOLVE, (void (*)(void))MatSolve_Transpose_LU)); 91 PetscCall(MatShellSetOperation(N, MATOP_SOLVE_ADD, (void (*)(void))MatSolveAdd_Transpose_LU)); 92 PetscCall(MatShellSetOperation(N, MATOP_SOLVE_TRANSPOSE, (void (*)(void))MatSolveTranspose_Transpose_LU)); 93 PetscCall(MatShellSetOperation(N, MATOP_SOLVE_TRANSPOSE_ADD, (void (*)(void))MatSolveTransposeAdd_Transpose_LU)); 94 PetscCall(MatShellSetOperation(N, MATOP_MAT_SOLVE, (void (*)(void))MatMatSolve_Transpose_LU)); 95 PetscCall(MatShellSetOperation(N, MATOP_MAT_SOLVE_TRANSPOSE, (void (*)(void))MatMatSolveTranspose_Transpose_LU)); 96 PetscFunctionReturn(PETSC_SUCCESS); 97 } 98 99 static PetscErrorCode MatSolve_Transpose_Cholesky(Mat N, Vec b, Vec x) 100 { 101 Mat A; 102 103 PetscFunctionBegin; 104 PetscCall(MatShellGetContext(N, &A)); 105 PetscCall(MatSolveTranspose(A, b, x)); 106 PetscFunctionReturn(PETSC_SUCCESS); 107 } 108 109 static PetscErrorCode MatSolveAdd_Transpose_Cholesky(Mat N, Vec b, Vec y, Vec x) 110 { 111 Mat A; 112 113 PetscFunctionBegin; 114 PetscCall(MatShellGetContext(N, &A)); 115 PetscCall(MatSolveTransposeAdd(A, b, y, x)); 116 PetscFunctionReturn(PETSC_SUCCESS); 117 } 118 119 static PetscErrorCode MatSolveTranspose_Transpose_Cholesky(Mat N, Vec b, Vec x) 120 { 121 Mat A; 122 123 PetscFunctionBegin; 124 PetscCall(MatShellGetContext(N, &A)); 125 PetscCall(MatSolve(A, b, x)); 126 PetscFunctionReturn(PETSC_SUCCESS); 127 } 128 129 static PetscErrorCode MatSolveTransposeAdd_Transpose_Cholesky(Mat N, Vec b, Vec y, Vec x) 130 { 131 Mat A; 132 133 PetscFunctionBegin; 134 PetscCall(MatShellGetContext(N, &A)); 135 PetscCall(MatSolveAdd(A, b, y, x)); 136 PetscFunctionReturn(PETSC_SUCCESS); 137 } 138 139 static PetscErrorCode MatMatSolve_Transpose_Cholesky(Mat N, Mat B, Mat X) 140 { 141 Mat A; 142 143 PetscFunctionBegin; 144 PetscCall(MatShellGetContext(N, &A)); 145 PetscCall(MatMatSolveTranspose(A, B, X)); 146 PetscFunctionReturn(PETSC_SUCCESS); 147 } 148 149 static PetscErrorCode MatMatSolveTranspose_Transpose_Cholesky(Mat N, Mat B, Mat X) 150 { 151 Mat A; 152 153 PetscFunctionBegin; 154 PetscCall(MatShellGetContext(N, &A)); 155 PetscCall(MatMatSolve(A, B, X)); 156 PetscFunctionReturn(PETSC_SUCCESS); 157 } 158 159 static PetscErrorCode MatCholeskyFactor_Transpose(Mat N, IS perm, const MatFactorInfo *minfo) 160 { 161 Mat A; 162 163 PetscFunctionBegin; 164 PetscCall(MatShellGetContext(N, &A)); 165 PetscCall(MatCholeskyFactor(A, perm, minfo)); 166 PetscCall(MatShellSetOperation(N, MATOP_SOLVE, (void (*)(void))MatSolve_Transpose_Cholesky)); 167 PetscCall(MatShellSetOperation(N, MATOP_SOLVE_ADD, (void (*)(void))MatSolveAdd_Transpose_Cholesky)); 168 PetscCall(MatShellSetOperation(N, MATOP_SOLVE_TRANSPOSE, (void (*)(void))MatSolveTranspose_Transpose_Cholesky)); 169 PetscCall(MatShellSetOperation(N, MATOP_SOLVE_TRANSPOSE_ADD, (void (*)(void))MatSolveTransposeAdd_Transpose_Cholesky)); 170 PetscCall(MatShellSetOperation(N, MATOP_MAT_SOLVE, (void (*)(void))MatMatSolve_Transpose_Cholesky)); 171 PetscCall(MatShellSetOperation(N, MATOP_MAT_SOLVE_TRANSPOSE, (void (*)(void))MatMatSolveTranspose_Transpose_Cholesky)); 172 PetscFunctionReturn(PETSC_SUCCESS); 173 } 174 175 static PetscErrorCode MatLUFactorNumeric_Transpose(Mat F, Mat N, const MatFactorInfo *info) 176 { 177 Mat A, FA; 178 179 PetscFunctionBegin; 180 PetscCall(MatShellGetContext(N, &A)); 181 PetscCall(MatShellGetContext(F, &FA)); 182 PetscCall(MatLUFactorNumeric(FA, A, info)); 183 PetscCall(MatShellSetOperation(F, MATOP_SOLVE, (void (*)(void))MatSolve_Transpose_LU)); 184 PetscCall(MatShellSetOperation(F, MATOP_SOLVE_ADD, (void (*)(void))MatSolveAdd_Transpose_LU)); 185 PetscCall(MatShellSetOperation(F, MATOP_SOLVE_TRANSPOSE, (void (*)(void))MatSolveTranspose_Transpose_LU)); 186 PetscCall(MatShellSetOperation(F, MATOP_SOLVE_TRANSPOSE_ADD, (void (*)(void))MatSolveTransposeAdd_Transpose_LU)); 187 PetscCall(MatShellSetOperation(F, MATOP_MAT_SOLVE, (void (*)(void))MatMatSolve_Transpose_LU)); 188 PetscCall(MatShellSetOperation(F, MATOP_MAT_SOLVE_TRANSPOSE, (void (*)(void))MatMatSolveTranspose_Transpose_LU)); 189 PetscFunctionReturn(PETSC_SUCCESS); 190 } 191 192 static PetscErrorCode MatLUFactorSymbolic_Transpose(Mat F, Mat N, IS row, IS col, const MatFactorInfo *info) 193 { 194 Mat A, FA; 195 196 PetscFunctionBegin; 197 PetscCall(MatShellGetContext(N, &A)); 198 PetscCall(MatShellGetContext(F, &FA)); 199 PetscCall(MatLUFactorSymbolic(FA, A, row, col, info)); 200 PetscCall(MatShellSetOperation(F, MATOP_LUFACTOR_NUMERIC, (void (*)(void))MatLUFactorNumeric_Transpose)); 201 PetscFunctionReturn(PETSC_SUCCESS); 202 } 203 204 static PetscErrorCode MatCholeskyFactorNumeric_Transpose(Mat F, Mat N, const MatFactorInfo *info) 205 { 206 Mat A, FA; 207 208 PetscFunctionBegin; 209 PetscCall(MatShellGetContext(N, &A)); 210 PetscCall(MatShellGetContext(F, &FA)); 211 PetscCall(MatCholeskyFactorNumeric(FA, A, info)); 212 PetscCall(MatShellSetOperation(F, MATOP_SOLVE, (void (*)(void))MatSolve_Transpose_Cholesky)); 213 PetscCall(MatShellSetOperation(F, MATOP_SOLVE_ADD, (void (*)(void))MatSolveAdd_Transpose_Cholesky)); 214 PetscCall(MatShellSetOperation(F, MATOP_SOLVE_TRANSPOSE, (void (*)(void))MatSolveTranspose_Transpose_Cholesky)); 215 PetscCall(MatShellSetOperation(F, MATOP_SOLVE_TRANSPOSE_ADD, (void (*)(void))MatSolveTransposeAdd_Transpose_Cholesky)); 216 PetscCall(MatShellSetOperation(F, MATOP_MAT_SOLVE, (void (*)(void))MatMatSolve_Transpose_Cholesky)); 217 PetscCall(MatShellSetOperation(F, MATOP_MAT_SOLVE_TRANSPOSE, (void (*)(void))MatMatSolveTranspose_Transpose_Cholesky)); 218 PetscFunctionReturn(PETSC_SUCCESS); 219 } 220 221 static PetscErrorCode MatCholeskyFactorSymbolic_Transpose(Mat F, Mat N, IS perm, const MatFactorInfo *info) 222 { 223 Mat A, FA; 224 225 PetscFunctionBegin; 226 PetscCall(MatShellGetContext(N, &A)); 227 PetscCall(MatShellGetContext(F, &FA)); 228 PetscCall(MatCholeskyFactorSymbolic(FA, A, perm, info)); 229 PetscCall(MatShellSetOperation(F, MATOP_CHOLESKY_FACTOR_NUMERIC, (void (*)(void))MatCholeskyFactorNumeric_Transpose)); 230 PetscFunctionReturn(PETSC_SUCCESS); 231 } 232 233 static PetscErrorCode MatGetFactor_Transpose(Mat N, MatSolverType type, MatFactorType ftype, Mat *F) 234 { 235 Mat A, FA; 236 237 PetscFunctionBegin; 238 PetscCall(MatShellGetContext(N, &A)); 239 PetscCall(MatGetFactor(A, type, ftype, &FA)); 240 PetscCall(MatCreateTranspose(FA, F)); 241 if (ftype == MAT_FACTOR_LU) PetscCall(MatShellSetOperation(*F, MATOP_LUFACTOR_SYMBOLIC, (void (*)(void))MatLUFactorSymbolic_Transpose)); 242 else if (ftype == MAT_FACTOR_CHOLESKY) { 243 PetscCall(MatShellSetOperation(*F, MATOP_CHOLESKY_FACTOR_SYMBOLIC, (void (*)(void))MatCholeskyFactorSymbolic_Transpose)); 244 PetscCall(MatPropagateSymmetryOptions(A, FA)); 245 } else SETERRQ(PetscObjectComm((PetscObject)N), PETSC_ERR_SUP, "Support for factor type %s not implemented in MATTRANSPOSEVIRTUAL", MatFactorTypes[ftype]); 246 (*F)->factortype = ftype; 247 PetscCall(MatDestroy(&FA)); 248 PetscFunctionReturn(PETSC_SUCCESS); 249 } 250 251 static PetscErrorCode MatDestroy_Transpose(Mat N) 252 { 253 Mat A; 254 255 PetscFunctionBegin; 256 PetscCall(MatShellGetContext(N, &A)); 257 PetscCall(MatDestroy(&A)); 258 PetscCall(PetscObjectComposeFunction((PetscObject)N, "MatTransposeGetMat_C", NULL)); 259 PetscCall(PetscObjectComposeFunction((PetscObject)N, "MatProductSetFromOptions_anytype_C", NULL)); 260 PetscCall(PetscObjectComposeFunction((PetscObject)N, "MatShellSetContext_C", NULL)); 261 PetscCall(PetscObjectComposeFunction((PetscObject)N, "MatFactorGetSolverType_C", NULL)); 262 PetscFunctionReturn(PETSC_SUCCESS); 263 } 264 265 static PetscErrorCode MatGetInfo_Transpose(Mat N, MatInfoType flag, MatInfo *info) 266 { 267 Mat A; 268 269 PetscFunctionBegin; 270 PetscCall(MatShellGetContext(N, &A)); 271 PetscCall(MatGetInfo(A, flag, info)); 272 PetscFunctionReturn(PETSC_SUCCESS); 273 } 274 275 static PetscErrorCode MatFactorGetSolverType_Transpose(Mat N, MatSolverType *type) 276 { 277 Mat A; 278 279 PetscFunctionBegin; 280 PetscCall(MatShellGetContext(N, &A)); 281 PetscCall(MatFactorGetSolverType(A, type)); 282 PetscFunctionReturn(PETSC_SUCCESS); 283 } 284 285 static PetscErrorCode MatDuplicate_Transpose(Mat N, MatDuplicateOption op, Mat *m) 286 { 287 Mat A, C; 288 289 PetscFunctionBegin; 290 PetscCall(MatShellGetContext(N, &A)); 291 PetscCall(MatDuplicate(A, op, &C)); 292 PetscCall(MatCreateTranspose(C, m)); 293 if (op == MAT_COPY_VALUES) PetscCall(MatCopy(N, *m, SAME_NONZERO_PATTERN)); 294 PetscCall(MatDestroy(&C)); 295 PetscFunctionReturn(PETSC_SUCCESS); 296 } 297 298 static PetscErrorCode MatHasOperation_Transpose(Mat mat, MatOperation op, PetscBool *has) 299 { 300 Mat A; 301 302 PetscFunctionBegin; 303 PetscCall(MatShellGetContext(mat, &A)); 304 *has = PETSC_FALSE; 305 if (op == MATOP_MULT || op == MATOP_MULT_ADD) { 306 PetscCall(MatHasOperation(A, MATOP_MULT_TRANSPOSE, has)); 307 } else if (op == MATOP_MULT_TRANSPOSE || op == MATOP_MULT_TRANSPOSE_ADD) { 308 PetscCall(MatHasOperation(A, MATOP_MULT, has)); 309 } else if (((void **)mat->ops)[op]) *has = PETSC_TRUE; 310 PetscFunctionReturn(PETSC_SUCCESS); 311 } 312 313 static PetscErrorCode MatProductSetFromOptions_Transpose(Mat D) 314 { 315 Mat A, B, C, Ain, Bin, Cin; 316 PetscBool Aistrans, Bistrans, Cistrans; 317 PetscInt Atrans, Btrans, Ctrans; 318 MatProductType ptype; 319 320 PetscFunctionBegin; 321 MatCheckProduct(D, 1); 322 A = D->product->A; 323 B = D->product->B; 324 C = D->product->C; 325 PetscCall(PetscObjectTypeCompare((PetscObject)A, MATTRANSPOSEVIRTUAL, &Aistrans)); 326 PetscCall(PetscObjectTypeCompare((PetscObject)B, MATTRANSPOSEVIRTUAL, &Bistrans)); 327 PetscCall(PetscObjectTypeCompare((PetscObject)C, MATTRANSPOSEVIRTUAL, &Cistrans)); 328 PetscCheck(Aistrans || Bistrans || Cistrans, PetscObjectComm((PetscObject)D), PETSC_ERR_PLIB, "This should not happen"); 329 Atrans = 0; 330 Ain = A; 331 while (Aistrans) { 332 Atrans++; 333 PetscCall(MatTransposeGetMat(Ain, &Ain)); 334 PetscCall(PetscObjectTypeCompare((PetscObject)Ain, MATTRANSPOSEVIRTUAL, &Aistrans)); 335 } 336 Btrans = 0; 337 Bin = B; 338 while (Bistrans) { 339 Btrans++; 340 PetscCall(MatTransposeGetMat(Bin, &Bin)); 341 PetscCall(PetscObjectTypeCompare((PetscObject)Bin, MATTRANSPOSEVIRTUAL, &Bistrans)); 342 } 343 Ctrans = 0; 344 Cin = C; 345 while (Cistrans) { 346 Ctrans++; 347 PetscCall(MatTransposeGetMat(Cin, &Cin)); 348 PetscCall(PetscObjectTypeCompare((PetscObject)Cin, MATTRANSPOSEVIRTUAL, &Cistrans)); 349 } 350 Atrans = Atrans % 2; 351 Btrans = Btrans % 2; 352 Ctrans = Ctrans % 2; 353 ptype = D->product->type; /* same product type by default */ 354 if (Ain->symmetric == PETSC_BOOL3_TRUE) Atrans = 0; 355 if (Bin->symmetric == PETSC_BOOL3_TRUE) Btrans = 0; 356 if (Cin && Cin->symmetric == PETSC_BOOL3_TRUE) Ctrans = 0; 357 358 if (Atrans || Btrans || Ctrans) { 359 ptype = MATPRODUCT_UNSPECIFIED; 360 switch (D->product->type) { 361 case MATPRODUCT_AB: 362 if (Atrans && Btrans) { /* At * Bt we do not have support for this */ 363 /* TODO custom implementation ? */ 364 } else if (Atrans) { /* At * B */ 365 ptype = MATPRODUCT_AtB; 366 } else { /* A * Bt */ 367 ptype = MATPRODUCT_ABt; 368 } 369 break; 370 case MATPRODUCT_AtB: 371 if (Atrans && Btrans) { /* A * Bt */ 372 ptype = MATPRODUCT_ABt; 373 } else if (Atrans) { /* A * B */ 374 ptype = MATPRODUCT_AB; 375 } else { /* At * Bt we do not have support for this */ 376 /* TODO custom implementation ? */ 377 } 378 break; 379 case MATPRODUCT_ABt: 380 if (Atrans && Btrans) { /* At * B */ 381 ptype = MATPRODUCT_AtB; 382 } else if (Atrans) { /* At * Bt we do not have support for this */ 383 /* TODO custom implementation ? */ 384 } else { /* A * B */ 385 ptype = MATPRODUCT_AB; 386 } 387 break; 388 case MATPRODUCT_PtAP: 389 if (Atrans) { /* PtAtP */ 390 /* TODO custom implementation ? */ 391 } else { /* RARt */ 392 ptype = MATPRODUCT_RARt; 393 } 394 break; 395 case MATPRODUCT_RARt: 396 if (Atrans) { /* RAtRt */ 397 /* TODO custom implementation ? */ 398 } else { /* PtAP */ 399 ptype = MATPRODUCT_PtAP; 400 } 401 break; 402 case MATPRODUCT_ABC: 403 /* TODO custom implementation ? */ 404 break; 405 default: 406 SETERRQ(PetscObjectComm((PetscObject)D), PETSC_ERR_SUP, "ProductType %s is not supported", MatProductTypes[D->product->type]); 407 } 408 } 409 PetscCall(MatProductReplaceMats(Ain, Bin, Cin, D)); 410 PetscCall(MatProductSetType(D, ptype)); 411 PetscCall(MatProductSetFromOptions(D)); 412 PetscFunctionReturn(PETSC_SUCCESS); 413 } 414 415 static PetscErrorCode MatGetDiagonal_Transpose(Mat N, Vec v) 416 { 417 Mat A; 418 419 PetscFunctionBegin; 420 PetscCall(MatShellGetContext(N, &A)); 421 PetscCall(MatGetDiagonal(A, v)); 422 PetscFunctionReturn(PETSC_SUCCESS); 423 } 424 425 static PetscErrorCode MatCopy_Transpose(Mat A, Mat B, MatStructure str) 426 { 427 Mat a, b; 428 429 PetscFunctionBegin; 430 PetscCall(MatShellGetContext(A, &a)); 431 PetscCall(MatShellGetContext(B, &b)); 432 PetscCall(MatCopy(a, b, str)); 433 PetscFunctionReturn(PETSC_SUCCESS); 434 } 435 436 static PetscErrorCode MatConvert_Transpose(Mat N, MatType newtype, MatReuse reuse, Mat *newmat) 437 { 438 Mat A; 439 PetscScalar vscale = 1.0, vshift = 0.0; 440 PetscBool flg; 441 442 PetscFunctionBegin; 443 PetscCall(MatShellGetContext(N, &A)); 444 PetscCall(MatHasOperation(A, MATOP_TRANSPOSE, &flg)); 445 if (flg || N->ops->getrow) { /* if this condition is false, MatConvert_Shell() will be called in MatConvert_Basic(), so the following checks are not needed */ 446 PetscCall(MatShellGetScalingShifts(N, &vshift, &vscale, (Vec *)MAT_SHELL_NOT_ALLOWED, (Vec *)MAT_SHELL_NOT_ALLOWED, (Vec *)MAT_SHELL_NOT_ALLOWED, (Mat *)MAT_SHELL_NOT_ALLOWED, (IS *)MAT_SHELL_NOT_ALLOWED, (IS *)MAT_SHELL_NOT_ALLOWED)); 447 } 448 if (flg) { 449 Mat B; 450 451 PetscCall(MatTranspose(A, MAT_INITIAL_MATRIX, &B)); 452 if (reuse != MAT_INPLACE_MATRIX) { 453 PetscCall(MatConvert(B, newtype, reuse, newmat)); 454 PetscCall(MatDestroy(&B)); 455 } else { 456 PetscCall(MatConvert(B, newtype, MAT_INPLACE_MATRIX, &B)); 457 PetscCall(MatHeaderReplace(N, &B)); 458 } 459 } else { /* use basic converter as fallback */ 460 flg = (PetscBool)(N->ops->getrow != NULL); 461 PetscCall(MatConvert_Basic(N, newtype, reuse, newmat)); 462 } 463 if (flg) { 464 PetscCall(MatScale(*newmat, vscale)); 465 PetscCall(MatShift(*newmat, vshift)); 466 } 467 PetscFunctionReturn(PETSC_SUCCESS); 468 } 469 470 static PetscErrorCode MatTransposeGetMat_Transpose(Mat N, Mat *M) 471 { 472 PetscFunctionBegin; 473 PetscCall(MatShellGetContext(N, M)); 474 PetscFunctionReturn(PETSC_SUCCESS); 475 } 476 477 /*@ 478 MatTransposeGetMat - Gets the `Mat` object stored inside a `MATTRANSPOSEVIRTUAL` 479 480 Logically Collective 481 482 Input Parameter: 483 . A - the `MATTRANSPOSEVIRTUAL` matrix 484 485 Output Parameter: 486 . M - the matrix object stored inside `A` 487 488 Level: intermediate 489 490 .seealso: [](ch_matrices), `Mat`, `MATTRANSPOSEVIRTUAL`, `MatCreateTranspose()` 491 @*/ 492 PetscErrorCode MatTransposeGetMat(Mat A, Mat *M) 493 { 494 PetscFunctionBegin; 495 PetscValidHeaderSpecific(A, MAT_CLASSID, 1); 496 PetscValidType(A, 1); 497 PetscAssertPointer(M, 2); 498 PetscUseMethod(A, "MatTransposeGetMat_C", (Mat, Mat *), (A, M)); 499 PetscFunctionReturn(PETSC_SUCCESS); 500 } 501 502 /*MC 503 MATTRANSPOSEVIRTUAL - "transpose" - A matrix type that represents a virtual transpose of a matrix 504 505 Level: advanced 506 507 Developer Notes: 508 This is implemented on top of `MATSHELL` to get support for scaling and shifting without requiring duplicate code 509 510 Users can not call `MatShellSetOperation()` operations on this class, there is some error checking for that incorrect usage 511 512 .seealso: [](ch_matrices), `Mat`, `MATHERMITIANTRANSPOSEVIRTUAL`, `Mat`, `MatCreateHermitianTranspose()`, `MatCreateTranspose()`, 513 `MATNORMALHERMITIAN`, `MATNORMAL` 514 M*/ 515 516 /*@ 517 MatCreateTranspose - Creates a new matrix `MATTRANSPOSEVIRTUAL` object that behaves like A' 518 519 Collective 520 521 Input Parameter: 522 . A - the (possibly rectangular) matrix 523 524 Output Parameter: 525 . N - the matrix that represents A' 526 527 Level: intermediate 528 529 Note: 530 The transpose A' is NOT actually formed! Rather the new matrix 531 object performs the matrix-vector product by using the `MatMultTranspose()` on 532 the original matrix 533 534 .seealso: [](ch_matrices), `Mat`, `MATTRANSPOSEVIRTUAL`, `MatCreateNormal()`, `MatMult()`, `MatMultTranspose()`, `MatCreate()`, 535 `MATNORMALHERMITIAN` 536 @*/ 537 PetscErrorCode MatCreateTranspose(Mat A, Mat *N) 538 { 539 VecType vtype; 540 541 PetscFunctionBegin; 542 PetscCall(MatCreate(PetscObjectComm((PetscObject)A), N)); 543 PetscCall(PetscLayoutReference(A->rmap, &((*N)->cmap))); 544 PetscCall(PetscLayoutReference(A->cmap, &((*N)->rmap))); 545 PetscCall(MatSetType(*N, MATSHELL)); 546 PetscCall(MatShellSetContext(*N, A)); 547 PetscCall(PetscObjectReference((PetscObject)A)); 548 549 PetscCall(MatSetBlockSizes(*N, PetscAbs(A->cmap->bs), PetscAbs(A->rmap->bs))); 550 PetscCall(MatGetVecType(A, &vtype)); 551 PetscCall(MatSetVecType(*N, vtype)); 552 #if defined(PETSC_HAVE_DEVICE) 553 PetscCall(MatBindToCPU(*N, A->boundtocpu)); 554 #endif 555 PetscCall(MatSetUp(*N)); 556 557 PetscCall(MatShellSetOperation(*N, MATOP_DESTROY, (void (*)(void))MatDestroy_Transpose)); 558 PetscCall(MatShellSetOperation(*N, MATOP_MULT, (void (*)(void))MatMult_Transpose)); 559 PetscCall(MatShellSetOperation(*N, MATOP_MULT_TRANSPOSE, (void (*)(void))MatMultTranspose_Transpose)); 560 PetscCall(MatShellSetOperation(*N, MATOP_LUFACTOR, (void (*)(void))MatLUFactor_Transpose)); 561 PetscCall(MatShellSetOperation(*N, MATOP_CHOLESKYFACTOR, (void (*)(void))MatCholeskyFactor_Transpose)); 562 PetscCall(MatShellSetOperation(*N, MATOP_GET_FACTOR, (void (*)(void))MatGetFactor_Transpose)); 563 PetscCall(MatShellSetOperation(*N, MATOP_GETINFO, (void (*)(void))MatGetInfo_Transpose)); 564 PetscCall(MatShellSetOperation(*N, MATOP_DUPLICATE, (void (*)(void))MatDuplicate_Transpose)); 565 PetscCall(MatShellSetOperation(*N, MATOP_HAS_OPERATION, (void (*)(void))MatHasOperation_Transpose)); 566 PetscCall(MatShellSetOperation(*N, MATOP_GET_DIAGONAL, (void (*)(void))MatGetDiagonal_Transpose)); 567 PetscCall(MatShellSetOperation(*N, MATOP_COPY, (void (*)(void))MatCopy_Transpose)); 568 PetscCall(MatShellSetOperation(*N, MATOP_CONVERT, (void (*)(void))MatConvert_Transpose)); 569 570 PetscCall(PetscObjectComposeFunction((PetscObject)*N, "MatTransposeGetMat_C", MatTransposeGetMat_Transpose)); 571 PetscCall(PetscObjectComposeFunction((PetscObject)*N, "MatProductSetFromOptions_anytype_C", MatProductSetFromOptions_Transpose)); 572 PetscCall(PetscObjectComposeFunction((PetscObject)*N, "MatFactorGetSolverType_C", MatFactorGetSolverType_Transpose)); 573 PetscCall(PetscObjectComposeFunction((PetscObject)*N, "MatShellSetContext_C", MatShellSetContext_Immutable)); 574 PetscCall(PetscObjectComposeFunction((PetscObject)*N, "MatShellSetContextDestroy_C", MatShellSetContextDestroy_Immutable)); 575 PetscCall(PetscObjectComposeFunction((PetscObject)*N, "MatShellSetManageScalingShifts_C", MatShellSetManageScalingShifts_Immutable)); 576 PetscCall(PetscObjectChangeTypeName((PetscObject)*N, MATTRANSPOSEVIRTUAL)); 577 PetscFunctionReturn(PETSC_SUCCESS); 578 } 579