Actual source code: transm.c

  1: #include <../src/mat/impls/shell/shell.h>

  3: static PetscErrorCode MatMult_Transpose(Mat N, Vec x, Vec y)
  4: {
  5:   Mat A;

  7:   PetscFunctionBegin;
  8:   PetscCall(MatShellGetContext(N, &A));
  9:   PetscCall(MatMultTranspose(A, x, y));
 10:   PetscFunctionReturn(PETSC_SUCCESS);
 11: }

 13: static PetscErrorCode MatMultTranspose_Transpose(Mat N, Vec x, Vec y)
 14: {
 15:   Mat A;

 17:   PetscFunctionBegin;
 18:   PetscCall(MatShellGetContext(N, &A));
 19:   PetscCall(MatMult(A, x, y));
 20:   PetscFunctionReturn(PETSC_SUCCESS);
 21: }

 23: static PetscErrorCode MatSolve_Transpose_LU(Mat N, Vec b, Vec x)
 24: {
 25:   Mat A;

 27:   PetscFunctionBegin;
 28:   PetscCall(MatShellGetContext(N, &A));
 29:   PetscCall(MatSolveTranspose(A, b, x));
 30:   PetscFunctionReturn(PETSC_SUCCESS);
 31: }

 33: static PetscErrorCode MatSolveAdd_Transpose_LU(Mat N, Vec b, Vec y, Vec x)
 34: {
 35:   Mat A;

 37:   PetscFunctionBegin;
 38:   PetscCall(MatShellGetContext(N, &A));
 39:   PetscCall(MatSolveTransposeAdd(A, b, y, x));
 40:   PetscFunctionReturn(PETSC_SUCCESS);
 41: }

 43: static PetscErrorCode MatSolveTranspose_Transpose_LU(Mat N, Vec b, Vec x)
 44: {
 45:   Mat A;

 47:   PetscFunctionBegin;
 48:   PetscCall(MatShellGetContext(N, &A));
 49:   PetscCall(MatSolve(A, b, x));
 50:   PetscFunctionReturn(PETSC_SUCCESS);
 51: }

 53: static PetscErrorCode MatSolveTransposeAdd_Transpose_LU(Mat N, Vec b, Vec y, Vec x)
 54: {
 55:   Mat A;

 57:   PetscFunctionBegin;
 58:   PetscCall(MatShellGetContext(N, &A));
 59:   PetscCall(MatSolveAdd(A, b, y, x));
 60:   PetscFunctionReturn(PETSC_SUCCESS);
 61: }

 63: static PetscErrorCode MatMatSolve_Transpose_LU(Mat N, Mat B, Mat X)
 64: {
 65:   Mat A;

 67:   PetscFunctionBegin;
 68:   PetscCall(MatShellGetContext(N, &A));
 69:   PetscCall(MatMatSolveTranspose(A, B, X));
 70:   PetscFunctionReturn(PETSC_SUCCESS);
 71: }

 73: static PetscErrorCode MatMatSolveTranspose_Transpose_LU(Mat N, Mat B, Mat X)
 74: {
 75:   Mat A;

 77:   PetscFunctionBegin;
 78:   PetscCall(MatShellGetContext(N, &A));
 79:   PetscCall(MatMatSolve(A, B, X));
 80:   PetscFunctionReturn(PETSC_SUCCESS);
 81: }

 83: static PetscErrorCode MatLUFactor_Transpose(Mat N, IS row, IS col, const MatFactorInfo *minfo)
 84: {
 85:   Mat A;

 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: }

 99: static PetscErrorCode MatSolve_Transpose_Cholesky(Mat N, Vec b, Vec x)
100: {
101:   Mat A;

103:   PetscFunctionBegin;
104:   PetscCall(MatShellGetContext(N, &A));
105:   PetscCall(MatSolveTranspose(A, b, x));
106:   PetscFunctionReturn(PETSC_SUCCESS);
107: }

109: static PetscErrorCode MatSolveAdd_Transpose_Cholesky(Mat N, Vec b, Vec y, Vec x)
110: {
111:   Mat A;

113:   PetscFunctionBegin;
114:   PetscCall(MatShellGetContext(N, &A));
115:   PetscCall(MatSolveTransposeAdd(A, b, y, x));
116:   PetscFunctionReturn(PETSC_SUCCESS);
117: }

119: static PetscErrorCode MatSolveTranspose_Transpose_Cholesky(Mat N, Vec b, Vec x)
120: {
121:   Mat A;

123:   PetscFunctionBegin;
124:   PetscCall(MatShellGetContext(N, &A));
125:   PetscCall(MatSolve(A, b, x));
126:   PetscFunctionReturn(PETSC_SUCCESS);
127: }

129: static PetscErrorCode MatSolveTransposeAdd_Transpose_Cholesky(Mat N, Vec b, Vec y, Vec x)
130: {
131:   Mat A;

133:   PetscFunctionBegin;
134:   PetscCall(MatShellGetContext(N, &A));
135:   PetscCall(MatSolveAdd(A, b, y, x));
136:   PetscFunctionReturn(PETSC_SUCCESS);
137: }

139: static PetscErrorCode MatMatSolve_Transpose_Cholesky(Mat N, Mat B, Mat X)
140: {
141:   Mat A;

143:   PetscFunctionBegin;
144:   PetscCall(MatShellGetContext(N, &A));
145:   PetscCall(MatMatSolveTranspose(A, B, X));
146:   PetscFunctionReturn(PETSC_SUCCESS);
147: }

149: static PetscErrorCode MatMatSolveTranspose_Transpose_Cholesky(Mat N, Mat B, Mat X)
150: {
151:   Mat A;

153:   PetscFunctionBegin;
154:   PetscCall(MatShellGetContext(N, &A));
155:   PetscCall(MatMatSolve(A, B, X));
156:   PetscFunctionReturn(PETSC_SUCCESS);
157: }

159: static PetscErrorCode MatCholeskyFactor_Transpose(Mat N, IS perm, const MatFactorInfo *minfo)
160: {
161:   Mat A;

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: }

175: static PetscErrorCode MatLUFactorNumeric_Transpose(Mat F, Mat N, const MatFactorInfo *info)
176: {
177:   Mat A, FA;

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: }

192: static PetscErrorCode MatLUFactorSymbolic_Transpose(Mat F, Mat N, IS row, IS col, const MatFactorInfo *info)
193: {
194:   Mat A, FA;

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: }

204: static PetscErrorCode MatCholeskyFactorNumeric_Transpose(Mat F, Mat N, const MatFactorInfo *info)
205: {
206:   Mat A, FA;

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: }

221: static PetscErrorCode MatCholeskyFactorSymbolic_Transpose(Mat F, Mat N, IS perm, const MatFactorInfo *info)
222: {
223:   Mat A, FA;

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: }

233: static PetscErrorCode MatGetFactor_Transpose(Mat N, MatSolverType type, MatFactorType ftype, Mat *F)
234: {
235:   Mat A, FA;

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: }

251: static PetscErrorCode MatDestroy_Transpose(Mat N)
252: {
253:   Mat A;

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: }

265: static PetscErrorCode MatGetInfo_Transpose(Mat N, MatInfoType flag, MatInfo *info)
266: {
267:   Mat A;

269:   PetscFunctionBegin;
270:   PetscCall(MatShellGetContext(N, &A));
271:   PetscCall(MatGetInfo(A, flag, info));
272:   PetscFunctionReturn(PETSC_SUCCESS);
273: }

275: static PetscErrorCode MatFactorGetSolverType_Transpose(Mat N, MatSolverType *type)
276: {
277:   Mat A;

279:   PetscFunctionBegin;
280:   PetscCall(MatShellGetContext(N, &A));
281:   PetscCall(MatFactorGetSolverType(A, type));
282:   PetscFunctionReturn(PETSC_SUCCESS);
283: }

285: static PetscErrorCode MatDuplicate_Transpose(Mat N, MatDuplicateOption op, Mat *m)
286: {
287:   Mat A, C;

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: }

298: static PetscErrorCode MatHasOperation_Transpose(Mat mat, MatOperation op, PetscBool *has)
299: {
300:   Mat A;

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: }

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;

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;

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: }

415: static PetscErrorCode MatGetDiagonal_Transpose(Mat N, Vec v)
416: {
417:   Mat A;

419:   PetscFunctionBegin;
420:   PetscCall(MatShellGetContext(N, &A));
421:   PetscCall(MatGetDiagonal(A, v));
422:   PetscFunctionReturn(PETSC_SUCCESS);
423: }

425: static PetscErrorCode MatCopy_Transpose(Mat A, Mat B, MatStructure str)
426: {
427:   Mat a, b;

429:   PetscFunctionBegin;
430:   PetscCall(MatShellGetContext(A, &a));
431:   PetscCall(MatShellGetContext(B, &b));
432:   PetscCall(MatCopy(a, b, str));
433:   PetscFunctionReturn(PETSC_SUCCESS);
434: }

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;

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;

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: }

470: static PetscErrorCode MatTransposeGetMat_Transpose(Mat N, Mat *M)
471: {
472:   PetscFunctionBegin;
473:   PetscCall(MatShellGetContext(N, M));
474:   PetscFunctionReturn(PETSC_SUCCESS);
475: }

477: /*@
478:   MatTransposeGetMat - Gets the `Mat` object stored inside a `MATTRANSPOSEVIRTUAL`

480:   Logically Collective

482:   Input Parameter:
483: . A - the `MATTRANSPOSEVIRTUAL` matrix

485:   Output Parameter:
486: . M - the matrix object stored inside `A`

488:   Level: intermediate

490: .seealso: [](ch_matrices), `Mat`, `MATTRANSPOSEVIRTUAL`, `MatCreateTranspose()`
491: @*/
492: PetscErrorCode MatTransposeGetMat(Mat A, Mat *M)
493: {
494:   PetscFunctionBegin;
497:   PetscAssertPointer(M, 2);
498:   PetscUseMethod(A, "MatTransposeGetMat_C", (Mat, Mat *), (A, M));
499:   PetscFunctionReturn(PETSC_SUCCESS);
500: }

502: /*MC
503:    MATTRANSPOSEVIRTUAL - "transpose" - A matrix type that represents a virtual transpose of a matrix

505:   Level: advanced

507:   Developer Notes:
508:   This is implemented on top of `MATSHELL` to get support for scaling and shifting without requiring duplicate code

510:   Users can not call `MatShellSetOperation()` operations on this class, there is some error checking for that incorrect usage

512: .seealso: [](ch_matrices), `Mat`, `MATHERMITIANTRANSPOSEVIRTUAL`, `Mat`, `MatCreateHermitianTranspose()`, `MatCreateTranspose()`,
513:           `MATNORMALHERMITIAN`, `MATNORMAL`
514: M*/

516: /*@
517:   MatCreateTranspose - Creates a new matrix `MATTRANSPOSEVIRTUAL` object that behaves like A'

519:   Collective

521:   Input Parameter:
522: . A - the (possibly rectangular) matrix

524:   Output Parameter:
525: . N - the matrix that represents A'

527:   Level: intermediate

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

534: .seealso: [](ch_matrices), `Mat`, `MATTRANSPOSEVIRTUAL`, `MatCreateNormal()`, `MatMult()`, `MatMultTranspose()`, `MatCreate()`,
535:           `MATNORMALHERMITIAN`
536: @*/
537: PetscErrorCode MatCreateTranspose(Mat A, Mat *N)
538: {
539:   VecType vtype;

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));

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));

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));

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: }