1 2 #include <petsc/private/matimpl.h> /*I "petscmat.h" I*/ 3 4 const char *const MatCompositeMergeTypes[] = {"left", "right", "MatCompositeMergeType", "MAT_COMPOSITE_", NULL}; 5 6 typedef struct _Mat_CompositeLink *Mat_CompositeLink; 7 struct _Mat_CompositeLink { 8 Mat mat; 9 Vec work; 10 Mat_CompositeLink next, prev; 11 }; 12 13 typedef struct { 14 MatCompositeType type; 15 Mat_CompositeLink head, tail; 16 Vec work; 17 PetscScalar scale; /* scale factor supplied with MatScale() */ 18 Vec left, right; /* left and right diagonal scaling provided with MatDiagonalScale() */ 19 Vec leftwork, rightwork, leftwork2, rightwork2; /* Two pairs of working vectors */ 20 PetscInt nmat; 21 PetscBool merge; 22 MatCompositeMergeType mergetype; 23 MatStructure structure; 24 25 PetscScalar *scalings; 26 PetscBool merge_mvctx; /* Whether need to merge mvctx of component matrices */ 27 Vec *lvecs; /* [nmat] Basically, they are Mvctx->lvec of each component matrix */ 28 PetscScalar *larray; /* [len] Data arrays of lvecs[] are stored consecutively in larray */ 29 PetscInt len; /* Length of larray[] */ 30 Vec gvec; /* Union of lvecs[] without duplicated entries */ 31 PetscInt *location; /* A map that maps entries in garray[] to larray[] */ 32 VecScatter Mvctx; 33 } Mat_Composite; 34 35 PetscErrorCode MatDestroy_Composite(Mat mat) 36 { 37 Mat_Composite *shell = (Mat_Composite *)mat->data; 38 Mat_CompositeLink next = shell->head, oldnext; 39 PetscInt i; 40 41 PetscFunctionBegin; 42 while (next) { 43 PetscCall(MatDestroy(&next->mat)); 44 if (next->work && (!next->next || next->work != next->next->work)) PetscCall(VecDestroy(&next->work)); 45 oldnext = next; 46 next = next->next; 47 PetscCall(PetscFree(oldnext)); 48 } 49 PetscCall(VecDestroy(&shell->work)); 50 PetscCall(VecDestroy(&shell->left)); 51 PetscCall(VecDestroy(&shell->right)); 52 PetscCall(VecDestroy(&shell->leftwork)); 53 PetscCall(VecDestroy(&shell->rightwork)); 54 PetscCall(VecDestroy(&shell->leftwork2)); 55 PetscCall(VecDestroy(&shell->rightwork2)); 56 57 if (shell->Mvctx) { 58 for (i = 0; i < shell->nmat; i++) PetscCall(VecDestroy(&shell->lvecs[i])); 59 PetscCall(PetscFree3(shell->location, shell->larray, shell->lvecs)); 60 PetscCall(PetscFree(shell->larray)); 61 PetscCall(VecDestroy(&shell->gvec)); 62 PetscCall(VecScatterDestroy(&shell->Mvctx)); 63 } 64 65 PetscCall(PetscFree(shell->scalings)); 66 PetscCall(PetscObjectComposeFunction((PetscObject)mat, "MatCompositeAddMat_C", NULL)); 67 PetscCall(PetscObjectComposeFunction((PetscObject)mat, "MatCompositeSetType_C", NULL)); 68 PetscCall(PetscObjectComposeFunction((PetscObject)mat, "MatCompositeGetType_C", NULL)); 69 PetscCall(PetscObjectComposeFunction((PetscObject)mat, "MatCompositeSetMergeType_C", NULL)); 70 PetscCall(PetscObjectComposeFunction((PetscObject)mat, "MatCompositeSetMatStructure_C", NULL)); 71 PetscCall(PetscObjectComposeFunction((PetscObject)mat, "MatCompositeGetMatStructure_C", NULL)); 72 PetscCall(PetscObjectComposeFunction((PetscObject)mat, "MatCompositeMerge_C", NULL)); 73 PetscCall(PetscObjectComposeFunction((PetscObject)mat, "MatCompositeGetNumberMat_C", NULL)); 74 PetscCall(PetscObjectComposeFunction((PetscObject)mat, "MatCompositeGetMat_C", NULL)); 75 PetscCall(PetscObjectComposeFunction((PetscObject)mat, "MatCompositeSetScalings_C", NULL)); 76 PetscCall(PetscFree(mat->data)); 77 PetscFunctionReturn(PETSC_SUCCESS); 78 } 79 80 PetscErrorCode MatMult_Composite_Multiplicative(Mat A, Vec x, Vec y) 81 { 82 Mat_Composite *shell = (Mat_Composite *)A->data; 83 Mat_CompositeLink next = shell->head; 84 Vec in, out; 85 PetscScalar scale; 86 PetscInt i; 87 88 PetscFunctionBegin; 89 PetscCheck(next, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Must provide at least one matrix with MatCompositeAddMat()"); 90 in = x; 91 if (shell->right) { 92 if (!shell->rightwork) PetscCall(VecDuplicate(shell->right, &shell->rightwork)); 93 PetscCall(VecPointwiseMult(shell->rightwork, shell->right, in)); 94 in = shell->rightwork; 95 } 96 while (next->next) { 97 if (!next->work) { /* should reuse previous work if the same size */ 98 PetscCall(MatCreateVecs(next->mat, NULL, &next->work)); 99 } 100 out = next->work; 101 PetscCall(MatMult(next->mat, in, out)); 102 in = out; 103 next = next->next; 104 } 105 PetscCall(MatMult(next->mat, in, y)); 106 if (shell->left) PetscCall(VecPointwiseMult(y, shell->left, y)); 107 scale = shell->scale; 108 if (shell->scalings) { 109 for (i = 0; i < shell->nmat; i++) scale *= shell->scalings[i]; 110 } 111 PetscCall(VecScale(y, scale)); 112 PetscFunctionReturn(PETSC_SUCCESS); 113 } 114 115 PetscErrorCode MatMultTranspose_Composite_Multiplicative(Mat A, Vec x, Vec y) 116 { 117 Mat_Composite *shell = (Mat_Composite *)A->data; 118 Mat_CompositeLink tail = shell->tail; 119 Vec in, out; 120 PetscScalar scale; 121 PetscInt i; 122 123 PetscFunctionBegin; 124 PetscCheck(tail, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Must provide at least one matrix with MatCompositeAddMat()"); 125 in = x; 126 if (shell->left) { 127 if (!shell->leftwork) PetscCall(VecDuplicate(shell->left, &shell->leftwork)); 128 PetscCall(VecPointwiseMult(shell->leftwork, shell->left, in)); 129 in = shell->leftwork; 130 } 131 while (tail->prev) { 132 if (!tail->prev->work) { /* should reuse previous work if the same size */ 133 PetscCall(MatCreateVecs(tail->mat, NULL, &tail->prev->work)); 134 } 135 out = tail->prev->work; 136 PetscCall(MatMultTranspose(tail->mat, in, out)); 137 in = out; 138 tail = tail->prev; 139 } 140 PetscCall(MatMultTranspose(tail->mat, in, y)); 141 if (shell->right) PetscCall(VecPointwiseMult(y, shell->right, y)); 142 143 scale = shell->scale; 144 if (shell->scalings) { 145 for (i = 0; i < shell->nmat; i++) scale *= shell->scalings[i]; 146 } 147 PetscCall(VecScale(y, scale)); 148 PetscFunctionReturn(PETSC_SUCCESS); 149 } 150 151 PetscErrorCode MatMult_Composite(Mat mat, Vec x, Vec y) 152 { 153 Mat_Composite *shell = (Mat_Composite *)mat->data; 154 Mat_CompositeLink cur = shell->head; 155 Vec in, y2, xin; 156 Mat A, B; 157 PetscInt i, j, k, n, nuniq, lo, hi, mid, *gindices, *buf, *tmp, tot; 158 const PetscScalar *vals; 159 const PetscInt *garray; 160 IS ix, iy; 161 PetscBool match; 162 163 PetscFunctionBegin; 164 PetscCheck(cur, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Must provide at least one matrix with MatCompositeAddMat()"); 165 in = x; 166 if (shell->right) { 167 if (!shell->rightwork) PetscCall(VecDuplicate(shell->right, &shell->rightwork)); 168 PetscCall(VecPointwiseMult(shell->rightwork, shell->right, in)); 169 in = shell->rightwork; 170 } 171 172 /* Try to merge Mvctx when instructed but not yet done. We did not do it in MatAssemblyEnd() since at that time 173 we did not know whether mat is ADDITIVE or MULTIPLICATIVE. Only now we are assured mat is ADDITIVE and 174 it is legal to merge Mvctx, because all component matrices have the same size. 175 */ 176 if (shell->merge_mvctx && !shell->Mvctx) { 177 /* Currently only implemented for MATMPIAIJ */ 178 for (cur = shell->head; cur; cur = cur->next) { 179 PetscCall(PetscObjectTypeCompare((PetscObject)cur->mat, MATMPIAIJ, &match)); 180 if (!match) { 181 shell->merge_mvctx = PETSC_FALSE; 182 goto skip_merge_mvctx; 183 } 184 } 185 186 /* Go through matrices first time to count total number of nonzero off-diag columns (may have dups) */ 187 tot = 0; 188 for (cur = shell->head; cur; cur = cur->next) { 189 PetscCall(MatMPIAIJGetSeqAIJ(cur->mat, NULL, &B, NULL)); 190 PetscCall(MatGetLocalSize(B, NULL, &n)); 191 tot += n; 192 } 193 PetscCall(PetscMalloc3(tot, &shell->location, tot, &shell->larray, shell->nmat, &shell->lvecs)); 194 shell->len = tot; 195 196 /* Go through matrices second time to sort off-diag columns and remove dups */ 197 PetscCall(PetscMalloc1(tot, &gindices)); /* No Malloc2() since we will give one to petsc and free the other */ 198 PetscCall(PetscMalloc1(tot, &buf)); 199 nuniq = 0; /* Number of unique nonzero columns */ 200 for (cur = shell->head; cur; cur = cur->next) { 201 PetscCall(MatMPIAIJGetSeqAIJ(cur->mat, NULL, &B, &garray)); 202 PetscCall(MatGetLocalSize(B, NULL, &n)); 203 /* Merge pre-sorted garray[0,n) and gindices[0,nuniq) to buf[] */ 204 i = j = k = 0; 205 while (i < n && j < nuniq) { 206 if (garray[i] < gindices[j]) buf[k++] = garray[i++]; 207 else if (garray[i] > gindices[j]) buf[k++] = gindices[j++]; 208 else { 209 buf[k++] = garray[i++]; 210 j++; 211 } 212 } 213 /* Copy leftover in garray[] or gindices[] */ 214 if (i < n) { 215 PetscCall(PetscArraycpy(buf + k, garray + i, n - i)); 216 nuniq = k + n - i; 217 } else if (j < nuniq) { 218 PetscCall(PetscArraycpy(buf + k, gindices + j, nuniq - j)); 219 nuniq = k + nuniq - j; 220 } else nuniq = k; 221 /* Swap gindices and buf to merge garray of the next matrix */ 222 tmp = gindices; 223 gindices = buf; 224 buf = tmp; 225 } 226 PetscCall(PetscFree(buf)); 227 228 /* Go through matrices third time to build a map from gindices[] to garray[] */ 229 tot = 0; 230 for (cur = shell->head, j = 0; cur; cur = cur->next, j++) { /* j-th matrix */ 231 PetscCall(MatMPIAIJGetSeqAIJ(cur->mat, NULL, &B, &garray)); 232 PetscCall(MatGetLocalSize(B, NULL, &n)); 233 PetscCall(VecCreateSeqWithArray(PETSC_COMM_SELF, 1, n, NULL, &shell->lvecs[j])); 234 /* This is an optimized PetscFindInt(garray[i],nuniq,gindices,&shell->location[tot+i]), using the fact that garray[] is also sorted */ 235 lo = 0; 236 for (i = 0; i < n; i++) { 237 hi = nuniq; 238 while (hi - lo > 1) { 239 mid = lo + (hi - lo) / 2; 240 if (garray[i] < gindices[mid]) hi = mid; 241 else lo = mid; 242 } 243 shell->location[tot + i] = lo; /* gindices[lo] = garray[i] */ 244 lo++; /* Since garray[i+1] > garray[i], we can safely advance lo */ 245 } 246 tot += n; 247 } 248 249 /* Build merged Mvctx */ 250 PetscCall(ISCreateGeneral(PETSC_COMM_SELF, nuniq, gindices, PETSC_OWN_POINTER, &ix)); 251 PetscCall(ISCreateStride(PETSC_COMM_SELF, nuniq, 0, 1, &iy)); 252 PetscCall(VecCreateMPIWithArray(PetscObjectComm((PetscObject)mat), 1, mat->cmap->n, mat->cmap->N, NULL, &xin)); 253 PetscCall(VecCreateSeq(PETSC_COMM_SELF, nuniq, &shell->gvec)); 254 PetscCall(VecScatterCreate(xin, ix, shell->gvec, iy, &shell->Mvctx)); 255 PetscCall(VecDestroy(&xin)); 256 PetscCall(ISDestroy(&ix)); 257 PetscCall(ISDestroy(&iy)); 258 } 259 260 skip_merge_mvctx: 261 PetscCall(VecSet(y, 0)); 262 if (!shell->leftwork2) PetscCall(VecDuplicate(y, &shell->leftwork2)); 263 y2 = shell->leftwork2; 264 265 if (shell->Mvctx) { /* Have a merged Mvctx */ 266 /* Suppose we want to compute y = sMx, where s is the scaling factor and A, B are matrix M's diagonal/off-diagonal part. We could do 267 in y = s(Ax1 + Bx2) or y = sAx1 + sBx2. The former incurs less FLOPS than the latter, but the latter provides an opportunity to 268 overlap communication/computation since we can do sAx1 while communicating x2. Here, we use the former approach. 269 */ 270 PetscCall(VecScatterBegin(shell->Mvctx, in, shell->gvec, INSERT_VALUES, SCATTER_FORWARD)); 271 PetscCall(VecScatterEnd(shell->Mvctx, in, shell->gvec, INSERT_VALUES, SCATTER_FORWARD)); 272 273 PetscCall(VecGetArrayRead(shell->gvec, &vals)); 274 for (i = 0; i < shell->len; i++) shell->larray[i] = vals[shell->location[i]]; 275 PetscCall(VecRestoreArrayRead(shell->gvec, &vals)); 276 277 for (cur = shell->head, tot = i = 0; cur; cur = cur->next, i++) { /* i-th matrix */ 278 PetscCall(MatMPIAIJGetSeqAIJ(cur->mat, &A, &B, NULL)); 279 PetscUseTypeMethod(A, mult, in, y2); 280 PetscCall(MatGetLocalSize(B, NULL, &n)); 281 PetscCall(VecPlaceArray(shell->lvecs[i], &shell->larray[tot])); 282 PetscCall((*B->ops->multadd)(B, shell->lvecs[i], y2, y2)); 283 PetscCall(VecResetArray(shell->lvecs[i])); 284 PetscCall(VecAXPY(y, (shell->scalings ? shell->scalings[i] : 1.0), y2)); 285 tot += n; 286 } 287 } else { 288 if (shell->scalings) { 289 for (cur = shell->head, i = 0; cur; cur = cur->next, i++) { 290 PetscCall(MatMult(cur->mat, in, y2)); 291 PetscCall(VecAXPY(y, shell->scalings[i], y2)); 292 } 293 } else { 294 for (cur = shell->head; cur; cur = cur->next) PetscCall(MatMultAdd(cur->mat, in, y, y)); 295 } 296 } 297 298 if (shell->left) PetscCall(VecPointwiseMult(y, shell->left, y)); 299 PetscCall(VecScale(y, shell->scale)); 300 PetscFunctionReturn(PETSC_SUCCESS); 301 } 302 303 PetscErrorCode MatMultTranspose_Composite(Mat A, Vec x, Vec y) 304 { 305 Mat_Composite *shell = (Mat_Composite *)A->data; 306 Mat_CompositeLink next = shell->head; 307 Vec in, y2 = NULL; 308 PetscInt i; 309 310 PetscFunctionBegin; 311 PetscCheck(next, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Must provide at least one matrix with MatCompositeAddMat()"); 312 in = x; 313 if (shell->left) { 314 if (!shell->leftwork) PetscCall(VecDuplicate(shell->left, &shell->leftwork)); 315 PetscCall(VecPointwiseMult(shell->leftwork, shell->left, in)); 316 in = shell->leftwork; 317 } 318 319 PetscCall(MatMultTranspose(next->mat, in, y)); 320 if (shell->scalings) { 321 PetscCall(VecScale(y, shell->scalings[0])); 322 if (!shell->rightwork2) PetscCall(VecDuplicate(y, &shell->rightwork2)); 323 y2 = shell->rightwork2; 324 } 325 i = 1; 326 while ((next = next->next)) { 327 if (!shell->scalings) PetscCall(MatMultTransposeAdd(next->mat, in, y, y)); 328 else { 329 PetscCall(MatMultTranspose(next->mat, in, y2)); 330 PetscCall(VecAXPY(y, shell->scalings[i++], y2)); 331 } 332 } 333 if (shell->right) PetscCall(VecPointwiseMult(y, shell->right, y)); 334 PetscCall(VecScale(y, shell->scale)); 335 PetscFunctionReturn(PETSC_SUCCESS); 336 } 337 338 PetscErrorCode MatMultAdd_Composite(Mat A, Vec x, Vec y, Vec z) 339 { 340 Mat_Composite *shell = (Mat_Composite *)A->data; 341 342 PetscFunctionBegin; 343 if (y != z) { 344 PetscCall(MatMult(A, x, z)); 345 PetscCall(VecAXPY(z, 1.0, y)); 346 } else { 347 if (!shell->leftwork) PetscCall(VecDuplicate(z, &shell->leftwork)); 348 PetscCall(MatMult(A, x, shell->leftwork)); 349 PetscCall(VecCopy(y, z)); 350 PetscCall(VecAXPY(z, 1.0, shell->leftwork)); 351 } 352 PetscFunctionReturn(PETSC_SUCCESS); 353 } 354 355 PetscErrorCode MatMultTransposeAdd_Composite(Mat A, Vec x, Vec y, Vec z) 356 { 357 Mat_Composite *shell = (Mat_Composite *)A->data; 358 359 PetscFunctionBegin; 360 if (y != z) { 361 PetscCall(MatMultTranspose(A, x, z)); 362 PetscCall(VecAXPY(z, 1.0, y)); 363 } else { 364 if (!shell->rightwork) PetscCall(VecDuplicate(z, &shell->rightwork)); 365 PetscCall(MatMultTranspose(A, x, shell->rightwork)); 366 PetscCall(VecCopy(y, z)); 367 PetscCall(VecAXPY(z, 1.0, shell->rightwork)); 368 } 369 PetscFunctionReturn(PETSC_SUCCESS); 370 } 371 372 PetscErrorCode MatGetDiagonal_Composite(Mat A, Vec v) 373 { 374 Mat_Composite *shell = (Mat_Composite *)A->data; 375 Mat_CompositeLink next = shell->head; 376 PetscInt i; 377 378 PetscFunctionBegin; 379 PetscCheck(next, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Must provide at least one matrix with MatCompositeAddMat()"); 380 PetscCheck(!shell->right && !shell->left, PETSC_COMM_SELF, PETSC_ERR_SUP, "Cannot get diagonal if left or right scaling"); 381 382 PetscCall(MatGetDiagonal(next->mat, v)); 383 if (shell->scalings) PetscCall(VecScale(v, shell->scalings[0])); 384 385 if (next->next && !shell->work) PetscCall(VecDuplicate(v, &shell->work)); 386 i = 1; 387 while ((next = next->next)) { 388 PetscCall(MatGetDiagonal(next->mat, shell->work)); 389 PetscCall(VecAXPY(v, (shell->scalings ? shell->scalings[i++] : 1.0), shell->work)); 390 } 391 PetscCall(VecScale(v, shell->scale)); 392 PetscFunctionReturn(PETSC_SUCCESS); 393 } 394 395 PetscErrorCode MatAssemblyEnd_Composite(Mat Y, MatAssemblyType t) 396 { 397 Mat_Composite *shell = (Mat_Composite *)Y->data; 398 399 PetscFunctionBegin; 400 if (shell->merge) PetscCall(MatCompositeMerge(Y)); 401 PetscFunctionReturn(PETSC_SUCCESS); 402 } 403 404 PetscErrorCode MatScale_Composite(Mat inA, PetscScalar alpha) 405 { 406 Mat_Composite *a = (Mat_Composite *)inA->data; 407 408 PetscFunctionBegin; 409 a->scale *= alpha; 410 PetscFunctionReturn(PETSC_SUCCESS); 411 } 412 413 PetscErrorCode MatDiagonalScale_Composite(Mat inA, Vec left, Vec right) 414 { 415 Mat_Composite *a = (Mat_Composite *)inA->data; 416 417 PetscFunctionBegin; 418 if (left) { 419 if (!a->left) { 420 PetscCall(VecDuplicate(left, &a->left)); 421 PetscCall(VecCopy(left, a->left)); 422 } else { 423 PetscCall(VecPointwiseMult(a->left, left, a->left)); 424 } 425 } 426 if (right) { 427 if (!a->right) { 428 PetscCall(VecDuplicate(right, &a->right)); 429 PetscCall(VecCopy(right, a->right)); 430 } else { 431 PetscCall(VecPointwiseMult(a->right, right, a->right)); 432 } 433 } 434 PetscFunctionReturn(PETSC_SUCCESS); 435 } 436 437 PetscErrorCode MatSetFromOptions_Composite(Mat A, PetscOptionItems *PetscOptionsObject) 438 { 439 Mat_Composite *a = (Mat_Composite *)A->data; 440 441 PetscFunctionBegin; 442 PetscOptionsHeadBegin(PetscOptionsObject, "MATCOMPOSITE options"); 443 PetscCall(PetscOptionsBool("-mat_composite_merge", "Merge at MatAssemblyEnd", "MatCompositeMerge", a->merge, &a->merge, NULL)); 444 PetscCall(PetscOptionsEnum("-mat_composite_merge_type", "Set composite merge direction", "MatCompositeSetMergeType", MatCompositeMergeTypes, (PetscEnum)a->mergetype, (PetscEnum *)&a->mergetype, NULL)); 445 PetscCall(PetscOptionsBool("-mat_composite_merge_mvctx", "Merge MatMult() vecscat contexts", "MatCreateComposite", a->merge_mvctx, &a->merge_mvctx, NULL)); 446 PetscOptionsHeadEnd(); 447 PetscFunctionReturn(PETSC_SUCCESS); 448 } 449 450 /*@ 451 MatCreateComposite - Creates a matrix as the sum or product of one or more matrices 452 453 Collective 454 455 Input Parameters: 456 + comm - MPI communicator 457 . nmat - number of matrices to put in 458 - mats - the matrices 459 460 Output Parameter: 461 . mat - the matrix 462 463 Options Database Keys: 464 + -mat_composite_merge - merge in `MatAssemblyEnd()` 465 . -mat_composite_merge_mvctx - merge Mvctx of component matrices to optimize communication in `MatMult()` for ADDITIVE matrices 466 - -mat_composite_merge_type - set merge direction 467 468 Level: advanced 469 470 Note: 471 Alternative construction 472 $ MatCreate(comm,&mat); 473 $ MatSetSizes(mat,m,n,M,N); 474 $ MatSetType(mat,MATCOMPOSITE); 475 $ MatCompositeAddMat(mat,mats[0]); 476 $ .... 477 $ MatCompositeAddMat(mat,mats[nmat-1]); 478 $ MatAssemblyBegin(mat,MAT_FINAL_ASSEMBLY); 479 $ MatAssemblyEnd(mat,MAT_FINAL_ASSEMBLY); 480 481 For the multiplicative form the product is mat[nmat-1]*mat[nmat-2]*....*mat[0] 482 483 .seealso: `MatDestroy()`, `MatMult()`, `MatCompositeAddMat()`, `MatCompositeGetMat()`, `MatCompositeMerge()`, `MatCompositeSetType()`, `MATCOMPOSITE` 484 @*/ 485 PetscErrorCode MatCreateComposite(MPI_Comm comm, PetscInt nmat, const Mat *mats, Mat *mat) 486 { 487 PetscInt m, n, M, N, i; 488 489 PetscFunctionBegin; 490 PetscCheck(nmat >= 1, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Must pass in at least one matrix"); 491 PetscValidPointer(mat, 4); 492 493 PetscCall(MatGetLocalSize(mats[0], PETSC_IGNORE, &n)); 494 PetscCall(MatGetLocalSize(mats[nmat - 1], &m, PETSC_IGNORE)); 495 PetscCall(MatGetSize(mats[0], PETSC_IGNORE, &N)); 496 PetscCall(MatGetSize(mats[nmat - 1], &M, PETSC_IGNORE)); 497 PetscCall(MatCreate(comm, mat)); 498 PetscCall(MatSetSizes(*mat, m, n, M, N)); 499 PetscCall(MatSetType(*mat, MATCOMPOSITE)); 500 for (i = 0; i < nmat; i++) PetscCall(MatCompositeAddMat(*mat, mats[i])); 501 PetscCall(MatAssemblyBegin(*mat, MAT_FINAL_ASSEMBLY)); 502 PetscCall(MatAssemblyEnd(*mat, MAT_FINAL_ASSEMBLY)); 503 PetscFunctionReturn(PETSC_SUCCESS); 504 } 505 506 static PetscErrorCode MatCompositeAddMat_Composite(Mat mat, Mat smat) 507 { 508 Mat_Composite *shell = (Mat_Composite *)mat->data; 509 Mat_CompositeLink ilink, next = shell->head; 510 511 PetscFunctionBegin; 512 PetscCall(PetscNew(&ilink)); 513 ilink->next = NULL; 514 PetscCall(PetscObjectReference((PetscObject)smat)); 515 ilink->mat = smat; 516 517 if (!next) shell->head = ilink; 518 else { 519 while (next->next) next = next->next; 520 next->next = ilink; 521 ilink->prev = next; 522 } 523 shell->tail = ilink; 524 shell->nmat += 1; 525 526 /* Retain the old scalings (if any) and expand it with a 1.0 for the newly added matrix */ 527 if (shell->scalings) { 528 PetscCall(PetscRealloc(sizeof(PetscScalar) * shell->nmat, &shell->scalings)); 529 shell->scalings[shell->nmat - 1] = 1.0; 530 } 531 PetscFunctionReturn(PETSC_SUCCESS); 532 } 533 534 /*@ 535 MatCompositeAddMat - Add another matrix to a composite matrix. 536 537 Collective 538 539 Input Parameters: 540 + mat - the composite matrix 541 - smat - the partial matrix 542 543 Level: advanced 544 545 .seealso: `Mat`, `MatCreateComposite()`, `MatCompositeGetMat()`, `MATCOMPOSITE` 546 @*/ 547 PetscErrorCode MatCompositeAddMat(Mat mat, Mat smat) 548 { 549 PetscFunctionBegin; 550 PetscValidHeaderSpecific(mat, MAT_CLASSID, 1); 551 PetscValidHeaderSpecific(smat, MAT_CLASSID, 2); 552 PetscUseMethod(mat, "MatCompositeAddMat_C", (Mat, Mat), (mat, smat)); 553 PetscFunctionReturn(PETSC_SUCCESS); 554 } 555 556 static PetscErrorCode MatCompositeSetType_Composite(Mat mat, MatCompositeType type) 557 { 558 Mat_Composite *b = (Mat_Composite *)mat->data; 559 560 PetscFunctionBegin; 561 b->type = type; 562 if (type == MAT_COMPOSITE_MULTIPLICATIVE) { 563 mat->ops->getdiagonal = NULL; 564 mat->ops->mult = MatMult_Composite_Multiplicative; 565 mat->ops->multtranspose = MatMultTranspose_Composite_Multiplicative; 566 b->merge_mvctx = PETSC_FALSE; 567 } else { 568 mat->ops->getdiagonal = MatGetDiagonal_Composite; 569 mat->ops->mult = MatMult_Composite; 570 mat->ops->multtranspose = MatMultTranspose_Composite; 571 } 572 PetscFunctionReturn(PETSC_SUCCESS); 573 } 574 575 /*@ 576 MatCompositeSetType - Indicates if the matrix is defined as the sum of a set of matrices or the product. 577 578 Logically Collective 579 580 Input Parameters: 581 . mat - the composite matrix 582 583 Level: advanced 584 585 .seealso: `Mat`, `MatDestroy()`, `MatMult()`, `MatCompositeAddMat()`, `MatCreateComposite()`, `MatCompositeGetType()`, `MATCOMPOSITE` 586 @*/ 587 PetscErrorCode MatCompositeSetType(Mat mat, MatCompositeType type) 588 { 589 PetscFunctionBegin; 590 PetscValidHeaderSpecific(mat, MAT_CLASSID, 1); 591 PetscValidLogicalCollectiveEnum(mat, type, 2); 592 PetscUseMethod(mat, "MatCompositeSetType_C", (Mat, MatCompositeType), (mat, type)); 593 PetscFunctionReturn(PETSC_SUCCESS); 594 } 595 596 static PetscErrorCode MatCompositeGetType_Composite(Mat mat, MatCompositeType *type) 597 { 598 Mat_Composite *b = (Mat_Composite *)mat->data; 599 600 PetscFunctionBegin; 601 *type = b->type; 602 PetscFunctionReturn(PETSC_SUCCESS); 603 } 604 605 /*@ 606 MatCompositeGetType - Returns type of composite. 607 608 Not Collective 609 610 Input Parameter: 611 . mat - the composite matrix 612 613 Output Parameter: 614 . type - type of composite 615 616 Level: advanced 617 618 .seealso: `MatCreateComposite()`, `MatCompositeSetType()`, `MATCOMPOSITE`, `MatCompositeType` 619 @*/ 620 PetscErrorCode MatCompositeGetType(Mat mat, MatCompositeType *type) 621 { 622 PetscFunctionBegin; 623 PetscValidHeaderSpecific(mat, MAT_CLASSID, 1); 624 PetscValidPointer(type, 2); 625 PetscUseMethod(mat, "MatCompositeGetType_C", (Mat, MatCompositeType *), (mat, type)); 626 PetscFunctionReturn(PETSC_SUCCESS); 627 } 628 629 static PetscErrorCode MatCompositeSetMatStructure_Composite(Mat mat, MatStructure str) 630 { 631 Mat_Composite *b = (Mat_Composite *)mat->data; 632 633 PetscFunctionBegin; 634 b->structure = str; 635 PetscFunctionReturn(PETSC_SUCCESS); 636 } 637 638 /*@ 639 MatCompositeSetMatStructure - Indicates structure of matrices in the composite matrix. 640 641 Not Collective 642 643 Input Parameters: 644 + mat - the composite matrix 645 - str - either `SAME_NONZERO_PATTERN`, `DIFFERENT_NONZERO_PATTERN` (default) or `SUBSET_NONZERO_PATTERN` 646 647 Level: advanced 648 649 Note: 650 Information about the matrices structure is used in `MatCompositeMerge()` for additive composite matrix. 651 652 .seealso: `MatAXPY()`, `MatCreateComposite()`, `MatCompositeMerge()` `MatCompositeGetMatStructure()`, `MATCOMPOSITE` 653 @*/ 654 PetscErrorCode MatCompositeSetMatStructure(Mat mat, MatStructure str) 655 { 656 PetscFunctionBegin; 657 PetscValidHeaderSpecific(mat, MAT_CLASSID, 1); 658 PetscUseMethod(mat, "MatCompositeSetMatStructure_C", (Mat, MatStructure), (mat, str)); 659 PetscFunctionReturn(PETSC_SUCCESS); 660 } 661 662 static PetscErrorCode MatCompositeGetMatStructure_Composite(Mat mat, MatStructure *str) 663 { 664 Mat_Composite *b = (Mat_Composite *)mat->data; 665 666 PetscFunctionBegin; 667 *str = b->structure; 668 PetscFunctionReturn(PETSC_SUCCESS); 669 } 670 671 /*@ 672 MatCompositeGetMatStructure - Returns the structure of matrices in the composite matrix. 673 674 Not Collective 675 676 Input Parameter: 677 . mat - the composite matrix 678 679 Output Parameter: 680 . str - structure of the matrices 681 682 Level: advanced 683 684 .seealso: `MatCreateComposite()`, `MatCompositeSetMatStructure()`, `MATCOMPOSITE` 685 @*/ 686 PetscErrorCode MatCompositeGetMatStructure(Mat mat, MatStructure *str) 687 { 688 PetscFunctionBegin; 689 PetscValidHeaderSpecific(mat, MAT_CLASSID, 1); 690 PetscValidPointer(str, 2); 691 PetscUseMethod(mat, "MatCompositeGetMatStructure_C", (Mat, MatStructure *), (mat, str)); 692 PetscFunctionReturn(PETSC_SUCCESS); 693 } 694 695 static PetscErrorCode MatCompositeSetMergeType_Composite(Mat mat, MatCompositeMergeType type) 696 { 697 Mat_Composite *shell = (Mat_Composite *)mat->data; 698 699 PetscFunctionBegin; 700 shell->mergetype = type; 701 PetscFunctionReturn(PETSC_SUCCESS); 702 } 703 704 /*@ 705 MatCompositeSetMergeType - Sets order of `MatCompositeMerge()`. 706 707 Logically Collective 708 709 Input Parameters: 710 + mat - the composite matrix 711 - type - `MAT_COMPOSITE_MERGE RIGHT` (default) to start merge from right with the first added matrix (mat[0]), 712 `MAT_COMPOSITE_MERGE_LEFT` to start merge from left with the last added matrix (mat[nmat-1]) 713 714 Level: advanced 715 716 Note: 717 The resulting matrix is the same regardless of the `MatCompositeMergeType`. Only the order of operation is changed. 718 If set to `MAT_COMPOSITE_MERGE_RIGHT` the order of the merge is mat[nmat-1]*(mat[nmat-2]*(...*(mat[1]*mat[0]))) 719 otherwise the order is (((mat[nmat-1]*mat[nmat-2])*mat[nmat-3])*...)*mat[0]. 720 721 .seealso: `MatCreateComposite()`, `MatCompositeMerge()`, `MATCOMPOSITE` 722 @*/ 723 PetscErrorCode MatCompositeSetMergeType(Mat mat, MatCompositeMergeType type) 724 { 725 PetscFunctionBegin; 726 PetscValidHeaderSpecific(mat, MAT_CLASSID, 1); 727 PetscValidLogicalCollectiveEnum(mat, type, 2); 728 PetscUseMethod(mat, "MatCompositeSetMergeType_C", (Mat, MatCompositeMergeType), (mat, type)); 729 PetscFunctionReturn(PETSC_SUCCESS); 730 } 731 732 static PetscErrorCode MatCompositeMerge_Composite(Mat mat) 733 { 734 Mat_Composite *shell = (Mat_Composite *)mat->data; 735 Mat_CompositeLink next = shell->head, prev = shell->tail; 736 Mat tmat, newmat; 737 Vec left, right; 738 PetscScalar scale; 739 PetscInt i; 740 741 PetscFunctionBegin; 742 PetscCheck(next, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Must provide at least one matrix with MatCompositeAddMat()"); 743 scale = shell->scale; 744 if (shell->type == MAT_COMPOSITE_ADDITIVE) { 745 if (shell->mergetype == MAT_COMPOSITE_MERGE_RIGHT) { 746 i = 0; 747 PetscCall(MatDuplicate(next->mat, MAT_COPY_VALUES, &tmat)); 748 if (shell->scalings) PetscCall(MatScale(tmat, shell->scalings[i++])); 749 while ((next = next->next)) PetscCall(MatAXPY(tmat, (shell->scalings ? shell->scalings[i++] : 1.0), next->mat, shell->structure)); 750 } else { 751 i = shell->nmat - 1; 752 PetscCall(MatDuplicate(prev->mat, MAT_COPY_VALUES, &tmat)); 753 if (shell->scalings) PetscCall(MatScale(tmat, shell->scalings[i--])); 754 while ((prev = prev->prev)) PetscCall(MatAXPY(tmat, (shell->scalings ? shell->scalings[i--] : 1.0), prev->mat, shell->structure)); 755 } 756 } else { 757 if (shell->mergetype == MAT_COMPOSITE_MERGE_RIGHT) { 758 PetscCall(MatDuplicate(next->mat, MAT_COPY_VALUES, &tmat)); 759 while ((next = next->next)) { 760 PetscCall(MatMatMult(next->mat, tmat, MAT_INITIAL_MATRIX, PETSC_DECIDE, &newmat)); 761 PetscCall(MatDestroy(&tmat)); 762 tmat = newmat; 763 } 764 } else { 765 PetscCall(MatDuplicate(prev->mat, MAT_COPY_VALUES, &tmat)); 766 while ((prev = prev->prev)) { 767 PetscCall(MatMatMult(tmat, prev->mat, MAT_INITIAL_MATRIX, PETSC_DECIDE, &newmat)); 768 PetscCall(MatDestroy(&tmat)); 769 tmat = newmat; 770 } 771 } 772 if (shell->scalings) { 773 for (i = 0; i < shell->nmat; i++) scale *= shell->scalings[i]; 774 } 775 } 776 777 if ((left = shell->left)) PetscCall(PetscObjectReference((PetscObject)left)); 778 if ((right = shell->right)) PetscCall(PetscObjectReference((PetscObject)right)); 779 780 PetscCall(MatHeaderReplace(mat, &tmat)); 781 782 PetscCall(MatDiagonalScale(mat, left, right)); 783 PetscCall(MatScale(mat, scale)); 784 PetscCall(VecDestroy(&left)); 785 PetscCall(VecDestroy(&right)); 786 PetscFunctionReturn(PETSC_SUCCESS); 787 } 788 789 /*@ 790 MatCompositeMerge - Given a composite matrix, replaces it with a "regular" matrix 791 by summing or computing the product of all the matrices inside the composite matrix. 792 793 Collective 794 795 Input Parameter: 796 . mat - the composite matrix 797 798 Options Database Keys: 799 + -mat_composite_merge - merge in `MatAssemblyEnd()` 800 - -mat_composite_merge_type - set merge direction 801 802 Level: advanced 803 804 Note: 805 The `MatType` of the resulting matrix will be the same as the `MatType` of the FIRST matrix in the composite matrix. 806 807 .seealso: `MatDestroy()`, `MatMult()`, `MatCompositeAddMat()`, `MatCreateComposite()`, `MatCompositeSetMatStructure()`, `MatCompositeSetMergeType()`, `MATCOMPOSITE` 808 @*/ 809 PetscErrorCode MatCompositeMerge(Mat mat) 810 { 811 PetscFunctionBegin; 812 PetscValidHeaderSpecific(mat, MAT_CLASSID, 1); 813 PetscUseMethod(mat, "MatCompositeMerge_C", (Mat), (mat)); 814 PetscFunctionReturn(PETSC_SUCCESS); 815 } 816 817 static PetscErrorCode MatCompositeGetNumberMat_Composite(Mat mat, PetscInt *nmat) 818 { 819 Mat_Composite *shell = (Mat_Composite *)mat->data; 820 821 PetscFunctionBegin; 822 *nmat = shell->nmat; 823 PetscFunctionReturn(PETSC_SUCCESS); 824 } 825 826 /*@ 827 MatCompositeGetNumberMat - Returns the number of matrices in the composite matrix. 828 829 Not Collective 830 831 Input Parameter: 832 . mat - the composite matrix 833 834 Output Parameter: 835 . nmat - number of matrices in the composite matrix 836 837 Level: advanced 838 839 .seealso: `MatCreateComposite()`, `MatCompositeGetMat()`, `MATCOMPOSITE` 840 @*/ 841 PetscErrorCode MatCompositeGetNumberMat(Mat mat, PetscInt *nmat) 842 { 843 PetscFunctionBegin; 844 PetscValidHeaderSpecific(mat, MAT_CLASSID, 1); 845 PetscValidIntPointer(nmat, 2); 846 PetscUseMethod(mat, "MatCompositeGetNumberMat_C", (Mat, PetscInt *), (mat, nmat)); 847 PetscFunctionReturn(PETSC_SUCCESS); 848 } 849 850 static PetscErrorCode MatCompositeGetMat_Composite(Mat mat, PetscInt i, Mat *Ai) 851 { 852 Mat_Composite *shell = (Mat_Composite *)mat->data; 853 Mat_CompositeLink ilink; 854 PetscInt k; 855 856 PetscFunctionBegin; 857 PetscCheck(i < shell->nmat, PetscObjectComm((PetscObject)mat), PETSC_ERR_ARG_OUTOFRANGE, "index out of range: %" PetscInt_FMT " >= %" PetscInt_FMT, i, shell->nmat); 858 ilink = shell->head; 859 for (k = 0; k < i; k++) ilink = ilink->next; 860 *Ai = ilink->mat; 861 PetscFunctionReturn(PETSC_SUCCESS); 862 } 863 864 /*@ 865 MatCompositeGetMat - Returns the ith matrix from the composite matrix. 866 867 Logically Collective 868 869 Input Parameters: 870 + mat - the composite matrix 871 - i - the number of requested matrix 872 873 Output Parameter: 874 . Ai - ith matrix in composite 875 876 Level: advanced 877 878 .seealso: `MatCreateComposite()`, `MatCompositeGetNumberMat()`, `MatCompositeAddMat()`, `MATCOMPOSITE` 879 @*/ 880 PetscErrorCode MatCompositeGetMat(Mat mat, PetscInt i, Mat *Ai) 881 { 882 PetscFunctionBegin; 883 PetscValidHeaderSpecific(mat, MAT_CLASSID, 1); 884 PetscValidLogicalCollectiveInt(mat, i, 2); 885 PetscValidPointer(Ai, 3); 886 PetscUseMethod(mat, "MatCompositeGetMat_C", (Mat, PetscInt, Mat *), (mat, i, Ai)); 887 PetscFunctionReturn(PETSC_SUCCESS); 888 } 889 890 PetscErrorCode MatCompositeSetScalings_Composite(Mat mat, const PetscScalar *scalings) 891 { 892 Mat_Composite *shell = (Mat_Composite *)mat->data; 893 PetscInt nmat; 894 895 PetscFunctionBegin; 896 PetscCall(MatCompositeGetNumberMat(mat, &nmat)); 897 if (!shell->scalings) PetscCall(PetscMalloc1(nmat, &shell->scalings)); 898 PetscCall(PetscArraycpy(shell->scalings, scalings, nmat)); 899 PetscFunctionReturn(PETSC_SUCCESS); 900 } 901 902 /*@ 903 MatCompositeSetScalings - Sets separate scaling factors for component matrices. 904 905 Logically Collective 906 907 Input Parameters: 908 + mat - the composite matrix 909 - scalings - array of scaling factors with scalings[i] being factor of i-th matrix, for i in [0, nmat) 910 911 Level: advanced 912 913 .seealso: `MatScale()`, `MatDiagonalScale()`, `MATCOMPOSITE` 914 @*/ 915 PetscErrorCode MatCompositeSetScalings(Mat mat, const PetscScalar *scalings) 916 { 917 PetscFunctionBegin; 918 PetscValidHeaderSpecific(mat, MAT_CLASSID, 1); 919 PetscValidScalarPointer(scalings, 2); 920 PetscValidLogicalCollectiveScalar(mat, *scalings, 2); 921 PetscUseMethod(mat, "MatCompositeSetScalings_C", (Mat, const PetscScalar *), (mat, scalings)); 922 PetscFunctionReturn(PETSC_SUCCESS); 923 } 924 925 static struct _MatOps MatOps_Values = {NULL, 926 NULL, 927 NULL, 928 MatMult_Composite, 929 MatMultAdd_Composite, 930 /* 5*/ MatMultTranspose_Composite, 931 MatMultTransposeAdd_Composite, 932 NULL, 933 NULL, 934 NULL, 935 /* 10*/ NULL, 936 NULL, 937 NULL, 938 NULL, 939 NULL, 940 /* 15*/ NULL, 941 NULL, 942 MatGetDiagonal_Composite, 943 MatDiagonalScale_Composite, 944 NULL, 945 /* 20*/ NULL, 946 MatAssemblyEnd_Composite, 947 NULL, 948 NULL, 949 /* 24*/ NULL, 950 NULL, 951 NULL, 952 NULL, 953 NULL, 954 /* 29*/ NULL, 955 NULL, 956 NULL, 957 NULL, 958 NULL, 959 /* 34*/ NULL, 960 NULL, 961 NULL, 962 NULL, 963 NULL, 964 /* 39*/ NULL, 965 NULL, 966 NULL, 967 NULL, 968 NULL, 969 /* 44*/ NULL, 970 MatScale_Composite, 971 MatShift_Basic, 972 NULL, 973 NULL, 974 /* 49*/ NULL, 975 NULL, 976 NULL, 977 NULL, 978 NULL, 979 /* 54*/ NULL, 980 NULL, 981 NULL, 982 NULL, 983 NULL, 984 /* 59*/ NULL, 985 MatDestroy_Composite, 986 NULL, 987 NULL, 988 NULL, 989 /* 64*/ NULL, 990 NULL, 991 NULL, 992 NULL, 993 NULL, 994 /* 69*/ NULL, 995 NULL, 996 NULL, 997 NULL, 998 NULL, 999 /* 74*/ NULL, 1000 NULL, 1001 MatSetFromOptions_Composite, 1002 NULL, 1003 NULL, 1004 /* 79*/ NULL, 1005 NULL, 1006 NULL, 1007 NULL, 1008 NULL, 1009 /* 84*/ NULL, 1010 NULL, 1011 NULL, 1012 NULL, 1013 NULL, 1014 /* 89*/ NULL, 1015 NULL, 1016 NULL, 1017 NULL, 1018 NULL, 1019 /* 94*/ NULL, 1020 NULL, 1021 NULL, 1022 NULL, 1023 NULL, 1024 /*99*/ NULL, 1025 NULL, 1026 NULL, 1027 NULL, 1028 NULL, 1029 /*104*/ NULL, 1030 NULL, 1031 NULL, 1032 NULL, 1033 NULL, 1034 /*109*/ NULL, 1035 NULL, 1036 NULL, 1037 NULL, 1038 NULL, 1039 /*114*/ NULL, 1040 NULL, 1041 NULL, 1042 NULL, 1043 NULL, 1044 /*119*/ NULL, 1045 NULL, 1046 NULL, 1047 NULL, 1048 NULL, 1049 /*124*/ NULL, 1050 NULL, 1051 NULL, 1052 NULL, 1053 NULL, 1054 /*129*/ NULL, 1055 NULL, 1056 NULL, 1057 NULL, 1058 NULL, 1059 /*134*/ NULL, 1060 NULL, 1061 NULL, 1062 NULL, 1063 NULL, 1064 /*139*/ NULL, 1065 NULL, 1066 NULL, 1067 NULL, 1068 NULL, 1069 /*144*/ NULL, 1070 NULL, 1071 NULL, 1072 NULL, 1073 NULL, 1074 NULL, 1075 /*150*/ NULL, 1076 NULL}; 1077 1078 /*MC 1079 MATCOMPOSITE - A matrix defined by the sum (or product) of one or more matrices. 1080 The matrices need to have a correct size and parallel layout for the sum or product to be valid. 1081 1082 Note: 1083 To use the product of the matrices call `MatCompositeSetType`(mat,`MAT_COMPOSITE_MULTIPLICATIVE`); 1084 1085 Level: advanced 1086 1087 .seealso: `MatCreateComposite()`, `MatCompositeSetScalings()`, `MatCompositeAddMat()`, `MatSetType()`, `MatCompositeSetType()`, `MatCompositeGetType()`, 1088 `MatCompositeSetMatStructure()`, `MatCompositeGetMatStructure()`, `MatCompositeMerge()`, `MatCompositeSetMergeType()`, `MatCompositeGetNumberMat()`, `MatCompositeGetMat()` 1089 M*/ 1090 1091 PETSC_EXTERN PetscErrorCode MatCreate_Composite(Mat A) 1092 { 1093 Mat_Composite *b; 1094 1095 PetscFunctionBegin; 1096 PetscCall(PetscNew(&b)); 1097 A->data = (void *)b; 1098 PetscCall(PetscMemcpy(A->ops, &MatOps_Values, sizeof(struct _MatOps))); 1099 1100 PetscCall(PetscLayoutSetUp(A->rmap)); 1101 PetscCall(PetscLayoutSetUp(A->cmap)); 1102 1103 A->assembled = PETSC_TRUE; 1104 A->preallocated = PETSC_TRUE; 1105 b->type = MAT_COMPOSITE_ADDITIVE; 1106 b->scale = 1.0; 1107 b->nmat = 0; 1108 b->merge = PETSC_FALSE; 1109 b->mergetype = MAT_COMPOSITE_MERGE_RIGHT; 1110 b->structure = DIFFERENT_NONZERO_PATTERN; 1111 b->merge_mvctx = PETSC_TRUE; 1112 1113 PetscCall(PetscObjectComposeFunction((PetscObject)A, "MatCompositeAddMat_C", MatCompositeAddMat_Composite)); 1114 PetscCall(PetscObjectComposeFunction((PetscObject)A, "MatCompositeSetType_C", MatCompositeSetType_Composite)); 1115 PetscCall(PetscObjectComposeFunction((PetscObject)A, "MatCompositeGetType_C", MatCompositeGetType_Composite)); 1116 PetscCall(PetscObjectComposeFunction((PetscObject)A, "MatCompositeSetMergeType_C", MatCompositeSetMergeType_Composite)); 1117 PetscCall(PetscObjectComposeFunction((PetscObject)A, "MatCompositeSetMatStructure_C", MatCompositeSetMatStructure_Composite)); 1118 PetscCall(PetscObjectComposeFunction((PetscObject)A, "MatCompositeGetMatStructure_C", MatCompositeGetMatStructure_Composite)); 1119 PetscCall(PetscObjectComposeFunction((PetscObject)A, "MatCompositeMerge_C", MatCompositeMerge_Composite)); 1120 PetscCall(PetscObjectComposeFunction((PetscObject)A, "MatCompositeGetNumberMat_C", MatCompositeGetNumberMat_Composite)); 1121 PetscCall(PetscObjectComposeFunction((PetscObject)A, "MatCompositeGetMat_C", MatCompositeGetMat_Composite)); 1122 PetscCall(PetscObjectComposeFunction((PetscObject)A, "MatCompositeSetScalings_C", MatCompositeSetScalings_Composite)); 1123 1124 PetscCall(PetscObjectChangeTypeName((PetscObject)A, MATCOMPOSITE)); 1125 PetscFunctionReturn(PETSC_SUCCESS); 1126 } 1127