1 2 #include <petsc/private/matimpl.h> /*I "petscmat.h" I*/ 3 4 typedef struct _Mat_CompositeLink *Mat_CompositeLink; 5 struct _Mat_CompositeLink { 6 Mat mat; 7 Vec work; 8 Mat_CompositeLink next,prev; 9 }; 10 11 typedef struct { 12 MatCompositeType type; 13 Mat_CompositeLink head,tail; 14 Vec work; 15 PetscScalar scale; /* scale factor supplied with MatScale() */ 16 Vec left,right; /* left and right diagonal scaling provided with MatDiagonalScale() */ 17 Vec leftwork,rightwork; 18 PetscInt nmat; 19 PetscBool mergefromright; 20 } Mat_Composite; 21 22 PetscErrorCode MatDestroy_Composite(Mat mat) 23 { 24 PetscErrorCode ierr; 25 Mat_Composite *shell = (Mat_Composite*)mat->data; 26 Mat_CompositeLink next = shell->head,oldnext; 27 28 PetscFunctionBegin; 29 while (next) { 30 ierr = MatDestroy(&next->mat);CHKERRQ(ierr); 31 if (next->work && (!next->next || next->work != next->next->work)) { 32 ierr = VecDestroy(&next->work);CHKERRQ(ierr); 33 } 34 oldnext = next; 35 next = next->next; 36 ierr = PetscFree(oldnext);CHKERRQ(ierr); 37 } 38 ierr = VecDestroy(&shell->work);CHKERRQ(ierr); 39 ierr = VecDestroy(&shell->left);CHKERRQ(ierr); 40 ierr = VecDestroy(&shell->right);CHKERRQ(ierr); 41 ierr = VecDestroy(&shell->leftwork);CHKERRQ(ierr); 42 ierr = VecDestroy(&shell->rightwork);CHKERRQ(ierr); 43 ierr = PetscFree(mat->data);CHKERRQ(ierr); 44 PetscFunctionReturn(0); 45 } 46 47 PetscErrorCode MatMult_Composite_Multiplicative(Mat A,Vec x,Vec y) 48 { 49 Mat_Composite *shell = (Mat_Composite*)A->data; 50 Mat_CompositeLink next = shell->head; 51 PetscErrorCode ierr; 52 Vec in,out; 53 54 PetscFunctionBegin; 55 if (!next) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must provide at least one matrix with MatCompositeAddMat()"); 56 in = x; 57 if (shell->right) { 58 if (!shell->rightwork) { 59 ierr = VecDuplicate(shell->right,&shell->rightwork);CHKERRQ(ierr); 60 } 61 ierr = VecPointwiseMult(shell->rightwork,shell->right,in);CHKERRQ(ierr); 62 in = shell->rightwork; 63 } 64 while (next->next) { 65 if (!next->work) { /* should reuse previous work if the same size */ 66 ierr = MatCreateVecs(next->mat,NULL,&next->work);CHKERRQ(ierr); 67 } 68 out = next->work; 69 ierr = MatMult(next->mat,in,out);CHKERRQ(ierr); 70 in = out; 71 next = next->next; 72 } 73 ierr = MatMult(next->mat,in,y);CHKERRQ(ierr); 74 if (shell->left) { 75 ierr = VecPointwiseMult(y,shell->left,y);CHKERRQ(ierr); 76 } 77 ierr = VecScale(y,shell->scale);CHKERRQ(ierr); 78 PetscFunctionReturn(0); 79 } 80 81 PetscErrorCode MatMultTranspose_Composite_Multiplicative(Mat A,Vec x,Vec y) 82 { 83 Mat_Composite *shell = (Mat_Composite*)A->data; 84 Mat_CompositeLink tail = shell->tail; 85 PetscErrorCode ierr; 86 Vec in,out; 87 88 PetscFunctionBegin; 89 if (!tail) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must provide at least one matrix with MatCompositeAddMat()"); 90 in = x; 91 if (shell->left) { 92 if (!shell->leftwork) { 93 ierr = VecDuplicate(shell->left,&shell->leftwork);CHKERRQ(ierr); 94 } 95 ierr = VecPointwiseMult(shell->leftwork,shell->left,in);CHKERRQ(ierr); 96 in = shell->leftwork; 97 } 98 while (tail->prev) { 99 if (!tail->prev->work) { /* should reuse previous work if the same size */ 100 ierr = MatCreateVecs(tail->mat,NULL,&tail->prev->work);CHKERRQ(ierr); 101 } 102 out = tail->prev->work; 103 ierr = MatMultTranspose(tail->mat,in,out);CHKERRQ(ierr); 104 in = out; 105 tail = tail->prev; 106 } 107 ierr = MatMultTranspose(tail->mat,in,y);CHKERRQ(ierr); 108 if (shell->right) { 109 ierr = VecPointwiseMult(y,shell->right,y);CHKERRQ(ierr); 110 } 111 ierr = VecScale(y,shell->scale);CHKERRQ(ierr); 112 PetscFunctionReturn(0); 113 } 114 115 PetscErrorCode MatMult_Composite(Mat A,Vec x,Vec y) 116 { 117 Mat_Composite *shell = (Mat_Composite*)A->data; 118 Mat_CompositeLink next = shell->head; 119 PetscErrorCode ierr; 120 Vec in; 121 122 PetscFunctionBegin; 123 if (!next) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must provide at least one matrix with MatCompositeAddMat()"); 124 in = x; 125 if (shell->right) { 126 if (!shell->rightwork) { 127 ierr = VecDuplicate(shell->right,&shell->rightwork);CHKERRQ(ierr); 128 } 129 ierr = VecPointwiseMult(shell->rightwork,shell->right,in);CHKERRQ(ierr); 130 in = shell->rightwork; 131 } 132 ierr = MatMult(next->mat,in,y);CHKERRQ(ierr); 133 while ((next = next->next)) { 134 ierr = MatMultAdd(next->mat,in,y,y);CHKERRQ(ierr); 135 } 136 if (shell->left) { 137 ierr = VecPointwiseMult(y,shell->left,y);CHKERRQ(ierr); 138 } 139 ierr = VecScale(y,shell->scale);CHKERRQ(ierr); 140 PetscFunctionReturn(0); 141 } 142 143 PetscErrorCode MatMultTranspose_Composite(Mat A,Vec x,Vec y) 144 { 145 Mat_Composite *shell = (Mat_Composite*)A->data; 146 Mat_CompositeLink next = shell->head; 147 PetscErrorCode ierr; 148 Vec in; 149 150 PetscFunctionBegin; 151 if (!next) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must provide at least one matrix with MatCompositeAddMat()"); 152 in = x; 153 if (shell->left) { 154 if (!shell->leftwork) { 155 ierr = VecDuplicate(shell->left,&shell->leftwork);CHKERRQ(ierr); 156 } 157 ierr = VecPointwiseMult(shell->leftwork,shell->left,in);CHKERRQ(ierr); 158 in = shell->leftwork; 159 } 160 ierr = MatMultTranspose(next->mat,in,y);CHKERRQ(ierr); 161 while ((next = next->next)) { 162 ierr = MatMultTransposeAdd(next->mat,in,y,y);CHKERRQ(ierr); 163 } 164 if (shell->right) { 165 ierr = VecPointwiseMult(y,shell->right,y);CHKERRQ(ierr); 166 } 167 ierr = VecScale(y,shell->scale);CHKERRQ(ierr); 168 PetscFunctionReturn(0); 169 } 170 171 PetscErrorCode MatGetDiagonal_Composite(Mat A,Vec v) 172 { 173 Mat_Composite *shell = (Mat_Composite*)A->data; 174 Mat_CompositeLink next = shell->head; 175 PetscErrorCode ierr; 176 177 PetscFunctionBegin; 178 if (!next) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must provide at least one matrix with MatCompositeAddMat()"); 179 if (shell->right || shell->left) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Cannot get diagonal if left or right scaling"); 180 181 ierr = MatGetDiagonal(next->mat,v);CHKERRQ(ierr); 182 if (next->next && !shell->work) { 183 ierr = VecDuplicate(v,&shell->work);CHKERRQ(ierr); 184 } 185 while ((next = next->next)) { 186 ierr = MatGetDiagonal(next->mat,shell->work);CHKERRQ(ierr); 187 ierr = VecAXPY(v,1.0,shell->work);CHKERRQ(ierr); 188 } 189 ierr = VecScale(v,shell->scale);CHKERRQ(ierr); 190 PetscFunctionReturn(0); 191 } 192 193 PetscErrorCode MatAssemblyEnd_Composite(Mat Y,MatAssemblyType t) 194 { 195 PetscErrorCode ierr; 196 PetscBool flg = PETSC_FALSE; 197 198 PetscFunctionBegin; 199 ierr = PetscOptionsGetBool(((PetscObject)Y)->options,((PetscObject)Y)->prefix,"-mat_composite_merge",&flg,NULL);CHKERRQ(ierr); 200 if (flg) { 201 ierr = MatCompositeMerge(Y);CHKERRQ(ierr); 202 } 203 PetscFunctionReturn(0); 204 } 205 206 PetscErrorCode MatScale_Composite(Mat inA,PetscScalar alpha) 207 { 208 Mat_Composite *a = (Mat_Composite*)inA->data; 209 210 PetscFunctionBegin; 211 a->scale *= alpha; 212 PetscFunctionReturn(0); 213 } 214 215 PetscErrorCode MatDiagonalScale_Composite(Mat inA,Vec left,Vec right) 216 { 217 Mat_Composite *a = (Mat_Composite*)inA->data; 218 PetscErrorCode ierr; 219 220 PetscFunctionBegin; 221 if (left) { 222 if (!a->left) { 223 ierr = VecDuplicate(left,&a->left);CHKERRQ(ierr); 224 ierr = VecCopy(left,a->left);CHKERRQ(ierr); 225 } else { 226 ierr = VecPointwiseMult(a->left,left,a->left);CHKERRQ(ierr); 227 } 228 } 229 if (right) { 230 if (!a->right) { 231 ierr = VecDuplicate(right,&a->right);CHKERRQ(ierr); 232 ierr = VecCopy(right,a->right);CHKERRQ(ierr); 233 } else { 234 ierr = VecPointwiseMult(a->right,right,a->right);CHKERRQ(ierr); 235 } 236 } 237 PetscFunctionReturn(0); 238 } 239 240 /*@ 241 MatCreateComposite - Creates a matrix as the sum of zero or more matrices 242 243 Collective on MPI_Comm 244 245 Input Parameters: 246 + comm - MPI communicator 247 . nmat - number of matrices to put in 248 - mats - the matrices 249 250 Output Parameter: 251 . mat - the matrix 252 253 Level: advanced 254 255 Notes: 256 Alternative construction 257 $ MatCreate(comm,&mat); 258 $ MatSetSizes(mat,m,n,M,N); 259 $ MatSetType(mat,MATCOMPOSITE); 260 $ MatCompositeAddMat(mat,mats[0]); 261 $ .... 262 $ MatCompositeAddMat(mat,mats[nmat-1]); 263 $ MatAssemblyBegin(mat,MAT_FINAL_ASSEMBLY); 264 $ MatAssemblyEnd(mat,MAT_FINAL_ASSEMBLY); 265 266 For the multiplicative form the product is mat[nmat-1]*mat[nmat-2]*....*mat[0] 267 268 .seealso: MatDestroy(), MatMult(), MatCompositeAddMat(), MatCompositeGetMat(), MatCompositeMerge(), MatCompositeSetType(), MATCOMPOSITE 269 270 @*/ 271 PetscErrorCode MatCreateComposite(MPI_Comm comm,PetscInt nmat,const Mat *mats,Mat *mat) 272 { 273 PetscErrorCode ierr; 274 PetscInt m,n,M,N,i; 275 276 PetscFunctionBegin; 277 if (nmat < 1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Must pass in at least one matrix"); 278 PetscValidPointer(mat,3); 279 280 ierr = MatGetLocalSize(mats[0],PETSC_IGNORE,&n);CHKERRQ(ierr); 281 ierr = MatGetLocalSize(mats[nmat-1],&m,PETSC_IGNORE);CHKERRQ(ierr); 282 ierr = MatGetSize(mats[0],PETSC_IGNORE,&N);CHKERRQ(ierr); 283 ierr = MatGetSize(mats[nmat-1],&M,PETSC_IGNORE);CHKERRQ(ierr); 284 ierr = MatCreate(comm,mat);CHKERRQ(ierr); 285 ierr = MatSetSizes(*mat,m,n,M,N);CHKERRQ(ierr); 286 ierr = MatSetType(*mat,MATCOMPOSITE);CHKERRQ(ierr); 287 for (i=0; i<nmat; i++) { 288 ierr = MatCompositeAddMat(*mat,mats[i]);CHKERRQ(ierr); 289 } 290 ierr = MatAssemblyBegin(*mat,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 291 ierr = MatAssemblyEnd(*mat,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 292 PetscFunctionReturn(0); 293 } 294 295 296 static PetscErrorCode MatCompositeAddMat_Composite(Mat mat,Mat smat) 297 { 298 Mat_Composite *shell = (Mat_Composite*)mat->data; 299 Mat_CompositeLink ilink,next = shell->head; 300 PetscErrorCode ierr; 301 302 PetscFunctionBegin; 303 ierr = PetscNewLog(mat,&ilink);CHKERRQ(ierr); 304 ilink->next = 0; 305 ierr = PetscObjectReference((PetscObject)smat);CHKERRQ(ierr); 306 ilink->mat = smat; 307 308 if (!next) shell->head = ilink; 309 else { 310 while (next->next) { 311 next = next->next; 312 } 313 next->next = ilink; 314 ilink->prev = next; 315 } 316 shell->tail = ilink; 317 shell->nmat += 1; 318 PetscFunctionReturn(0); 319 } 320 321 /*@ 322 MatCompositeAddMat - add another matrix to a composite matrix 323 324 Collective on Mat 325 326 Input Parameters: 327 + mat - the composite matrix 328 - smat - the partial matrix 329 330 Level: advanced 331 332 .seealso: MatCreateComposite() 333 @*/ 334 PetscErrorCode MatCompositeAddMat(Mat mat,Mat smat) 335 { 336 PetscErrorCode ierr; 337 338 PetscFunctionBegin; 339 PetscValidHeaderSpecific(mat,MAT_CLASSID,1); 340 PetscValidHeaderSpecific(smat,MAT_CLASSID,2); 341 ierr = PetscUseMethod(mat,"MatCompositeAddMat_C",(Mat,Mat),(mat,smat));CHKERRQ(ierr); 342 PetscFunctionReturn(0); 343 } 344 345 static PetscErrorCode MatCompositeSetType_Composite(Mat mat,MatCompositeType type) 346 { 347 Mat_Composite *b = (Mat_Composite*)mat->data; 348 349 PetscFunctionBegin; 350 if (type == MAT_COMPOSITE_MULTIPLICATIVE) { 351 mat->ops->getdiagonal = 0; 352 mat->ops->mult = MatMult_Composite_Multiplicative; 353 mat->ops->multtranspose = MatMultTranspose_Composite_Multiplicative; 354 b->type = MAT_COMPOSITE_MULTIPLICATIVE; 355 } else { 356 mat->ops->getdiagonal = MatGetDiagonal_Composite; 357 mat->ops->mult = MatMult_Composite; 358 mat->ops->multtranspose = MatMultTranspose_Composite; 359 b->type = MAT_COMPOSITE_ADDITIVE; 360 } 361 PetscFunctionReturn(0); 362 } 363 364 /*@ 365 MatCompositeSetType - Indicates if the matrix is defined as the sum of a set of matrices or the product 366 367 Collective on MPI_Comm 368 369 Input Parameters: 370 . mat - the composite matrix 371 372 373 Level: advanced 374 375 Notes: 376 The MatType of the resulting matrix will be the same as the MatType of the FIRST 377 matrix in the composite matrix. 378 379 .seealso: MatDestroy(), MatMult(), MatCompositeAddMat(), MatCreateComposite(), MatCompositeGetType(), MATCOMPOSITE 380 381 @*/ 382 PetscErrorCode MatCompositeSetType(Mat mat,MatCompositeType type) 383 { 384 PetscErrorCode ierr; 385 386 PetscFunctionBegin; 387 PetscValidHeaderSpecific(mat,MAT_CLASSID,1); 388 ierr = PetscUseMethod(mat,"MatCompositeSetType_C",(Mat,MatCompositeType),(mat,type));CHKERRQ(ierr); 389 PetscFunctionReturn(0); 390 } 391 392 static PetscErrorCode MatCompositeGetType_Composite(Mat mat,MatCompositeType *type) 393 { 394 Mat_Composite *b = (Mat_Composite*)mat->data; 395 396 PetscFunctionBegin; 397 *type = b->type; 398 PetscFunctionReturn(0); 399 } 400 401 /*@ 402 MatCompositeGetType - Returns type of composite. 403 404 Not Collective 405 406 Input Parameter: 407 . mat - the composite matrix 408 409 Output Parameter: 410 . type - type of composite 411 412 Level: advanced 413 414 .seealso: MatCreateComposite(), MatCompositeSetType(), MATCOMPOSITE 415 416 @*/ 417 PetscErrorCode MatCompositeGetType(Mat mat,MatCompositeType *type) 418 { 419 PetscErrorCode ierr; 420 421 PetscFunctionBegin; 422 PetscValidHeaderSpecific(mat,MAT_CLASSID,1); 423 PetscValidPointer(type,2); 424 ierr = PetscUseMethod(mat,"MatCompositeGetType_C",(Mat,MatCompositeType*),(mat,type));CHKERRQ(ierr); 425 PetscFunctionReturn(0); 426 } 427 428 static PetscErrorCode MatCompositeMerge_Composite(Mat mat) 429 { 430 Mat_Composite *shell = (Mat_Composite*)mat->data; 431 Mat_CompositeLink next = shell->head, prev = shell->tail; 432 PetscErrorCode ierr; 433 Mat tmat,newmat; 434 Vec left,right; 435 PetscScalar scale; 436 437 PetscFunctionBegin; 438 if (!next) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must provide at least one matrix with MatCompositeAddMat()"); 439 440 PetscFunctionBegin; 441 if (shell->type == MAT_COMPOSITE_ADDITIVE) { 442 ierr = MatDuplicate(next->mat,MAT_COPY_VALUES,&tmat);CHKERRQ(ierr); 443 while ((next = next->next)) { 444 ierr = MatAXPY(tmat,1.0,next->mat,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr); 445 } 446 } else { 447 if (shell->mergefromright) { 448 ierr = MatDuplicate(next->mat,MAT_COPY_VALUES,&tmat);CHKERRQ(ierr); 449 while ((next = next->next)) { 450 ierr = MatMatMult(next->mat,tmat,MAT_INITIAL_MATRIX,PETSC_DECIDE,&newmat);CHKERRQ(ierr); 451 ierr = MatDestroy(&tmat);CHKERRQ(ierr); 452 tmat = newmat; 453 } 454 } else { 455 ierr = MatDuplicate(prev->mat,MAT_COPY_VALUES,&tmat);CHKERRQ(ierr); 456 while ((prev = prev->prev)) { 457 ierr = MatMatMult(tmat,prev->mat,MAT_INITIAL_MATRIX,PETSC_DECIDE,&newmat);CHKERRQ(ierr); 458 ierr = MatDestroy(&tmat);CHKERRQ(ierr); 459 tmat = newmat; 460 } 461 } 462 } 463 464 scale = shell->scale; 465 if ((left = shell->left)) {ierr = PetscObjectReference((PetscObject)left);CHKERRQ(ierr);} 466 if ((right = shell->right)) {ierr = PetscObjectReference((PetscObject)right);CHKERRQ(ierr);} 467 468 ierr = MatHeaderReplace(mat,&tmat);CHKERRQ(ierr); 469 470 ierr = MatDiagonalScale(mat,left,right);CHKERRQ(ierr); 471 ierr = MatScale(mat,scale);CHKERRQ(ierr); 472 ierr = VecDestroy(&left);CHKERRQ(ierr); 473 ierr = VecDestroy(&right);CHKERRQ(ierr); 474 PetscFunctionReturn(0); 475 } 476 477 /*@ 478 MatCompositeMerge - Given a composite matrix, replaces it with a "regular" matrix 479 by summing all the matrices inside the composite matrix. 480 481 Collective on MPI_Comm 482 483 Input Parameters: 484 . mat - the composite matrix 485 486 487 Options Database: 488 . -mat_composite_merge (you must call MatAssemblyBegin()/MatAssemblyEnd() to have this checked) 489 490 Level: advanced 491 492 Notes: 493 The MatType of the resulting matrix will be the same as the MatType of the FIRST 494 matrix in the composite matrix. 495 496 .seealso: MatDestroy(), MatMult(), MatCompositeAddMat(), MatCreateComposite(), MatCompositeSetMergeFromRight(), MATCOMPOSITE 497 498 @*/ 499 PetscErrorCode MatCompositeMerge(Mat mat) 500 { 501 PetscErrorCode ierr; 502 503 PetscFunctionBegin; 504 PetscValidHeaderSpecific(mat,MAT_CLASSID,1); 505 ierr = PetscUseMethod(mat,"MatCompositeMerge_C",(Mat),(mat));CHKERRQ(ierr); 506 PetscFunctionReturn(0); 507 } 508 509 static PetscErrorCode MatCompositeGetNMat_Composite(Mat mat,PetscInt *nmat) 510 { 511 Mat_Composite *shell = (Mat_Composite*)mat->data; 512 513 PetscFunctionBegin; 514 *nmat = shell->nmat; 515 PetscFunctionReturn(0); 516 } 517 518 /*@ 519 MatCompositeGetNMat - Returns the number of matrices in composite. 520 521 Not Collective 522 523 Input Parameter: 524 . mat - the composite matrix 525 526 Output Parameter: 527 + size - the local size 528 - nmat - number of matrices in composite 529 530 Level: advanced 531 532 .seealso: MatCreateComposite(), MatCompositeGetMat() 533 534 @*/ 535 PetscErrorCode MatCompositeGetNMat(Mat mat,PetscInt *nmat) 536 { 537 PetscErrorCode ierr; 538 539 PetscFunctionBegin; 540 PetscValidHeaderSpecific(mat,MAT_CLASSID,1); 541 PetscValidPointer(nmat,2); 542 ierr = PetscUseMethod(mat,"MatCompositeGetNMat_C",(Mat,PetscInt*),(mat,nmat));CHKERRQ(ierr); 543 PetscFunctionReturn(0); 544 } 545 546 static PetscErrorCode MatCompositeGetMat_Composite(Mat mat,PetscInt i,Mat *Ai) 547 { 548 Mat_Composite *shell = (Mat_Composite*)mat->data; 549 Mat_CompositeLink ilink; 550 PetscInt k; 551 552 PetscFunctionBegin; 553 if (i >= shell->nmat) SETERRQ2(PetscObjectComm((PetscObject)mat),PETSC_ERR_ARG_OUTOFRANGE,"index out of range: %d >= %d",i,shell->nmat); 554 ilink = shell->head; 555 for (k=0; k<i; k++) { 556 ilink = ilink->next; 557 } 558 *Ai = ilink->mat; 559 PetscFunctionReturn(0); 560 } 561 562 /*@ 563 MatCompositeGetMat - Returns the ith matrix from composite. 564 565 Logically Collective on Mat 566 567 Input Parameter: 568 + mat - the composite matrix 569 - i - the number of requested matrix 570 571 Output Parameter: 572 . Ai - ith matrix in composite 573 574 Level: advanced 575 576 .seealso: MatCreateComposite(), MatCompositeGetNMat() 577 578 @*/ 579 PetscErrorCode MatCompositeGetMat(Mat mat,PetscInt i,Mat *Ai) 580 { 581 PetscErrorCode ierr; 582 583 PetscFunctionBegin; 584 PetscValidHeaderSpecific(mat,MAT_CLASSID,1); 585 PetscValidLogicalCollectiveInt(mat,i,2); 586 PetscValidPointer(Ai,3); 587 ierr = PetscUseMethod(mat,"MatCompositeGetMat_C",(Mat,PetscInt,Mat*),(mat,i,Ai));CHKERRQ(ierr); 588 PetscFunctionReturn(0); 589 } 590 591 static PetscErrorCode MatCompositeSetMergeFromRight_Composite(Mat mat,PetscBool flg) 592 { 593 Mat_Composite *shell = (Mat_Composite*)mat->data; 594 595 PetscFunctionBegin; 596 shell->mergefromright = flg; 597 PetscFunctionReturn(0); 598 } 599 600 /*@ 601 MatCompositeSetMergeFromRight - Sets MatCompositeMerge() to start from right. 602 603 Logically Collective on Mat 604 605 Input Parameter: 606 + mat - the composite matrix 607 - flg - if true (default) the merge starts with the first added matrix (mat[0]) 608 609 Level: advanced 610 611 .seealso: MatCreateComposite(), MatCompositeMerge() 612 613 @*/ 614 PetscErrorCode MatCompositeSetMergeFromRight(Mat mat,PetscBool flg) 615 { 616 PetscErrorCode ierr; 617 618 PetscFunctionBegin; 619 PetscValidHeaderSpecific(mat,MAT_CLASSID,1); 620 PetscValidLogicalCollectiveBool(mat,flg,2); 621 ierr = PetscUseMethod(mat,"MatCompositeSetMergeFromRight_C",(Mat,PetscBool),(mat,flg));CHKERRQ(ierr); 622 PetscFunctionReturn(0); 623 } 624 625 static struct _MatOps MatOps_Values = {0, 626 0, 627 0, 628 MatMult_Composite, 629 0, 630 /* 5*/ MatMultTranspose_Composite, 631 0, 632 0, 633 0, 634 0, 635 /* 10*/ 0, 636 0, 637 0, 638 0, 639 0, 640 /* 15*/ 0, 641 0, 642 MatGetDiagonal_Composite, 643 MatDiagonalScale_Composite, 644 0, 645 /* 20*/ 0, 646 MatAssemblyEnd_Composite, 647 0, 648 0, 649 /* 24*/ 0, 650 0, 651 0, 652 0, 653 0, 654 /* 29*/ 0, 655 0, 656 0, 657 0, 658 0, 659 /* 34*/ 0, 660 0, 661 0, 662 0, 663 0, 664 /* 39*/ 0, 665 0, 666 0, 667 0, 668 0, 669 /* 44*/ 0, 670 MatScale_Composite, 671 MatShift_Basic, 672 0, 673 0, 674 /* 49*/ 0, 675 0, 676 0, 677 0, 678 0, 679 /* 54*/ 0, 680 0, 681 0, 682 0, 683 0, 684 /* 59*/ 0, 685 MatDestroy_Composite, 686 0, 687 0, 688 0, 689 /* 64*/ 0, 690 0, 691 0, 692 0, 693 0, 694 /* 69*/ 0, 695 0, 696 0, 697 0, 698 0, 699 /* 74*/ 0, 700 0, 701 0, 702 0, 703 0, 704 /* 79*/ 0, 705 0, 706 0, 707 0, 708 0, 709 /* 84*/ 0, 710 0, 711 0, 712 0, 713 0, 714 /* 89*/ 0, 715 0, 716 0, 717 0, 718 0, 719 /* 94*/ 0, 720 0, 721 0, 722 0, 723 0, 724 /*99*/ 0, 725 0, 726 0, 727 0, 728 0, 729 /*104*/ 0, 730 0, 731 0, 732 0, 733 0, 734 /*109*/ 0, 735 0, 736 0, 737 0, 738 0, 739 /*114*/ 0, 740 0, 741 0, 742 0, 743 0, 744 /*119*/ 0, 745 0, 746 0, 747 0, 748 0, 749 /*124*/ 0, 750 0, 751 0, 752 0, 753 0, 754 /*129*/ 0, 755 0, 756 0, 757 0, 758 0, 759 /*134*/ 0, 760 0, 761 0, 762 0, 763 0, 764 /*139*/ 0, 765 0, 766 0 767 }; 768 769 /*MC 770 MATCOMPOSITE - A matrix defined by the sum (or product) of one or more matrices. 771 The matrices need to have a correct size and parallel layout for the sum or product to be valid. 772 773 Notes: 774 to use the product of the matrices call MatCompositeSetType(mat,MAT_COMPOSITE_MULTIPLICATIVE); 775 776 Level: advanced 777 778 .seealso: MatCreateComposite(), MatCompositeAddMat(), MatSetType(), MatCompositeSetType(), MatCompositeGetType(), MatCompositeMerge(), MatCompositeSetMergeFromRight(), MatCompositeGetNmat(), MatCompositeGetMat() 779 M*/ 780 781 PETSC_EXTERN PetscErrorCode MatCreate_Composite(Mat A) 782 { 783 Mat_Composite *b; 784 PetscErrorCode ierr; 785 786 PetscFunctionBegin; 787 ierr = PetscNewLog(A,&b);CHKERRQ(ierr); 788 A->data = (void*)b; 789 ierr = PetscMemcpy(A->ops,&MatOps_Values,sizeof(struct _MatOps));CHKERRQ(ierr); 790 791 ierr = PetscLayoutSetUp(A->rmap);CHKERRQ(ierr); 792 ierr = PetscLayoutSetUp(A->cmap);CHKERRQ(ierr); 793 794 A->assembled = PETSC_TRUE; 795 A->preallocated = PETSC_TRUE; 796 b->type = MAT_COMPOSITE_ADDITIVE; 797 b->scale = 1.0; 798 b->nmat = 0; 799 b->mergefromright = PETSC_TRUE; 800 801 ierr = PetscObjectComposeFunction((PetscObject)A,"MatCompositeAddMat_C",MatCompositeAddMat_Composite);CHKERRQ(ierr); 802 ierr = PetscObjectComposeFunction((PetscObject)A,"MatCompositeSetType_C",MatCompositeSetType_Composite);CHKERRQ(ierr); 803 ierr = PetscObjectComposeFunction((PetscObject)A,"MatCompositeGetType_C",MatCompositeGetType_Composite);CHKERRQ(ierr); 804 ierr = PetscObjectComposeFunction((PetscObject)A,"MatCompositeMerge_C",MatCompositeMerge_Composite);CHKERRQ(ierr); 805 ierr = PetscObjectComposeFunction((PetscObject)A,"MatCompositeGetNMat_C",MatCompositeGetNMat_Composite);CHKERRQ(ierr); 806 ierr = PetscObjectComposeFunction((PetscObject)A,"MatCompositeGetMat_C",MatCompositeGetMat_Composite);CHKERRQ(ierr); 807 ierr = PetscObjectComposeFunction((PetscObject)A,"MatCompositeSetMergeFromRight_C",MatCompositeSetMergeFromRight_Composite);CHKERRQ(ierr); 808 809 ierr = PetscObjectChangeTypeName((PetscObject)A,MATCOMPOSITE);CHKERRQ(ierr); 810 PetscFunctionReturn(0); 811 } 812 813