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