1 #ifdef PETSC_RCS_HEADER 2 static char vcid[] = "$Id: matrix.c,v 1.316 1999/01/27 19:47:00 bsmith Exp curfman $"; 3 #endif 4 5 /* 6 This is where the abstract matrix operations are defined 7 */ 8 9 #include "src/mat/matimpl.h" /*I "mat.h" I*/ 10 #include "src/vec/vecimpl.h" 11 12 #undef __FUNC__ 13 #define __FUNC__ "MatGetRow" 14 /*@C 15 MatGetRow - Gets a row of a matrix. You MUST call MatRestoreRow() 16 for each row that you get to ensure that your application does 17 not bleed memory. 18 19 Not Collective 20 21 Input Parameters: 22 + mat - the matrix 23 - row - the row to get 24 25 Output Parameters: 26 + ncols - the number of nonzeros in the row 27 . cols - if nonzero, the column numbers 28 - vals - if nonzero, the values 29 30 Notes: 31 This routine is provided for people who need to have direct access 32 to the structure of a matrix. We hope that we provide enough 33 high-level matrix routines that few users will need it. 34 35 MatGetRow() always returns 0-based column indices, regardless of 36 whether the internal representation is 0-based (default) or 1-based. 37 38 For better efficiency, set cols and/or vals to PETSC_NULL if you do 39 not wish to extract these quantities. 40 41 The user can only examine the values extracted with MatGetRow(); 42 the values cannot be altered. To change the matrix entries, one 43 must use MatSetValues(). 44 45 You can only have one call to MatGetRow() outstanding for a particular 46 matrix at a time. 47 48 Fortran Notes: 49 The calling sequence from Fortran is 50 .vb 51 MatGetRow(matrix,row,ncols,cols,values,ierr) 52 Mat matrix (input) 53 integer row (input) 54 integer ncols (output) 55 integer cols(maxcols) (output) 56 double precision (or double complex) values(maxcols) output 57 .ve 58 where maxcols >= maximum nonzeros in any row of the matrix. 59 60 Caution: 61 Do not try to change the contents of the output arrays (cols and vals). 62 In some cases, this may corrupt the matrix. 63 64 .keywords: matrix, row, get, extract 65 66 .seealso: MatRestoreRow(), MatSetValues() 67 @*/ 68 int MatGetRow(Mat mat,int row,int *ncols,int **cols,Scalar **vals) 69 { 70 int ierr; 71 72 PetscFunctionBegin; 73 PetscValidHeaderSpecific(mat,MAT_COOKIE); 74 PetscValidIntPointer(ncols); 75 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 76 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 77 if (!mat->ops->getrow) SETERRQ(PETSC_ERR_SUP,0,""); 78 PLogEventBegin(MAT_GetRow,mat,0,0,0); 79 ierr = (*mat->ops->getrow)(mat,row,ncols,cols,vals); CHKERRQ(ierr); 80 PLogEventEnd(MAT_GetRow,mat,0,0,0); 81 PetscFunctionReturn(0); 82 } 83 84 #undef __FUNC__ 85 #define __FUNC__ "MatRestoreRow" 86 /*@C 87 MatRestoreRow - Frees any temporary space allocated by MatGetRow(). 88 89 Not Collective 90 91 Input Parameters: 92 + mat - the matrix 93 . row - the row to get 94 . ncols, cols - the number of nonzeros and their columns 95 - vals - if nonzero the column values 96 97 Notes: 98 This routine should be called after you have finished examining the entries. 99 100 Fortran Notes: 101 The calling sequence from Fortran is 102 .vb 103 MatRestoreRow(matrix,row,ncols,cols,values,ierr) 104 Mat matrix (input) 105 integer row (input) 106 integer ncols (output) 107 integer cols(maxcols) (output) 108 double precision (or double complex) values(maxcols) output 109 .ve 110 Where maxcols >= maximum nonzeros in any row of the matrix. 111 112 In Fortran MatRestoreRow() MUST be called after MatGetRow() 113 before another call to MatGetRow() can be made. 114 115 .keywords: matrix, row, restore 116 117 .seealso: MatGetRow() 118 @*/ 119 int MatRestoreRow(Mat mat,int row,int *ncols,int **cols,Scalar **vals) 120 { 121 int ierr; 122 123 PetscFunctionBegin; 124 PetscValidHeaderSpecific(mat,MAT_COOKIE); 125 PetscValidIntPointer(ncols); 126 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 127 if (!mat->ops->restorerow) PetscFunctionReturn(0); 128 ierr = (*mat->ops->restorerow)(mat,row,ncols,cols,vals);CHKERRQ(ierr); 129 PetscFunctionReturn(0); 130 } 131 132 #undef __FUNC__ 133 #define __FUNC__ "MatView" 134 /*@C 135 MatView - Visualizes a matrix object. 136 137 Collective on Mat unless Viewer is VIEWER_STDOUT_SELF 138 139 Input Parameters: 140 + mat - the matrix 141 - ptr - visualization context 142 143 Notes: 144 The available visualization contexts include 145 + VIEWER_STDOUT_SELF - standard output (default) 146 . VIEWER_STDOUT_WORLD - synchronized standard 147 output where only the first processor opens 148 the file. All other processors send their 149 data to the first processor to print. 150 - VIEWER_DRAW_WORLD - graphical display of nonzero structure 151 152 The user can open alternative visualization contexts with 153 + ViewerASCIIOpen() - Outputs matrix to a specified file 154 . ViewerBinaryOpen() - Outputs matrix in binary to a 155 specified file; corresponding input uses MatLoad() 156 . ViewerDrawOpen() - Outputs nonzero matrix structure to 157 an X window display 158 - ViewerSocketOpen() - Outputs matrix to Socket viewer. 159 Currently only the sequential dense and AIJ 160 matrix types support the Socket viewer. 161 162 The user can call ViewerSetFormat() to specify the output 163 format of ASCII printed objects (when using VIEWER_STDOUT_SELF, 164 VIEWER_STDOUT_WORLD and ViewerASCIIOpen). Available formats include 165 + VIEWER_FORMAT_ASCII_DEFAULT - default, prints matrix contents 166 . VIEWER_FORMAT_ASCII_MATLAB - prints matrix contents in Matlab format 167 . VIEWER_FORMAT_ASCII_DENSE - prints entire matrix including zeros 168 . VIEWER_FORMAT_ASCII_COMMON - prints matrix contents, using a sparse 169 format common among all matrix types 170 . VIEWER_FORMAT_ASCII_IMPL - prints matrix contents, using an implementation-specific 171 format (which is in many cases the same as the default) 172 . VIEWER_FORMAT_ASCII_INFO - prints basic information about the matrix 173 size and structure (not the matrix entries) 174 - VIEWER_FORMAT_ASCII_INFO_LONG - prints more detailed information about 175 the matrix structure 176 177 .keywords: matrix, view, visualize, output, print, write, draw 178 179 .seealso: ViewerSetFormat(), ViewerASCIIOpen(), ViewerDrawOpen(), 180 ViewerSocketOpen(), ViewerBinaryOpen(), MatLoad() 181 @*/ 182 int MatView(Mat mat,Viewer viewer) 183 { 184 int format, ierr, rows, cols; 185 char *cstr; 186 ViewerType vtype; 187 188 PetscFunctionBegin; 189 PetscValidHeaderSpecific(mat,MAT_COOKIE); 190 if (viewer) PetscValidHeaderSpecific(viewer,VIEWER_COOKIE); 191 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 192 193 if (!viewer) { 194 viewer = VIEWER_STDOUT_SELF; 195 } 196 197 ierr = ViewerGetType(viewer,&vtype); 198 if (PetscTypeCompare(vtype,ASCII_VIEWER)) { 199 ierr = ViewerGetFormat(viewer,&format); CHKERRQ(ierr); 200 if (format == VIEWER_FORMAT_ASCII_INFO || format == VIEWER_FORMAT_ASCII_INFO_LONG) { 201 ViewerASCIIPrintf(viewer,"Matrix Object:\n"); 202 ierr = MatGetType(mat,PETSC_NULL,&cstr); CHKERRQ(ierr); 203 ierr = MatGetSize(mat,&rows,&cols); CHKERRQ(ierr); 204 ViewerASCIIPrintf(viewer," type=%s, rows=%d, cols=%d\n",cstr,rows,cols); 205 if (mat->ops->getinfo) { 206 MatInfo info; 207 ierr = MatGetInfo(mat,MAT_GLOBAL_SUM,&info); CHKERRQ(ierr); 208 ViewerASCIIPrintf(viewer," total: nonzeros=%d, allocated nonzeros=%d\n", 209 (int)info.nz_used,(int)info.nz_allocated); 210 } 211 } 212 } 213 if (mat->ops->view) { 214 ierr = ViewerASCIIPushTab(viewer);CHKERRQ(ierr); 215 ierr = (*mat->ops->view)(mat,viewer); CHKERRQ(ierr); 216 ierr = ViewerASCIIPopTab(viewer);CHKERRQ(ierr); 217 } 218 PetscFunctionReturn(0); 219 } 220 221 #undef __FUNC__ 222 #define __FUNC__ "MatScaleSystem" 223 /*@C 224 MatScaleSystem - Scale a vector solution and right hand side to 225 match the scaling of a scaled matrix. 226 227 Collective on Mat 228 229 Input Parameter: 230 + mat - the matrix 231 . x - solution vector (or PETSC_NULL) 232 + b - right hand side vector (or PETSC_NULL) 233 234 235 Notes: 236 For AIJ, BAIJ, and BDiag matrices the matrices are not 237 internally scaled, so this does nothing. For MPIRowBS it 238 permutes and diagonally scales. 239 240 The SLES methods automatically call this routine when required 241 (via PCPreSolve()) so it is rarely used directly. 242 243 Level: Developer 244 245 .keywords: matrix, scale 246 247 .seealso: MatUseScaledForm(), MatUnScaleSystem() 248 @*/ 249 int MatScaleSystem(Mat mat,Vec x,Vec b) 250 { 251 int ierr; 252 253 PetscFunctionBegin; 254 PetscValidHeaderSpecific(mat,MAT_COOKIE); 255 if (x) {PetscValidHeaderSpecific(x,VEC_COOKIE);} 256 if (b) {PetscValidHeaderSpecific(b,VEC_COOKIE);} 257 if (mat->ops->scalesystem) { 258 ierr = (*mat->ops->scalesystem)(mat,x,b); CHKERRQ(ierr); 259 } 260 PetscFunctionReturn(0); 261 } 262 263 #undef __FUNC__ 264 #define __FUNC__ "MatUnScaleSystem" 265 /*@C 266 MatUnScaleSystem - Unscales a vector solution and right hand side to 267 match the original scaling of a scaled matrix. 268 269 Collective on Mat 270 271 Input Parameter: 272 + mat - the matrix 273 . x - solution vector (or PETSC_NULL) 274 + b - right hand side vector (or PETSC_NULL) 275 276 277 Notes: 278 For AIJ, BAIJ, and BDiag matrices the matrices are not 279 internally scaled, so this does nothing. For MPIRowBS it 280 permutes and diagonally scales. 281 282 The SLES methods automatically call this routine when required 283 (via PCPostSolve()) so it is rarely used directly. 284 285 Level: Developer 286 287 .keywords: matrix, scale 288 289 .seealso: MatUseScaledForm(), MatScaleSystem() 290 @*/ 291 int MatUnScaleSystem(Mat mat,Vec x,Vec b) 292 { 293 int ierr; 294 295 PetscFunctionBegin; 296 PetscValidHeaderSpecific(mat,MAT_COOKIE); 297 if (x) {PetscValidHeaderSpecific(x,VEC_COOKIE);} 298 if (b) {PetscValidHeaderSpecific(b,VEC_COOKIE);} 299 if (mat->ops->unscalesystem) { 300 ierr = (*mat->ops->unscalesystem)(mat,x,b); CHKERRQ(ierr); 301 } 302 PetscFunctionReturn(0); 303 } 304 305 #undef __FUNC__ 306 #define __FUNC__ "MatUseScaledForm" 307 /*@C 308 MatUseScaledForm - For matrix storage formats that scale the 309 matrix (for example MPIRowbs matrices are diagonally scaled on 310 assembly) indicates matrix operations (MatMult() etc) are 311 applied using the scaled matrix. 312 313 Collective on Mat 314 315 Input Parameter: 316 + mat - the matrix 317 - scaled - PETSC_TRUE for applying the scaled, PETSC_FALSE for 318 applying the original matrix 319 320 Notes: 321 For scaled formats applying the original, unscaled matrix 322 will be slightly more expensive 323 324 Level: Developer 325 326 .keywords: matrix, scale 327 328 .seealso: MatScaleSystem(), MatUnScaleSystem() 329 @*/ 330 int MatUseScaledForm(Mat mat,PetscTruth scaled) 331 { 332 int ierr; 333 334 PetscFunctionBegin; 335 PetscValidHeaderSpecific(mat,MAT_COOKIE); 336 if (mat->ops->usescaledform) { 337 ierr = (*mat->ops->usescaledform)(mat,scaled); CHKERRQ(ierr); 338 } 339 PetscFunctionReturn(0); 340 } 341 342 #undef __FUNC__ 343 #define __FUNC__ "MatDestroy" 344 /*@C 345 MatDestroy - Frees space taken by a matrix. 346 347 Collective on Mat 348 349 Input Parameter: 350 . mat - the matrix 351 352 .keywords: matrix, destroy 353 @*/ 354 int MatDestroy(Mat mat) 355 { 356 int ierr; 357 358 PetscFunctionBegin; 359 PetscValidHeaderSpecific(mat,MAT_COOKIE); 360 ierr = (*mat->ops->destroy)(mat); CHKERRQ(ierr); 361 PetscFunctionReturn(0); 362 } 363 364 #undef __FUNC__ 365 #define __FUNC__ "MatValid" 366 /*@ 367 MatValid - Checks whether a matrix object is valid. 368 369 Collective on Mat 370 371 Input Parameter: 372 . m - the matrix to check 373 374 Output Parameter: 375 flg - flag indicating matrix status, either 376 PETSC_TRUE if matrix is valid, or PETSC_FALSE otherwise. 377 378 .keywords: matrix, valid 379 @*/ 380 int MatValid(Mat m,PetscTruth *flg) 381 { 382 PetscFunctionBegin; 383 PetscValidIntPointer(flg); 384 if (!m) *flg = PETSC_FALSE; 385 else if (m->cookie != MAT_COOKIE) *flg = PETSC_FALSE; 386 else *flg = PETSC_TRUE; 387 PetscFunctionReturn(0); 388 } 389 390 #undef __FUNC__ 391 #define __FUNC__ "MatSetValues" 392 /*@ 393 MatSetValues - Inserts or adds a block of values into a matrix. 394 These values may be cached, so MatAssemblyBegin() and MatAssemblyEnd() 395 MUST be called after all calls to MatSetValues() have been completed. 396 397 Not Collective 398 399 Input Parameters: 400 + mat - the matrix 401 . v - a logically two-dimensional array of values 402 . m, idxm - the number of rows and their global indices 403 . n, idxn - the number of columns and their global indices 404 - addv - either ADD_VALUES or INSERT_VALUES, where 405 ADD_VALUES adds values to any existing entries, and 406 INSERT_VALUES replaces existing entries with new values 407 408 Notes: 409 By default the values, v, are row-oriented and unsorted. 410 See MatSetOption() for other options. 411 412 Calls to MatSetValues() with the INSERT_VALUES and ADD_VALUES 413 options cannot be mixed without intervening calls to the assembly 414 routines. 415 416 MatSetValues() uses 0-based row and column numbers in Fortran 417 as well as in C. 418 419 Negative indices may be passed in idxm and idxn, these rows and columns are 420 simply ignored. This allows easily inserting element stiffness matrices 421 with homogeneous Dirchlet boundary conditions that you don't want represented 422 in the matrix. 423 424 Efficiency Alert: 425 The routine MatSetValuesBlocked() may offer much better efficiency 426 for users of block sparse formats (MATSEQBAIJ and MATMPIBAIJ). 427 428 Level: beginner 429 430 .keywords: matrix, insert, add, set, values 431 432 .seealso: MatSetOption(), MatAssemblyBegin(), MatAssemblyEnd(), MatSetValuesBlocked(), MatSetValuesLocal() 433 @*/ 434 int MatSetValues(Mat mat,int m,int *idxm,int n,int *idxn,Scalar *v,InsertMode addv) 435 { 436 int ierr; 437 438 PetscFunctionBegin; 439 if (!m || !n) PetscFunctionReturn(0); /* no values to insert */ 440 PetscValidHeaderSpecific(mat,MAT_COOKIE); 441 PetscValidIntPointer(idxm); 442 PetscValidIntPointer(idxn); 443 PetscValidScalarPointer(v); 444 if (mat->insertmode == NOT_SET_VALUES) { 445 mat->insertmode = addv; 446 } 447 #if defined(USE_PETSC_BOPT_g) 448 else if (mat->insertmode != addv) { 449 SETERRQ(PETSC_ERR_ARG_WRONGSTATE,1,"Cannot mix add values and insert values"); 450 } 451 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 452 #endif 453 454 if (mat->assembled) { 455 mat->was_assembled = PETSC_TRUE; 456 mat->assembled = PETSC_FALSE; 457 } 458 PLogEventBegin(MAT_SetValues,mat,0,0,0); 459 if (!mat->ops->setvalues) SETERRQ(PETSC_ERR_SUP,1,"Not supported for this matrix type"); 460 ierr = (*mat->ops->setvalues)(mat,m,idxm,n,idxn,v,addv);CHKERRQ(ierr); 461 PLogEventEnd(MAT_SetValues,mat,0,0,0); 462 PetscFunctionReturn(0); 463 } 464 465 #undef __FUNC__ 466 #define __FUNC__ "MatSetValuesBlocked" 467 /*@ 468 MatSetValuesBlocked - Inserts or adds a block of values into a matrix. 469 470 Not Collective 471 472 Input Parameters: 473 + mat - the matrix 474 . v - a logically two-dimensional array of values 475 . m, idxm - the number of block rows and their global block indices 476 . n, idxn - the number of block columns and their global block indices 477 - addv - either ADD_VALUES or INSERT_VALUES, where 478 ADD_VALUES adds values to any existing entries, and 479 INSERT_VALUES replaces existing entries with new values 480 481 Notes: 482 By default the values, v, are row-oriented and unsorted. So the layout of 483 v is the same as for MatSetValues(). See MatSetOption() for other options. 484 485 Calls to MatSetValuesBlocked() with the INSERT_VALUES and ADD_VALUES 486 options cannot be mixed without intervening calls to the assembly 487 routines. 488 489 MatSetValuesBlocked() uses 0-based row and column numbers in Fortran 490 as well as in C. 491 492 Negative indices may be passed in idxm and idxn, these rows and columns are 493 simply ignored. This allows easily inserting element stiffness matrices 494 with homogeneous Dirchlet boundary conditions that you don't want represented 495 in the matrix. 496 497 Each time an entry is set within a sparse matrix via MatSetValues(), 498 internal searching must be done to determine where to place the the 499 data in the matrix storage space. By instead inserting blocks of 500 entries via MatSetValuesBlocked(), the overhead of matrix assembly is 501 reduced. 502 503 Restrictions: 504 MatSetValuesBlocked() is currently supported only for the block AIJ 505 matrix format (MATSEQBAIJ and MATMPIBAIJ, which are created via 506 MatCreateSeqBAIJ() and MatCreateMPIBAIJ()). 507 508 .keywords: matrix, insert, add, set, values 509 510 .seealso: MatSetOption(), MatAssemblyBegin(), MatAssemblyEnd(), MatSetValues(), MatSetValuesBlockedLocal() 511 @*/ 512 int MatSetValuesBlocked(Mat mat,int m,int *idxm,int n,int *idxn,Scalar *v,InsertMode addv) 513 { 514 int ierr; 515 516 PetscFunctionBegin; 517 if (!m || !n) PetscFunctionReturn(0); /* no values to insert */ 518 PetscValidHeaderSpecific(mat,MAT_COOKIE); 519 PetscValidIntPointer(idxm); 520 PetscValidIntPointer(idxn); 521 PetscValidScalarPointer(v); 522 if (mat->insertmode == NOT_SET_VALUES) { 523 mat->insertmode = addv; 524 } 525 #if defined(USE_PETSC_BOPT_g) 526 else if (mat->insertmode != addv) { 527 SETERRQ(PETSC_ERR_ARG_WRONGSTATE,1,"Cannot mix add values and insert values"); 528 } 529 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 530 #endif 531 532 if (mat->assembled) { 533 mat->was_assembled = PETSC_TRUE; 534 mat->assembled = PETSC_FALSE; 535 } 536 PLogEventBegin(MAT_SetValues,mat,0,0,0); 537 if (!mat->ops->setvaluesblocked) SETERRQ(PETSC_ERR_SUP,1,"Not supported for this matrix type"); 538 ierr = (*mat->ops->setvaluesblocked)(mat,m,idxm,n,idxn,v,addv);CHKERRQ(ierr); 539 PLogEventEnd(MAT_SetValues,mat,0,0,0); 540 PetscFunctionReturn(0); 541 } 542 543 /*MC 544 MatSetValue - Set a single entry into a matrix. 545 546 Synopsis: 547 void MatSetValue(Mat m,int row,int col,Scalar value,InsertMode mode); 548 549 Not collective 550 551 Input Parameters: 552 + m - the matrix 553 . row - the row location of the entry 554 . col - the column location of the entry 555 . value - the value to insert 556 - mode - either INSERT_VALUES or ADD_VALUES 557 558 Notes: 559 For efficiency one should use MatSetValues() and set several or many 560 values simultaneously if possible. 561 562 Note that VecSetValue() does NOT return an error code (since this 563 is checked internally). 564 565 .seealso: MatSetValues() 566 M*/ 567 568 #undef __FUNC__ 569 #define __FUNC__ "MatGetValues" 570 /*@ 571 MatGetValues - Gets a block of values from a matrix. 572 573 Not Collective; currently only returns a local block 574 575 Input Parameters: 576 + mat - the matrix 577 . v - a logically two-dimensional array for storing the values 578 . m, idxm - the number of rows and their global indices 579 - n, idxn - the number of columns and their global indices 580 581 Notes: 582 The user must allocate space (m*n Scalars) for the values, v. 583 The values, v, are then returned in a row-oriented format, 584 analogous to that used by default in MatSetValues(). 585 586 MatGetValues() uses 0-based row and column numbers in 587 Fortran as well as in C. 588 589 MatGetValues() requires that the matrix has been assembled 590 with MatAssemblyBegin()/MatAssemblyEnd(). Thus, calls to 591 MatSetValues() and MatGetValues() CANNOT be made in succession 592 without intermediate matrix assembly. 593 594 .keywords: matrix, get, values 595 596 .seealso: MatGetRow(), MatGetSubmatrices(), MatSetValues() 597 @*/ 598 int MatGetValues(Mat mat,int m,int *idxm,int n,int *idxn,Scalar *v) 599 { 600 int ierr; 601 602 PetscFunctionBegin; 603 PetscValidHeaderSpecific(mat,MAT_COOKIE); 604 PetscValidIntPointer(idxm); 605 PetscValidIntPointer(idxn); 606 PetscValidScalarPointer(v); 607 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 608 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 609 if (!mat->ops->getvalues) SETERRQ(PETSC_ERR_SUP,0,""); 610 611 PLogEventBegin(MAT_GetValues,mat,0,0,0); 612 ierr = (*mat->ops->getvalues)(mat,m,idxm,n,idxn,v); CHKERRQ(ierr); 613 PLogEventEnd(MAT_GetValues,mat,0,0,0); 614 PetscFunctionReturn(0); 615 } 616 617 #undef __FUNC__ 618 #define __FUNC__ "MatSetLocalToGlobalMapping" 619 /*@ 620 MatSetLocalToGlobalMapping - Sets a local-to-global numbering for use by 621 the routine MatSetValuesLocal() to allow users to insert matrix entries 622 using a local (per-processor) numbering. 623 624 Not Collective 625 626 Input Parameters: 627 + x - the matrix 628 - mapping - mapping created with ISLocalToGlobalMappingCreate() 629 or ISLocalToGlobalMappingCreateIS() 630 631 .keywords: matrix, set, values, local, global, mapping 632 633 .seealso: MatAssemblyBegin(), MatAssemblyEnd(), MatSetValues(), MatSetValuesLocal() 634 @*/ 635 int MatSetLocalToGlobalMapping(Mat x,ISLocalToGlobalMapping mapping) 636 { 637 PetscFunctionBegin; 638 PetscValidHeaderSpecific(x,MAT_COOKIE); 639 PetscValidHeaderSpecific(mapping,IS_LTOGM_COOKIE); 640 641 if (x->mapping) { 642 SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Mapping already set for matrix"); 643 } 644 645 x->mapping = mapping; 646 PetscObjectReference((PetscObject)mapping); 647 PetscFunctionReturn(0); 648 } 649 650 #undef __FUNC__ 651 #define __FUNC__ "MatSetLocalToGlobalMappingBlocked" 652 /*@ 653 MatSetLocalToGlobalMappingBlocked - Sets a local-to-global numbering for use 654 by the routine MatSetValuesBlockedLocal() to allow users to insert matrix 655 entries using a local (per-processor) numbering. 656 657 Not Collective 658 659 Input Parameters: 660 + x - the matrix 661 - mapping - mapping created with ISLocalToGlobalMappingCreate() or 662 ISLocalToGlobalMappingCreateIS() 663 664 .keywords: matrix, set, values, local ordering 665 666 .seealso: MatAssemblyBegin(), MatAssemblyEnd(), MatSetValues(), MatSetValuesBlockedLocal(), 667 MatSetValuesBlocked(), MatSetValuesLocal() 668 @*/ 669 int MatSetLocalToGlobalMappingBlocked(Mat x,ISLocalToGlobalMapping mapping) 670 { 671 PetscFunctionBegin; 672 PetscValidHeaderSpecific(x,MAT_COOKIE); 673 PetscValidHeaderSpecific(mapping,IS_LTOGM_COOKIE); 674 675 if (x->bmapping) { 676 SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Mapping already set for matrix"); 677 } 678 679 x->bmapping = mapping; 680 PetscObjectReference((PetscObject)mapping); 681 PetscFunctionReturn(0); 682 } 683 684 #undef __FUNC__ 685 #define __FUNC__ "MatSetValuesLocal" 686 /*@ 687 MatSetValuesLocal - Inserts or adds values into certain locations of a matrix, 688 using a local ordering of the nodes. 689 690 Not Collective 691 692 Input Parameters: 693 + x - the matrix 694 . nrow, irow - number of rows and their local indices 695 . ncol, icol - number of columns and their local indices 696 . y - a logically two-dimensional array of values 697 - addv - either INSERT_VALUES or ADD_VALUES, where 698 ADD_VALUES adds values to any existing entries, and 699 INSERT_VALUES replaces existing entries with new values 700 701 Notes: 702 Before calling MatSetValuesLocal(), the user must first set the 703 local-to-global mapping by calling MatSetLocalToGlobalMapping(). 704 705 Calls to MatSetValuesLocal() with the INSERT_VALUES and ADD_VALUES 706 options cannot be mixed without intervening calls to the assembly 707 routines. 708 709 These values may be cached, so MatAssemblyBegin() and MatAssemblyEnd() 710 MUST be called after all calls to MatSetValuesLocal() have been completed. 711 712 .keywords: matrix, set, values, local ordering 713 714 .seealso: MatAssemblyBegin(), MatAssemblyEnd(), MatSetValues(), MatSetLocalToGlobalMapping() 715 @*/ 716 int MatSetValuesLocal(Mat mat,int nrow,int *irow,int ncol, int *icol,Scalar *y,InsertMode addv) 717 { 718 int ierr,irowm[2048],icolm[2048]; 719 720 PetscFunctionBegin; 721 PetscValidHeaderSpecific(mat,MAT_COOKIE); 722 PetscValidIntPointer(irow); 723 PetscValidIntPointer(icol); 724 PetscValidScalarPointer(y); 725 726 if (mat->insertmode == NOT_SET_VALUES) { 727 mat->insertmode = addv; 728 } 729 #if defined(USE_PETSC_BOPT_g) 730 else if (mat->insertmode != addv) { 731 SETERRQ(PETSC_ERR_ARG_WRONGSTATE,1,"Cannot mix add values and insert values"); 732 } 733 if (!mat->mapping) { 734 SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Local to global never set with MatSetLocalToGlobalMapping()"); 735 } 736 if (nrow > 2048 || ncol > 2048) { 737 SETERRQ2(PETSC_ERR_SUP,0,"Number column/row indices must be <= 2048: are %d %d",nrow,ncol); 738 } 739 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 740 #endif 741 742 if (mat->assembled) { 743 mat->was_assembled = PETSC_TRUE; 744 mat->assembled = PETSC_FALSE; 745 } 746 PLogEventBegin(MAT_SetValues,mat,0,0,0); 747 ierr = ISLocalToGlobalMappingApply(mat->mapping,nrow,irow,irowm); CHKERRQ(ierr); 748 ierr = ISLocalToGlobalMappingApply(mat->mapping,ncol,icol,icolm);CHKERRQ(ierr); 749 ierr = (*mat->ops->setvalues)(mat,nrow,irowm,ncol,icolm,y,addv);CHKERRQ(ierr); 750 PLogEventEnd(MAT_SetValues,mat,0,0,0); 751 PetscFunctionReturn(0); 752 } 753 754 #undef __FUNC__ 755 #define __FUNC__ "MatSetValuesBlockedLocal" 756 /*@ 757 MatSetValuesBlockedLocal - Inserts or adds values into certain locations of a matrix, 758 using a local ordering of the nodes a block at a time. 759 760 Not Collective 761 762 Input Parameters: 763 + x - the matrix 764 . nrow, irow - number of rows and their local indices 765 . ncol, icol - number of columns and their local indices 766 . y - a logically two-dimensional array of values 767 - addv - either INSERT_VALUES or ADD_VALUES, where 768 ADD_VALUES adds values to any existing entries, and 769 INSERT_VALUES replaces existing entries with new values 770 771 Notes: 772 Before calling MatSetValuesBlockedLocal(), the user must first set the 773 local-to-global mapping by calling MatSetLocalToGlobalMappingBlocked(), 774 where the mapping MUST be set for matrix blocks, not for matrix elements. 775 776 Calls to MatSetValuesBlockedLocal() with the INSERT_VALUES and ADD_VALUES 777 options cannot be mixed without intervening calls to the assembly 778 routines. 779 780 These values may be cached, so MatAssemblyBegin() and MatAssemblyEnd() 781 MUST be called after all calls to MatSetValuesBlockedLocal() have been completed. 782 783 .keywords: matrix, set, values, blocked, local 784 785 .seealso: MatAssemblyBegin(), MatAssemblyEnd(), MatSetValuesLocal(), MatSetLocalToGlobalMappingBlocked(), MatSetValuesBlocked() 786 @*/ 787 int MatSetValuesBlockedLocal(Mat mat,int nrow,int *irow,int ncol,int *icol,Scalar *y,InsertMode addv) 788 { 789 int ierr,irowm[2048],icolm[2048]; 790 791 PetscFunctionBegin; 792 PetscValidHeaderSpecific(mat,MAT_COOKIE); 793 PetscValidIntPointer(irow); 794 PetscValidIntPointer(icol); 795 PetscValidScalarPointer(y); 796 if (mat->insertmode == NOT_SET_VALUES) { 797 mat->insertmode = addv; 798 } 799 #if defined(USE_PETSC_BOPT_g) 800 else if (mat->insertmode != addv) { 801 SETERRQ(PETSC_ERR_ARG_WRONGSTATE,1,"Cannot mix add values and insert values"); 802 } 803 if (!mat->bmapping) { 804 SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Local to global never set with MatSetLocalToGlobalMappingBlocked()"); 805 } 806 if (nrow > 2048 || ncol > 2048) { 807 SETERRQ2(PETSC_ERR_SUP,0,"Number column/row indices must be <= 2048: are %d %d",nrow,ncol); 808 } 809 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 810 #endif 811 812 if (mat->assembled) { 813 mat->was_assembled = PETSC_TRUE; 814 mat->assembled = PETSC_FALSE; 815 } 816 PLogEventBegin(MAT_SetValues,mat,0,0,0); 817 ierr = ISLocalToGlobalMappingApply(mat->bmapping,nrow,irow,irowm); CHKERRQ(ierr); 818 ierr = ISLocalToGlobalMappingApply(mat->bmapping,ncol,icol,icolm); CHKERRQ(ierr); 819 ierr = (*mat->ops->setvaluesblocked)(mat,nrow,irowm,ncol,icolm,y,addv);CHKERRQ(ierr); 820 PLogEventEnd(MAT_SetValues,mat,0,0,0); 821 PetscFunctionReturn(0); 822 } 823 824 /* --------------------------------------------------------*/ 825 #undef __FUNC__ 826 #define __FUNC__ "MatMult" 827 /*@ 828 MatMult - Computes the matrix-vector product, y = Ax. 829 830 Collective on Mat and Vec 831 832 Input Parameters: 833 + mat - the matrix 834 - x - the vector to be multilplied 835 836 Output Parameters: 837 . y - the result 838 839 Notes: 840 The vectors x and y cannot be the same. I.e., one cannot 841 call MatMult(A,y,y). 842 843 .keywords: matrix, multiply, matrix-vector product 844 845 .seealso: MatMultTrans(), MatMultAdd(), MatMultTransAdd() 846 @*/ 847 int MatMult(Mat mat,Vec x,Vec y) 848 { 849 int ierr; 850 851 PetscFunctionBegin; 852 PetscValidHeaderSpecific(mat,MAT_COOKIE); 853 PetscValidHeaderSpecific(x,VEC_COOKIE);PetscValidHeaderSpecific(y,VEC_COOKIE); 854 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 855 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 856 if (x == y) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"x and y must be different vectors"); 857 if (mat->N != x->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim %d %d",mat->N,x->N); 858 if (mat->M != y->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec y: global dim %d %d",mat->M,y->N); 859 if (mat->m != y->n) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec y: local dim %d %d",mat->m,y->n); 860 861 PLogEventBegin(MAT_Mult,mat,x,y,0); 862 ierr = (*mat->ops->mult)(mat,x,y); CHKERRQ(ierr); 863 PLogEventEnd(MAT_Mult,mat,x,y,0); 864 865 PetscFunctionReturn(0); 866 } 867 868 #undef __FUNC__ 869 #define __FUNC__ "MatMultTrans" 870 /*@ 871 MatMultTrans - Computes matrix transpose times a vector. 872 873 Collective on Mat and Vec 874 875 Input Parameters: 876 + mat - the matrix 877 - x - the vector to be multilplied 878 879 Output Parameters: 880 . y - the result 881 882 Notes: 883 The vectors x and y cannot be the same. I.e., one cannot 884 call MatMultTrans(A,y,y). 885 886 .keywords: matrix, multiply, matrix-vector product, transpose 887 888 .seealso: MatMult(), MatMultAdd(), MatMultTransAdd() 889 @*/ 890 int MatMultTrans(Mat mat,Vec x,Vec y) 891 { 892 int ierr; 893 894 PetscFunctionBegin; 895 PetscValidHeaderSpecific(mat,MAT_COOKIE); 896 PetscValidHeaderSpecific(x,VEC_COOKIE); PetscValidHeaderSpecific(y,VEC_COOKIE); 897 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 898 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 899 if (x == y) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"x and y must be different vectors"); 900 if (mat->M != x->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim %d %d",mat->M,x->N); 901 if (mat->N != y->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec y: global dim %d %d",mat->N,y->N); 902 PLogEventBegin(MAT_MultTrans,mat,x,y,0); 903 ierr = (*mat->ops->multtrans)(mat,x,y); CHKERRQ(ierr); 904 PLogEventEnd(MAT_MultTrans,mat,x,y,0); 905 PetscFunctionReturn(0); 906 } 907 908 #undef __FUNC__ 909 #define __FUNC__ "MatMultAdd" 910 /*@ 911 MatMultAdd - Computes v3 = v2 + A * v1. 912 913 Collective on Mat and Vec 914 915 Input Parameters: 916 + mat - the matrix 917 - v1, v2 - the vectors 918 919 Output Parameters: 920 . v3 - the result 921 922 Notes: 923 The vectors v1 and v3 cannot be the same. I.e., one cannot 924 call MatMultAdd(A,v1,v2,v1). 925 926 .keywords: matrix, multiply, matrix-vector product, add 927 928 .seealso: MatMultTrans(), MatMult(), MatMultTransAdd() 929 @*/ 930 int MatMultAdd(Mat mat,Vec v1,Vec v2,Vec v3) 931 { 932 int ierr; 933 934 PetscFunctionBegin; 935 PetscValidHeaderSpecific(mat,MAT_COOKIE);PetscValidHeaderSpecific(v1,VEC_COOKIE); 936 PetscValidHeaderSpecific(v2,VEC_COOKIE); PetscValidHeaderSpecific(v3,VEC_COOKIE); 937 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 938 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 939 if (mat->N != v1->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec v1: global dim %d %d",mat->N,v1->N); 940 if (mat->M != v2->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec v2: global dim %d %d",mat->M,v2->N); 941 if (mat->M != v3->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec v3: global dim %d %d",mat->M,v3->N); 942 if (mat->m != v3->n) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec v3: local dim %d %d",mat->m,v3->n); 943 if (mat->m != v2->n) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec v2: local dim %d %d",mat->m,v2->n); 944 if (v1 == v3) SETERRQ(PETSC_ERR_ARG_IDN,0,"v1 and v3 must be different vectors"); 945 946 PLogEventBegin(MAT_MultAdd,mat,v1,v2,v3); 947 ierr = (*mat->ops->multadd)(mat,v1,v2,v3); CHKERRQ(ierr); 948 PLogEventEnd(MAT_MultAdd,mat,v1,v2,v3); 949 PetscFunctionReturn(0); 950 } 951 952 #undef __FUNC__ 953 #define __FUNC__ "MatMultTransAdd" 954 /*@ 955 MatMultTransAdd - Computes v3 = v2 + A' * v1. 956 957 Collective on Mat and Vec 958 959 Input Parameters: 960 + mat - the matrix 961 - v1, v2 - the vectors 962 963 Output Parameters: 964 . v3 - the result 965 966 Notes: 967 The vectors v1 and v3 cannot be the same. I.e., one cannot 968 call MatMultTransAdd(A,v1,v2,v1). 969 970 .keywords: matrix, multiply, matrix-vector product, transpose, add 971 972 .seealso: MatMultTrans(), MatMultAdd(), MatMult() 973 @*/ 974 int MatMultTransAdd(Mat mat,Vec v1,Vec v2,Vec v3) 975 { 976 int ierr; 977 978 PetscFunctionBegin; 979 PetscValidHeaderSpecific(mat,MAT_COOKIE);PetscValidHeaderSpecific(v1,VEC_COOKIE); 980 PetscValidHeaderSpecific(v2,VEC_COOKIE);PetscValidHeaderSpecific(v3,VEC_COOKIE); 981 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 982 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 983 if (!mat->ops->multtransadd) SETERRQ(PETSC_ERR_SUP,0,""); 984 if (v1 == v3) SETERRQ(PETSC_ERR_ARG_IDN,0,"v1 and v3 must be different vectors"); 985 if (mat->M != v1->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec v1: global dim %d %d",mat->M,v1->N); 986 if (mat->N != v2->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec v2: global dim %d %d",mat->N,v2->N); 987 if (mat->N != v3->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec v3: global dim %d %d",mat->N,v3->N); 988 989 PLogEventBegin(MAT_MultTransAdd,mat,v1,v2,v3); 990 ierr = (*mat->ops->multtransadd)(mat,v1,v2,v3); CHKERRQ(ierr); 991 PLogEventEnd(MAT_MultTransAdd,mat,v1,v2,v3); 992 PetscFunctionReturn(0); 993 } 994 /* ------------------------------------------------------------*/ 995 #undef __FUNC__ 996 #define __FUNC__ "MatGetInfo" 997 /*@C 998 MatGetInfo - Returns information about matrix storage (number of 999 nonzeros, memory, etc.). 1000 1001 Collective on Mat if MAT_GLOBAL_MAX or MAT_GLOBAL_SUM is used 1002 as the flag 1003 1004 Input Parameters: 1005 . mat - the matrix 1006 1007 Output Parameters: 1008 + flag - flag indicating the type of parameters to be returned 1009 (MAT_LOCAL - local matrix, MAT_GLOBAL_MAX - maximum over all processors, 1010 MAT_GLOBAL_SUM - sum over all processors) 1011 - info - matrix information context 1012 1013 Notes: 1014 The MatInfo context contains a variety of matrix data, including 1015 number of nonzeros allocated and used, number of mallocs during 1016 matrix assembly, etc. Additional information for factored matrices 1017 is provided (such as the fill ratio, number of mallocs during 1018 factorization, etc.). Much of this info is printed to STDOUT 1019 when using the runtime options 1020 $ -log_info -mat_view_info 1021 1022 Example for C/C++ Users: 1023 See the file ${PETSC_DIR}/include/mat.h for a complete list of 1024 data within the MatInfo context. For example, 1025 .vb 1026 MatInfo info; 1027 Mat A; 1028 double mal, nz_a, nz_u; 1029 1030 MatGetInfo(A,MAT_LOCAL,&info); 1031 mal = info.mallocs; 1032 nz_a = info.nz_allocated; 1033 .ve 1034 1035 Example for Fortran Users: 1036 Fortran users should declare info as a double precision 1037 array of dimension MAT_INFO_SIZE, and then extract the parameters 1038 of interest. See the file ${PETSC_DIR}/include/finclude/mat.h 1039 a complete list of parameter names. 1040 .vb 1041 double precision info(MAT_INFO_SIZE) 1042 double precision mal, nz_a 1043 Mat A 1044 integer ierr 1045 1046 call MatGetInfo(A,MAT_LOCAL,info,ierr) 1047 mal = info(MAT_INFO_MALLOCS) 1048 nz_a = info(MAT_INFO_NZ_ALLOCATED) 1049 .ve 1050 1051 .keywords: matrix, get, info, storage, nonzeros, memory, fill 1052 @*/ 1053 int MatGetInfo(Mat mat,MatInfoType flag,MatInfo *info) 1054 { 1055 int ierr; 1056 1057 PetscFunctionBegin; 1058 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1059 PetscValidPointer(info); 1060 if (!mat->ops->getinfo) SETERRQ(PETSC_ERR_SUP,0,""); 1061 ierr = (*mat->ops->getinfo)(mat,flag,info);CHKERRQ(ierr); 1062 PetscFunctionReturn(0); 1063 } 1064 1065 /* ----------------------------------------------------------*/ 1066 #undef __FUNC__ 1067 #define __FUNC__ "MatILUDTFactor" 1068 /*@ 1069 MatILUDTFactor - Performs a drop tolerance ILU factorization. 1070 1071 Collective on Mat 1072 1073 Input Parameters: 1074 + mat - the matrix 1075 . dt - the drop tolerance 1076 . maxnz - the maximum number of nonzeros per row allowed? 1077 . row - row permutation 1078 - col - column permutation 1079 1080 Output Parameters: 1081 . fact - the factored matrix 1082 1083 .keywords: matrix, factor, LU, in-place 1084 1085 .seealso: MatLUFactorSymbolic(), MatLUFactorNumeric(), MatCholeskyFactor() 1086 @*/ 1087 int MatILUDTFactor(Mat mat,double dt,int maxnz,IS row,IS col,Mat *fact) 1088 { 1089 int ierr; 1090 1091 PetscFunctionBegin; 1092 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1093 PetscValidPointer(fact); 1094 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 1095 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 1096 if (!mat->ops->iludtfactor) SETERRQ(PETSC_ERR_SUP,0,""); 1097 1098 PLogEventBegin(MAT_ILUFactor,mat,row,col,0); 1099 ierr = (*mat->ops->iludtfactor)(mat,dt,maxnz,row,col,fact); CHKERRQ(ierr); 1100 PLogEventEnd(MAT_ILUFactor,mat,row,col,0); 1101 1102 PetscFunctionReturn(0); 1103 } 1104 1105 #undef __FUNC__ 1106 #define __FUNC__ "MatLUFactor" 1107 /*@ 1108 MatLUFactor - Performs in-place LU factorization of matrix. 1109 1110 Collective on Mat 1111 1112 Input Parameters: 1113 + mat - the matrix 1114 . row - row permutation 1115 . col - column permutation 1116 - f - expected fill as ratio of original fill. 1117 1118 .keywords: matrix, factor, LU, in-place 1119 1120 .seealso: MatLUFactorSymbolic(), MatLUFactorNumeric(), MatCholeskyFactor() 1121 @*/ 1122 int MatLUFactor(Mat mat,IS row,IS col,double f) 1123 { 1124 int ierr; 1125 1126 PetscFunctionBegin; 1127 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1128 if (mat->M != mat->N) SETERRQ(PETSC_ERR_ARG_WRONG,0,"matrix must be square"); 1129 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 1130 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 1131 if (!mat->ops->lufactor) SETERRQ(PETSC_ERR_SUP,0,""); 1132 1133 PLogEventBegin(MAT_LUFactor,mat,row,col,0); 1134 ierr = (*mat->ops->lufactor)(mat,row,col,f); CHKERRQ(ierr); 1135 PLogEventEnd(MAT_LUFactor,mat,row,col,0); 1136 PetscFunctionReturn(0); 1137 } 1138 1139 #undef __FUNC__ 1140 #define __FUNC__ "MatILUFactor" 1141 /*@ 1142 MatILUFactor - Performs in-place ILU factorization of matrix. 1143 1144 Collective on Mat 1145 1146 Input Parameters: 1147 + mat - the matrix 1148 . row - row permutation 1149 . col - column permutation 1150 - info - structure containing 1151 $ levels - number of levels of fill. 1152 $ expected fill - as ratio of original fill. 1153 $ 1 or 0 - indicating force fill on diagonal (improves robustness for matrices 1154 missing diagonal entries) 1155 1156 Notes: 1157 Probably really in-place only when level of fill is zero, otherwise allocates 1158 new space to store factored matrix and deletes previous memory. 1159 1160 .keywords: matrix, factor, ILU, in-place 1161 1162 .seealso: MatILUFactorSymbolic(), MatLUFactorNumeric(), MatCholeskyFactor() 1163 @*/ 1164 int MatILUFactor(Mat mat,IS row,IS col,MatILUInfo *info) 1165 { 1166 int ierr; 1167 1168 PetscFunctionBegin; 1169 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1170 if (mat->M != mat->N) SETERRQ(PETSC_ERR_ARG_WRONG,0,"matrix must be square"); 1171 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 1172 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 1173 if (!mat->ops->ilufactor) SETERRQ(PETSC_ERR_SUP,0,""); 1174 1175 PLogEventBegin(MAT_ILUFactor,mat,row,col,0); 1176 ierr = (*mat->ops->ilufactor)(mat,row,col,info); CHKERRQ(ierr); 1177 PLogEventEnd(MAT_ILUFactor,mat,row,col,0); 1178 PetscFunctionReturn(0); 1179 } 1180 1181 #undef __FUNC__ 1182 #define __FUNC__ "MatLUFactorSymbolic" 1183 /*@ 1184 MatLUFactorSymbolic - Performs symbolic LU factorization of matrix. 1185 Call this routine before calling MatLUFactorNumeric(). 1186 1187 Collective on Mat 1188 1189 Input Parameters: 1190 + mat - the matrix 1191 . row, col - row and column permutations 1192 - f - expected fill as ratio of the original number of nonzeros, 1193 for example 3.0; choosing this parameter well can result in 1194 more efficient use of time and space. Run with the option -log_info 1195 to determine an optimal value to use 1196 1197 Output Parameter: 1198 . fact - new matrix that has been symbolically factored 1199 1200 Notes: 1201 See the file ${PETSC_DIR}/Performance for additional information about 1202 choosing the fill factor for better efficiency. 1203 1204 .keywords: matrix, factor, LU, symbolic, fill 1205 1206 .seealso: MatLUFactor(), MatLUFactorNumeric(), MatCholeskyFactor() 1207 @*/ 1208 int MatLUFactorSymbolic(Mat mat,IS row,IS col,double f,Mat *fact) 1209 { 1210 int ierr; 1211 1212 PetscFunctionBegin; 1213 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1214 if (mat->M != mat->N) SETERRQ(PETSC_ERR_ARG_WRONG,0,"matrix must be square"); 1215 PetscValidPointer(fact); 1216 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 1217 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 1218 if (!mat->ops->lufactorsymbolic) SETERRQ(PETSC_ERR_SUP,0,""); 1219 1220 PLogEventBegin(MAT_LUFactorSymbolic,mat,row,col,0); 1221 ierr = (*mat->ops->lufactorsymbolic)(mat,row,col,f,fact); CHKERRQ(ierr); 1222 PLogEventEnd(MAT_LUFactorSymbolic,mat,row,col,0); 1223 PetscFunctionReturn(0); 1224 } 1225 1226 #undef __FUNC__ 1227 #define __FUNC__ "MatLUFactorNumeric" 1228 /*@ 1229 MatLUFactorNumeric - Performs numeric LU factorization of a matrix. 1230 Call this routine after first calling MatLUFactorSymbolic(). 1231 1232 Collective on Mat 1233 1234 Input Parameters: 1235 + mat - the matrix 1236 - row, col - row and column permutations 1237 1238 Output Parameters: 1239 . fact - symbolically factored matrix that must have been generated 1240 by MatLUFactorSymbolic() 1241 1242 Notes: 1243 See MatLUFactor() for in-place factorization. See 1244 MatCholeskyFactorNumeric() for the symmetric, positive definite case. 1245 1246 .keywords: matrix, factor, LU, numeric 1247 1248 .seealso: MatLUFactorSymbolic(), MatLUFactor(), MatCholeskyFactor() 1249 @*/ 1250 int MatLUFactorNumeric(Mat mat,Mat *fact) 1251 { 1252 int ierr,flg; 1253 1254 PetscFunctionBegin; 1255 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1256 PetscValidPointer(fact); 1257 PetscValidHeaderSpecific(*fact,MAT_COOKIE); 1258 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 1259 if (mat->M != (*fact)->M || mat->N != (*fact)->N) { 1260 SETERRQ4(PETSC_ERR_ARG_SIZ,0,"Mat mat,Mat *fact: global dimensions are different %d should = %d %d should = %d", 1261 mat->M,(*fact)->M,mat->N,(*fact)->N); 1262 } 1263 if (!(*fact)->ops->lufactornumeric) SETERRQ(PETSC_ERR_SUP,0,""); 1264 1265 PLogEventBegin(MAT_LUFactorNumeric,mat,*fact,0,0); 1266 ierr = (*(*fact)->ops->lufactornumeric)(mat,fact); CHKERRQ(ierr); 1267 PLogEventEnd(MAT_LUFactorNumeric,mat,*fact,0,0); 1268 ierr = OptionsHasName(PETSC_NULL,"-mat_view_draw",&flg); CHKERRQ(ierr); 1269 if (flg) { 1270 ierr = OptionsHasName(0,"-mat_view_contour",&flg); CHKERRQ(ierr); 1271 if (flg) { 1272 ViewerPushFormat(VIEWER_DRAW_(mat->comm),VIEWER_FORMAT_DRAW_CONTOUR,0);CHKERRQ(ierr); 1273 } 1274 ierr = MatView(*fact,VIEWER_DRAW_(mat->comm)); CHKERRQ(ierr); 1275 ierr = ViewerFlush(VIEWER_DRAW_(mat->comm)); CHKERRQ(ierr); 1276 if (flg) { 1277 ViewerPopFormat(VIEWER_DRAW_(mat->comm));CHKERRQ(ierr); 1278 } 1279 } 1280 PetscFunctionReturn(0); 1281 } 1282 1283 #undef __FUNC__ 1284 #define __FUNC__ "MatCholeskyFactor" 1285 /*@ 1286 MatCholeskyFactor - Performs in-place Cholesky factorization of a 1287 symmetric matrix. 1288 1289 Collective on Mat 1290 1291 Input Parameters: 1292 + mat - the matrix 1293 . perm - row and column permutations 1294 - f - expected fill as ratio of original fill 1295 1296 Notes: 1297 See MatLUFactor() for the nonsymmetric case. See also 1298 MatCholeskyFactorSymbolic(), and MatCholeskyFactorNumeric(). 1299 1300 .keywords: matrix, factor, in-place, Cholesky 1301 1302 .seealso: MatLUFactor(), MatCholeskyFactorSymbolic(), MatCholeskyFactorNumeric() 1303 @*/ 1304 int MatCholeskyFactor(Mat mat,IS perm,double f) 1305 { 1306 int ierr; 1307 1308 PetscFunctionBegin; 1309 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1310 if (mat->M != mat->N) SETERRQ(PETSC_ERR_ARG_WRONG,0,"Matrix must be square"); 1311 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 1312 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 1313 if (!mat->ops->choleskyfactor) SETERRQ(PETSC_ERR_SUP,0,""); 1314 1315 PLogEventBegin(MAT_CholeskyFactor,mat,perm,0,0); 1316 ierr = (*mat->ops->choleskyfactor)(mat,perm,f); CHKERRQ(ierr); 1317 PLogEventEnd(MAT_CholeskyFactor,mat,perm,0,0); 1318 PetscFunctionReturn(0); 1319 } 1320 1321 #undef __FUNC__ 1322 #define __FUNC__ "MatCholeskyFactorSymbolic" 1323 /*@ 1324 MatCholeskyFactorSymbolic - Performs symbolic Cholesky factorization 1325 of a symmetric matrix. 1326 1327 Collective on Mat 1328 1329 Input Parameters: 1330 + mat - the matrix 1331 . perm - row and column permutations 1332 - f - expected fill as ratio of original 1333 1334 Output Parameter: 1335 . fact - the factored matrix 1336 1337 Notes: 1338 See MatLUFactorSymbolic() for the nonsymmetric case. See also 1339 MatCholeskyFactor() and MatCholeskyFactorNumeric(). 1340 1341 .keywords: matrix, factor, factorization, symbolic, Cholesky 1342 1343 .seealso: MatLUFactorSymbolic(), MatCholeskyFactor(), MatCholeskyFactorNumeric() 1344 @*/ 1345 int MatCholeskyFactorSymbolic(Mat mat,IS perm,double f,Mat *fact) 1346 { 1347 int ierr; 1348 1349 PetscFunctionBegin; 1350 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1351 PetscValidPointer(fact); 1352 if (mat->M != mat->N) SETERRQ(PETSC_ERR_ARG_WRONG,0,"Matrix must be square"); 1353 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 1354 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 1355 if (!mat->ops->choleskyfactorsymbolic) SETERRQ(PETSC_ERR_SUP,0,""); 1356 1357 PLogEventBegin(MAT_CholeskyFactorSymbolic,mat,perm,0,0); 1358 ierr = (*mat->ops->choleskyfactorsymbolic)(mat,perm,f,fact); CHKERRQ(ierr); 1359 PLogEventEnd(MAT_CholeskyFactorSymbolic,mat,perm,0,0); 1360 PetscFunctionReturn(0); 1361 } 1362 1363 #undef __FUNC__ 1364 #define __FUNC__ "MatCholeskyFactorNumeric" 1365 /*@ 1366 MatCholeskyFactorNumeric - Performs numeric Cholesky factorization 1367 of a symmetric matrix. Call this routine after first calling 1368 MatCholeskyFactorSymbolic(). 1369 1370 Collective on Mat 1371 1372 Input Parameter: 1373 . mat - the initial matrix 1374 1375 Output Parameter: 1376 . fact - the factored matrix 1377 1378 .keywords: matrix, factor, numeric, Cholesky 1379 1380 .seealso: MatCholeskyFactorSymbolic(), MatCholeskyFactor(), MatLUFactorNumeric() 1381 @*/ 1382 int MatCholeskyFactorNumeric(Mat mat,Mat *fact) 1383 { 1384 int ierr; 1385 1386 PetscFunctionBegin; 1387 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1388 PetscValidPointer(fact); 1389 if (!mat->ops->choleskyfactornumeric) SETERRQ(PETSC_ERR_SUP,0,""); 1390 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 1391 if (mat->M != (*fact)->M || mat->N != (*fact)->N) { 1392 SETERRQ4(PETSC_ERR_ARG_SIZ,0,"Mat mat,Mat *fact: global dim %d should = %d %d should = %d", 1393 mat->M,(*fact)->M,mat->N,(*fact)->N); 1394 } 1395 1396 PLogEventBegin(MAT_CholeskyFactorNumeric,mat,*fact,0,0); 1397 ierr = (*mat->ops->choleskyfactornumeric)(mat,fact); CHKERRQ(ierr); 1398 PLogEventEnd(MAT_CholeskyFactorNumeric,mat,*fact,0,0); 1399 PetscFunctionReturn(0); 1400 } 1401 1402 /* ----------------------------------------------------------------*/ 1403 #undef __FUNC__ 1404 #define __FUNC__ "MatSolve" 1405 /*@ 1406 MatSolve - Solves A x = b, given a factored matrix. 1407 1408 Collective on Mat and Vec 1409 1410 Input Parameters: 1411 + mat - the factored matrix 1412 - b - the right-hand-side vector 1413 1414 Output Parameter: 1415 . x - the result vector 1416 1417 Notes: 1418 The vectors b and x cannot be the same. I.e., one cannot 1419 call MatSolve(A,x,x). 1420 1421 .keywords: matrix, linear system, solve, LU, Cholesky, triangular solve 1422 1423 .seealso: MatSolveAdd(), MatSolveTrans(), MatSolveTransAdd() 1424 @*/ 1425 int MatSolve(Mat mat,Vec b,Vec x) 1426 { 1427 int ierr; 1428 1429 PetscFunctionBegin; 1430 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1431 PetscValidHeaderSpecific(b,VEC_COOKIE); PetscValidHeaderSpecific(x,VEC_COOKIE); 1432 if (x == b) SETERRQ(PETSC_ERR_ARG_IDN,0,"x and b must be different vectors"); 1433 if (!mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Unfactored matrix"); 1434 if (mat->N != x->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim %d %d",mat->N,x->N); 1435 if (mat->M != b->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: global dim %d %d",mat->M,b->N); 1436 if (mat->m != b->n) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: local dim %d %d",mat->m,b->n); 1437 if (mat->M == 0 && mat->N == 0) PetscFunctionReturn(0); 1438 1439 if (!mat->ops->solve) SETERRQ(PETSC_ERR_SUP,0,""); 1440 PLogEventBegin(MAT_Solve,mat,b,x,0); 1441 ierr = (*mat->ops->solve)(mat,b,x); CHKERRQ(ierr); 1442 PLogEventEnd(MAT_Solve,mat,b,x,0); 1443 PetscFunctionReturn(0); 1444 } 1445 1446 #undef __FUNC__ 1447 #define __FUNC__ "MatForwardSolve" 1448 /* @ 1449 MatForwardSolve - Solves L x = b, given a factored matrix, A = LU. 1450 1451 Collective on Mat and Vec 1452 1453 Input Parameters: 1454 + mat - the factored matrix 1455 - b - the right-hand-side vector 1456 1457 Output Parameter: 1458 . x - the result vector 1459 1460 Notes: 1461 MatSolve() should be used for most applications, as it performs 1462 a forward solve followed by a backward solve. 1463 1464 The vectors b and x cannot be the same. I.e., one cannot 1465 call MatForwardSolve(A,x,x). 1466 1467 .keywords: matrix, forward, LU, Cholesky, triangular solve 1468 1469 .seealso: MatSolve(), MatBackwardSolve() 1470 @ */ 1471 int MatForwardSolve(Mat mat,Vec b,Vec x) 1472 { 1473 int ierr; 1474 1475 PetscFunctionBegin; 1476 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1477 PetscValidHeaderSpecific(b,VEC_COOKIE); PetscValidHeaderSpecific(x,VEC_COOKIE); 1478 if (x == b) SETERRQ(PETSC_ERR_ARG_IDN,0,"x and b must be different vectors"); 1479 if (!mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Unfactored matrix"); 1480 if (!mat->ops->forwardsolve) SETERRQ(PETSC_ERR_SUP,0,""); 1481 if (mat->N != x->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim %d %d",mat->N,x->N); 1482 if (mat->M != b->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: global dim %d %d",mat->M,b->N); 1483 if (mat->m != b->n) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: local dim %d %d",mat->m,b->n); 1484 1485 PLogEventBegin(MAT_ForwardSolve,mat,b,x,0); 1486 ierr = (*mat->ops->forwardsolve)(mat,b,x); CHKERRQ(ierr); 1487 PLogEventEnd(MAT_ForwardSolve,mat,b,x,0); 1488 PetscFunctionReturn(0); 1489 } 1490 1491 #undef __FUNC__ 1492 #define __FUNC__ "MatBackwardSolve" 1493 /* @ 1494 MatBackwardSolve - Solves U x = b, given a factored matrix, A = LU. 1495 1496 Collective on Mat and Vec 1497 1498 Input Parameters: 1499 + mat - the factored matrix 1500 - b - the right-hand-side vector 1501 1502 Output Parameter: 1503 . x - the result vector 1504 1505 Notes: 1506 MatSolve() should be used for most applications, as it performs 1507 a forward solve followed by a backward solve. 1508 1509 The vectors b and x cannot be the same. I.e., one cannot 1510 call MatBackwardSolve(A,x,x). 1511 1512 .keywords: matrix, backward, LU, Cholesky, triangular solve 1513 1514 .seealso: MatSolve(), MatForwardSolve() 1515 @ */ 1516 int MatBackwardSolve(Mat mat,Vec b,Vec x) 1517 { 1518 int ierr; 1519 1520 PetscFunctionBegin; 1521 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1522 PetscValidHeaderSpecific(b,VEC_COOKIE); PetscValidHeaderSpecific(x,VEC_COOKIE); 1523 if (x == b) SETERRQ(PETSC_ERR_ARG_IDN,0,"x and b must be different vectors"); 1524 if (!mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Unfactored matrix"); 1525 if (!mat->ops->backwardsolve) SETERRQ(PETSC_ERR_SUP,0,""); 1526 if (mat->N != x->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim %d %d",mat->N,x->N); 1527 if (mat->M != b->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: global dim %d %d",mat->M,b->N); 1528 if (mat->m != b->n) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: local dim %d %d",mat->m,b->n); 1529 1530 PLogEventBegin(MAT_BackwardSolve,mat,b,x,0); 1531 ierr = (*mat->ops->backwardsolve)(mat,b,x); CHKERRQ(ierr); 1532 PLogEventEnd(MAT_BackwardSolve,mat,b,x,0); 1533 PetscFunctionReturn(0); 1534 } 1535 1536 #undef __FUNC__ 1537 #define __FUNC__ "MatSolveAdd" 1538 /*@ 1539 MatSolveAdd - Computes x = y + inv(A)*b, given a factored matrix. 1540 1541 Collective on Mat and Vec 1542 1543 Input Parameters: 1544 + mat - the factored matrix 1545 . b - the right-hand-side vector 1546 - y - the vector to be added to 1547 1548 Output Parameter: 1549 . x - the result vector 1550 1551 Notes: 1552 The vectors b and x cannot be the same. I.e., one cannot 1553 call MatSolveAdd(A,x,y,x). 1554 1555 .keywords: matrix, linear system, solve, LU, Cholesky, add 1556 1557 .seealso: MatSolve(), MatSolveTrans(), MatSolveTransAdd() 1558 @*/ 1559 int MatSolveAdd(Mat mat,Vec b,Vec y,Vec x) 1560 { 1561 Scalar one = 1.0; 1562 Vec tmp; 1563 int ierr; 1564 1565 PetscFunctionBegin; 1566 PetscValidHeaderSpecific(mat,MAT_COOKIE);PetscValidHeaderSpecific(y,VEC_COOKIE); 1567 PetscValidHeaderSpecific(b,VEC_COOKIE); PetscValidHeaderSpecific(x,VEC_COOKIE); 1568 if (x == b) SETERRQ(PETSC_ERR_ARG_IDN,0,"x and b must be different vectors"); 1569 if (!mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Unfactored matrix"); 1570 if (mat->N != x->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim %d %d",mat->N,x->N); 1571 if (mat->M != b->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: global dim %d %d",mat->M,b->N); 1572 if (mat->M != y->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec y: global dim %d %d",mat->M,y->N); 1573 if (mat->m != b->n) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: local dim %d %d",mat->m,b->n); 1574 if (x->n != y->n) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Vec x,Vec y: local dim %d %d",x->n,y->n); 1575 1576 PLogEventBegin(MAT_SolveAdd,mat,b,x,y); 1577 if (mat->ops->solveadd) { 1578 ierr = (*mat->ops->solveadd)(mat,b,y,x); CHKERRQ(ierr); 1579 } else { 1580 /* do the solve then the add manually */ 1581 if (x != y) { 1582 ierr = MatSolve(mat,b,x); CHKERRQ(ierr); 1583 ierr = VecAXPY(&one,y,x); CHKERRQ(ierr); 1584 } else { 1585 ierr = VecDuplicate(x,&tmp); CHKERRQ(ierr); 1586 PLogObjectParent(mat,tmp); 1587 ierr = VecCopy(x,tmp); CHKERRQ(ierr); 1588 ierr = MatSolve(mat,b,x); CHKERRQ(ierr); 1589 ierr = VecAXPY(&one,tmp,x); CHKERRQ(ierr); 1590 ierr = VecDestroy(tmp); CHKERRQ(ierr); 1591 } 1592 } 1593 PLogEventEnd(MAT_SolveAdd,mat,b,x,y); 1594 PetscFunctionReturn(0); 1595 } 1596 1597 #undef __FUNC__ 1598 #define __FUNC__ "MatSolveTrans" 1599 /*@ 1600 MatSolveTrans - Solves A' x = b, given a factored matrix. 1601 1602 Collective on Mat and Vec 1603 1604 Input Parameters: 1605 + mat - the factored matrix 1606 - b - the right-hand-side vector 1607 1608 Output Parameter: 1609 . x - the result vector 1610 1611 Notes: 1612 The vectors b and x cannot be the same. I.e., one cannot 1613 call MatSolveTrans(A,x,x). 1614 1615 .keywords: matrix, linear system, solve, LU, Cholesky, transpose 1616 1617 .seealso: MatSolve(), MatSolveAdd(), MatSolveTransAdd() 1618 @*/ 1619 int MatSolveTrans(Mat mat,Vec b,Vec x) 1620 { 1621 int ierr; 1622 1623 PetscFunctionBegin; 1624 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1625 PetscValidHeaderSpecific(b,VEC_COOKIE); PetscValidHeaderSpecific(x,VEC_COOKIE); 1626 if (!mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Unfactored matrix"); 1627 if (x == b) SETERRQ(PETSC_ERR_ARG_IDN,0,"x and b must be different vectors"); 1628 if (!mat->ops->solvetrans) SETERRQ(PETSC_ERR_SUP,0,""); 1629 if (mat->M != x->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim %d %d",mat->M,x->N); 1630 if (mat->N != b->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: global dim %d %d",mat->N,b->N); 1631 1632 PLogEventBegin(MAT_SolveTrans,mat,b,x,0); 1633 ierr = (*mat->ops->solvetrans)(mat,b,x); CHKERRQ(ierr); 1634 PLogEventEnd(MAT_SolveTrans,mat,b,x,0); 1635 PetscFunctionReturn(0); 1636 } 1637 1638 #undef __FUNC__ 1639 #define __FUNC__ "MatSolveTransAdd" 1640 /*@ 1641 MatSolveTransAdd - Computes x = y + inv(trans(A)) b, given a 1642 factored matrix. 1643 1644 Collective on Mat and Vec 1645 1646 Input Parameters: 1647 + mat - the factored matrix 1648 . b - the right-hand-side vector 1649 - y - the vector to be added to 1650 1651 Output Parameter: 1652 . x - the result vector 1653 1654 Notes: 1655 The vectors b and x cannot be the same. I.e., one cannot 1656 call MatSolveTransAdd(A,x,y,x). 1657 1658 .keywords: matrix, linear system, solve, LU, Cholesky, transpose, add 1659 1660 .seealso: MatSolve(), MatSolveAdd(), MatSolveTrans() 1661 @*/ 1662 int MatSolveTransAdd(Mat mat,Vec b,Vec y,Vec x) 1663 { 1664 Scalar one = 1.0; 1665 int ierr; 1666 Vec tmp; 1667 1668 PetscFunctionBegin; 1669 PetscValidHeaderSpecific(mat,MAT_COOKIE);PetscValidHeaderSpecific(y,VEC_COOKIE); 1670 PetscValidHeaderSpecific(b,VEC_COOKIE); PetscValidHeaderSpecific(x,VEC_COOKIE); 1671 if (x == b) SETERRQ(PETSC_ERR_ARG_IDN,0,"x and b must be different vectors"); 1672 if (!mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Unfactored matrix"); 1673 if (mat->M != x->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim %d %d",mat->M,x->N); 1674 if (mat->N != b->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: global dim %d %d",mat->N,b->N); 1675 if (mat->N != y->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec y: global dim %d %d",mat->N,y->N); 1676 if (x->n != y->n) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Vec x,Vec y: local dim %d %d",x->n,y->n); 1677 1678 PLogEventBegin(MAT_SolveTransAdd,mat,b,x,y); 1679 if (mat->ops->solvetransadd) { 1680 ierr = (*mat->ops->solvetransadd)(mat,b,y,x); CHKERRQ(ierr); 1681 } else { 1682 /* do the solve then the add manually */ 1683 if (x != y) { 1684 ierr = MatSolveTrans(mat,b,x); CHKERRQ(ierr); 1685 ierr = VecAXPY(&one,y,x); CHKERRQ(ierr); 1686 } else { 1687 ierr = VecDuplicate(x,&tmp); CHKERRQ(ierr); 1688 PLogObjectParent(mat,tmp); 1689 ierr = VecCopy(x,tmp); CHKERRQ(ierr); 1690 ierr = MatSolveTrans(mat,b,x); CHKERRQ(ierr); 1691 ierr = VecAXPY(&one,tmp,x); CHKERRQ(ierr); 1692 ierr = VecDestroy(tmp); CHKERRQ(ierr); 1693 } 1694 } 1695 PLogEventEnd(MAT_SolveTransAdd,mat,b,x,y); 1696 PetscFunctionReturn(0); 1697 } 1698 /* ----------------------------------------------------------------*/ 1699 1700 #undef __FUNC__ 1701 #define __FUNC__ "MatRelax" 1702 /*@ 1703 MatRelax - Computes one relaxation sweep. 1704 1705 Collective on Mat and Vec 1706 1707 Input Parameters: 1708 + mat - the matrix 1709 . b - the right hand side 1710 . omega - the relaxation factor 1711 . flag - flag indicating the type of SOR (see below) 1712 . shift - diagonal shift 1713 - its - the number of iterations 1714 1715 Output Parameters: 1716 . x - the solution (can contain an initial guess) 1717 1718 SOR Flags: 1719 . SOR_FORWARD_SWEEP - forward SOR 1720 . SOR_BACKWARD_SWEEP - backward SOR 1721 . SOR_SYMMETRIC_SWEEP - SSOR (symmetric SOR) 1722 . SOR_LOCAL_FORWARD_SWEEP - local forward SOR 1723 . SOR_LOCAL_BACKWARD_SWEEP - local forward SOR 1724 . SOR_LOCAL_SYMMETRIC_SWEEP - local SSOR 1725 . SOR_APPLY_UPPER, SOR_APPLY_LOWER - applies 1726 upper/lower triangular part of matrix to 1727 vector (with omega) 1728 . SOR_ZERO_INITIAL_GUESS - zero initial guess 1729 1730 Notes: 1731 SOR_LOCAL_FORWARD_SWEEP, SOR_LOCAL_BACKWARD_SWEEP, and 1732 SOR_LOCAL_SYMMETRIC_SWEEP perform seperate independent smoothings 1733 on each processor. 1734 1735 Application programmers will not generally use MatRelax() directly, 1736 but instead will employ the SLES/PC interface. 1737 1738 Notes for Advanced Users: 1739 The flags are implemented as bitwise inclusive or operations. 1740 For example, use (SOR_ZERO_INITIAL_GUESS | SOR_SYMMETRIC_SWEEP) 1741 to specify a zero initial guess for SSOR. 1742 1743 .keywords: matrix, relax, relaxation, sweep 1744 @*/ 1745 int MatRelax(Mat mat,Vec b,double omega,MatSORType flag,double shift, 1746 int its,Vec x) 1747 { 1748 int ierr; 1749 1750 PetscFunctionBegin; 1751 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1752 PetscValidHeaderSpecific(b,VEC_COOKIE); PetscValidHeaderSpecific(x,VEC_COOKIE); 1753 if (!mat->ops->relax) SETERRQ(PETSC_ERR_SUP,0,""); 1754 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 1755 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 1756 if (mat->N != x->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim %d %d",mat->N,x->N); 1757 if (mat->M != b->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: global dim %d %d",mat->M,b->N); 1758 if (mat->m != b->n) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: local dim %d %d",mat->m,b->n); 1759 1760 PLogEventBegin(MAT_Relax,mat,b,x,0); 1761 ierr =(*mat->ops->relax)(mat,b,omega,flag,shift,its,x); CHKERRQ(ierr); 1762 PLogEventEnd(MAT_Relax,mat,b,x,0); 1763 PetscFunctionReturn(0); 1764 } 1765 1766 #undef __FUNC__ 1767 #define __FUNC__ "MatCopy_Basic" 1768 /* 1769 Default matrix copy routine. 1770 */ 1771 int MatCopy_Basic(Mat A,Mat B,MatStructure str) 1772 { 1773 int ierr,i,rstart,rend,nz,*cwork; 1774 Scalar *vwork; 1775 1776 PetscFunctionBegin; 1777 ierr = MatZeroEntries(B); CHKERRQ(ierr); 1778 ierr = MatGetOwnershipRange(A,&rstart,&rend); CHKERRQ(ierr); 1779 for (i=rstart; i<rend; i++) { 1780 ierr = MatGetRow(A,i,&nz,&cwork,&vwork); CHKERRQ(ierr); 1781 ierr = MatSetValues(B,1,&i,nz,cwork,vwork,INSERT_VALUES); CHKERRQ(ierr); 1782 ierr = MatRestoreRow(A,i,&nz,&cwork,&vwork); CHKERRQ(ierr); 1783 } 1784 ierr = MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY); CHKERRQ(ierr); 1785 ierr = MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY); CHKERRQ(ierr); 1786 PetscFunctionReturn(0); 1787 } 1788 1789 #undef __FUNC__ 1790 #define __FUNC__ "MatCopy" 1791 /*@C 1792 MatCopy - Copys a matrix to another matrix. 1793 1794 Collective on Mat 1795 1796 Input Parameters: 1797 + A - the matrix 1798 - str - SAME_NONZERO_PATTERN or DIFFERENT_NONZERO_PATTERN 1799 1800 Output Parameter: 1801 . B - where the copy is put 1802 1803 Notes: 1804 MatCopy() copies the matrix entries of a matrix to another existing 1805 matrix (after first zeroing the second matrix). A related routine is 1806 MatConvert(), which first creates a new matrix and then copies the data. 1807 1808 .keywords: matrix, copy, convert 1809 1810 .seealso: MatConvert() 1811 @*/ 1812 int MatCopy(Mat A,Mat B,MatStructure str) 1813 { 1814 int ierr; 1815 1816 PetscFunctionBegin; 1817 PetscValidHeaderSpecific(A,MAT_COOKIE); PetscValidHeaderSpecific(B,MAT_COOKIE); 1818 if (!A->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 1819 if (A->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 1820 if (A->M != B->M || A->N != B->N) SETERRQ4(PETSC_ERR_ARG_SIZ,0,"Mat A,Mat B: global dim %d %d",A->M,B->M, 1821 A->N,B->N); 1822 1823 PLogEventBegin(MAT_Copy,A,B,0,0); 1824 if (A->ops->copy) { 1825 ierr = (*A->ops->copy)(A,B,str); CHKERRQ(ierr); 1826 } else { /* generic conversion */ 1827 ierr = MatCopy_Basic(A,B,str); CHKERRQ(ierr); 1828 } 1829 PLogEventEnd(MAT_Copy,A,B,0,0); 1830 PetscFunctionReturn(0); 1831 } 1832 1833 static int MatConvertersSet = 0; 1834 static int (*MatConverters[MAX_MATRIX_TYPES][MAX_MATRIX_TYPES])(Mat,MatType,Mat*) = 1835 {{0,0,0,0,0,0,0,0,0,0,0,0},{0,0,0,0,0,0,0,0,0,0,0,0}, 1836 {0,0,0,0,0,0,0,0,0,0,0,0},{0,0,0,0,0,0,0,0,0,0,0,0}, 1837 {0,0,0,0,0,0,0,0,0,0,0,0},{0,0,0,0,0,0,0,0,0,0,0,0}, 1838 {0,0,0,0,0,0,0,0,0,0,0,0},{0,0,0,0,0,0,0,0,0,0,0,0}, 1839 {0,0,0,0,0,0,0,0,0,0,0,0},{0,0,0,0,0,0,0,0,0,0,0,0}, 1840 {0,0,0,0,0,0,0,0,0,0,0,0},{0,0,0,0,0,0,0,0,0,0,0,0}}; 1841 1842 #undef __FUNC__ 1843 #define __FUNC__ "MatConvertRegister" 1844 /*@C 1845 MatConvertRegister - Allows one to register a routine that converts between 1846 two matrix types. 1847 1848 Not Collective 1849 1850 Input Parameters: 1851 . intype - the type of matrix (defined in include/mat.h), for example, MATSEQAIJ. 1852 . outtype - new matrix type, or MATSAME 1853 1854 .seealso: MatConvertRegisterAll() 1855 @*/ 1856 int MatConvertRegister(MatType intype,MatType outtype,int (*converter)(Mat,MatType,Mat*)) 1857 { 1858 PetscFunctionBegin; 1859 MatConverters[intype][outtype] = converter; 1860 MatConvertersSet = 1; 1861 PetscFunctionReturn(0); 1862 } 1863 1864 #undef __FUNC__ 1865 #define __FUNC__ "MatConvert" 1866 /*@C 1867 MatConvert - Converts a matrix to another matrix, either of the same 1868 or different type. 1869 1870 Collective on Mat 1871 1872 Input Parameters: 1873 + mat - the matrix 1874 - newtype - new matrix type. Use MATSAME to create a new matrix of the 1875 same type as the original matrix. 1876 1877 Output Parameter: 1878 . M - pointer to place new matrix 1879 1880 Notes: 1881 MatConvert() first creates a new matrix and then copies the data from 1882 the first matrix. A related routine is MatCopy(), which copies the matrix 1883 entries of one matrix to another already existing matrix context. 1884 1885 .keywords: matrix, copy, convert 1886 1887 .seealso: MatCopy(), MatDuplicate() 1888 @*/ 1889 int MatConvert(Mat mat,MatType newtype,Mat *M) 1890 { 1891 int ierr; 1892 1893 PetscFunctionBegin; 1894 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1895 PetscValidPointer(M); 1896 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 1897 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 1898 1899 if (newtype > MAX_MATRIX_TYPES || newtype < -1) { 1900 SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,1,"Not a valid matrix type"); 1901 } 1902 *M = 0; 1903 1904 if (!MatConvertersSet) { 1905 ierr = MatLoadRegisterAll(); CHKERRQ(ierr); 1906 } 1907 1908 PLogEventBegin(MAT_Convert,mat,0,0,0); 1909 if ((newtype == mat->type || newtype == MATSAME) && mat->ops->duplicate) { 1910 ierr = (*mat->ops->duplicate)(mat,MAT_COPY_VALUES,M); CHKERRQ(ierr); 1911 } else { 1912 if (!MatConvertersSet) { 1913 ierr = MatConvertRegisterAll(); CHKERRQ(ierr); 1914 } 1915 if (!MatConverters[mat->type][newtype]) { 1916 SETERRQ(PETSC_ERR_ARG_WRONG,1,"Invalid matrix type, or matrix converter not registered"); 1917 } 1918 ierr = (*MatConverters[mat->type][newtype])(mat,newtype,M); CHKERRQ(ierr); 1919 } 1920 PLogEventEnd(MAT_Convert,mat,0,0,0); 1921 PetscFunctionReturn(0); 1922 } 1923 1924 #undef __FUNC__ 1925 #define __FUNC__ "MatDuplicate" 1926 /*@C 1927 MatDuplicate - Duplicates a matrix including the non-zero structure. 1928 1929 Collective on Mat 1930 1931 Input Parameters: 1932 + mat - the matrix 1933 - op - either MAT_DO_NO_COPY_VALUES or MAT_COPY_VALUES, cause it to copy nonzero 1934 values as well or not 1935 1936 Output Parameter: 1937 . M - pointer to place new matrix 1938 1939 .keywords: matrix, copy, convert, duplicate 1940 1941 .seealso: MatCopy(), MatConvert() 1942 @*/ 1943 int MatDuplicate(Mat mat,MatDuplicateOption op,Mat *M) 1944 { 1945 int ierr; 1946 1947 PetscFunctionBegin; 1948 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1949 PetscValidPointer(M); 1950 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 1951 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 1952 1953 *M = 0; 1954 PLogEventBegin(MAT_Convert,mat,0,0,0); 1955 if (!mat->ops->duplicate) { 1956 SETERRQ(PETSC_ERR_SUP,1,"Not written for this matrix type"); 1957 } 1958 ierr = (*mat->ops->duplicate)(mat,op,M); CHKERRQ(ierr); 1959 PLogEventEnd(MAT_Convert,mat,0,0,0); 1960 PetscFunctionReturn(0); 1961 } 1962 1963 #undef __FUNC__ 1964 #define __FUNC__ "MatGetDiagonal" 1965 /*@ 1966 MatGetDiagonal - Gets the diagonal of a matrix. 1967 1968 Collective on Mat and Vec 1969 1970 Input Parameters: 1971 + mat - the matrix 1972 - v - the vector for storing the diagonal 1973 1974 Output Parameter: 1975 . v - the diagonal of the matrix 1976 1977 Notes: 1978 For the SeqAIJ matrix format, this routine may also be called 1979 on a LU factored matrix; in that case it routines the reciprocal of 1980 the diagonal entries in U. It returns the entries permuted by the 1981 row and column permutation used during the symbolic factorization. 1982 1983 .keywords: matrix, get, diagonal 1984 @*/ 1985 int MatGetDiagonal(Mat mat,Vec v) 1986 { 1987 int ierr; 1988 1989 PetscFunctionBegin; 1990 PetscValidHeaderSpecific(mat,MAT_COOKIE);PetscValidHeaderSpecific(v,VEC_COOKIE); 1991 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 1992 if (!mat->ops->getdiagonal) SETERRQ(PETSC_ERR_SUP,0,""); 1993 ierr = (*mat->ops->getdiagonal)(mat,v);CHKERRQ(ierr); 1994 PetscFunctionReturn(0); 1995 } 1996 1997 #undef __FUNC__ 1998 #define __FUNC__ "MatTranspose" 1999 /*@C 2000 MatTranspose - Computes an in-place or out-of-place transpose of a matrix. 2001 2002 Collective on Mat 2003 2004 Input Parameter: 2005 . mat - the matrix to transpose 2006 2007 Output Parameters: 2008 . B - the transpose (or pass in PETSC_NULL for an in-place transpose) 2009 2010 .keywords: matrix, transpose 2011 2012 .seealso: MatMultTrans(), MatMultTransAdd() 2013 @*/ 2014 int MatTranspose(Mat mat,Mat *B) 2015 { 2016 int ierr; 2017 2018 PetscFunctionBegin; 2019 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2020 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2021 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2022 if (!mat->ops->transpose) SETERRQ(PETSC_ERR_SUP,0,""); 2023 ierr = (*mat->ops->transpose)(mat,B);CHKERRQ(ierr); 2024 PetscFunctionReturn(0); 2025 } 2026 2027 #undef __FUNC__ 2028 #define __FUNC__ "MatPermute" 2029 /*@C 2030 MatPermute - Creates a new matrix with rows and columns permuted from the 2031 original. 2032 2033 Collective on Mat 2034 2035 Input Parameters: 2036 + mat - the matrix to permute 2037 . row - row permutation 2038 - col - column permutation 2039 2040 Output Parameters: 2041 . B - the permuted matrix 2042 2043 .keywords: matrix, transpose 2044 2045 .seealso: MatGetReordering() 2046 @*/ 2047 int MatPermute(Mat mat,IS row,IS col,Mat *B) 2048 { 2049 int ierr; 2050 2051 PetscFunctionBegin; 2052 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2053 PetscValidHeaderSpecific(row,IS_COOKIE); 2054 PetscValidHeaderSpecific(col,IS_COOKIE); 2055 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2056 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2057 if (!mat->ops->permute) SETERRQ(PETSC_ERR_SUP,0,""); 2058 ierr = (*mat->ops->permute)(mat,row,col,B);CHKERRQ(ierr); 2059 PetscFunctionReturn(0); 2060 } 2061 2062 #undef __FUNC__ 2063 #define __FUNC__ "MatEqual" 2064 /*@ 2065 MatEqual - Compares two matrices. 2066 2067 Collective on Mat 2068 2069 Input Parameters: 2070 + A - the first matrix 2071 - B - the second matrix 2072 2073 Output Parameter: 2074 . flg - PETSC_TRUE if the matrices are equal; PETSC_FALSE otherwise. 2075 2076 .keywords: matrix, equal, equivalent 2077 @*/ 2078 int MatEqual(Mat A,Mat B,PetscTruth *flg) 2079 { 2080 int ierr; 2081 2082 PetscFunctionBegin; 2083 PetscValidHeaderSpecific(A,MAT_COOKIE); PetscValidHeaderSpecific(B,MAT_COOKIE); 2084 PetscValidIntPointer(flg); 2085 if (!A->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2086 if (!B->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2087 if (A->M != B->M || A->N != B->N) SETERRQ4(PETSC_ERR_ARG_SIZ,0,"Mat A,Mat B: global dim %d %d %d %d", 2088 A->M,B->M,A->N,B->N); 2089 if (!A->ops->equal) SETERRQ(PETSC_ERR_SUP,0,""); 2090 ierr = (*A->ops->equal)(A,B,flg);CHKERRQ(ierr); 2091 PetscFunctionReturn(0); 2092 } 2093 2094 #undef __FUNC__ 2095 #define __FUNC__ "MatDiagonalScale" 2096 /*@ 2097 MatDiagonalScale - Scales a matrix on the left and right by diagonal 2098 matrices that are stored as vectors. Either of the two scaling 2099 matrices can be PETSC_NULL. 2100 2101 Collective on Mat 2102 2103 Input Parameters: 2104 + mat - the matrix to be scaled 2105 . l - the left scaling vector (or PETSC_NULL) 2106 - r - the right scaling vector (or PETSC_NULL) 2107 2108 Notes: 2109 MatDiagonalScale() computes A = LAR, where 2110 L = a diagonal matrix, R = a diagonal matrix 2111 2112 .keywords: matrix, diagonal, scale 2113 2114 .seealso: MatDiagonalScale() 2115 @*/ 2116 int MatDiagonalScale(Mat mat,Vec l,Vec r) 2117 { 2118 int ierr; 2119 2120 PetscFunctionBegin; 2121 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2122 if (!mat->ops->diagonalscale) SETERRQ(PETSC_ERR_SUP,0,""); 2123 if (l) PetscValidHeaderSpecific(l,VEC_COOKIE); 2124 if (r) PetscValidHeaderSpecific(r,VEC_COOKIE); 2125 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2126 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2127 2128 PLogEventBegin(MAT_Scale,mat,0,0,0); 2129 ierr = (*mat->ops->diagonalscale)(mat,l,r); CHKERRQ(ierr); 2130 PLogEventEnd(MAT_Scale,mat,0,0,0); 2131 PetscFunctionReturn(0); 2132 } 2133 2134 #undef __FUNC__ 2135 #define __FUNC__ "MatScale" 2136 /*@ 2137 MatScale - Scales all elements of a matrix by a given number. 2138 2139 Collective on Mat 2140 2141 Input Parameters: 2142 + mat - the matrix to be scaled 2143 - a - the scaling value 2144 2145 Output Parameter: 2146 . mat - the scaled matrix 2147 2148 .keywords: matrix, scale 2149 2150 .seealso: MatDiagonalScale() 2151 @*/ 2152 int MatScale(Scalar *a,Mat mat) 2153 { 2154 int ierr; 2155 2156 PetscFunctionBegin; 2157 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2158 PetscValidScalarPointer(a); 2159 if (!mat->ops->scale) SETERRQ(PETSC_ERR_SUP,0,""); 2160 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2161 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2162 2163 PLogEventBegin(MAT_Scale,mat,0,0,0); 2164 ierr = (*mat->ops->scale)(a,mat); CHKERRQ(ierr); 2165 PLogEventEnd(MAT_Scale,mat,0,0,0); 2166 PetscFunctionReturn(0); 2167 } 2168 2169 #undef __FUNC__ 2170 #define __FUNC__ "MatNorm" 2171 /*@ 2172 MatNorm - Calculates various norms of a matrix. 2173 2174 Collective on Mat 2175 2176 Input Parameters: 2177 + mat - the matrix 2178 - type - the type of norm, NORM_1, NORM_2, NORM_FROBENIUS, NORM_INFINITY 2179 2180 Output Parameters: 2181 . norm - the resulting norm 2182 2183 .keywords: matrix, norm, Frobenius 2184 @*/ 2185 int MatNorm(Mat mat,NormType type,double *norm) 2186 { 2187 int ierr; 2188 2189 PetscFunctionBegin; 2190 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2191 PetscValidScalarPointer(norm); 2192 2193 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2194 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2195 if (!mat->ops->norm) SETERRQ(PETSC_ERR_SUP,0,"Not for this matrix type"); 2196 ierr = (*mat->ops->norm)(mat,type,norm);CHKERRQ(ierr); 2197 PetscFunctionReturn(0); 2198 } 2199 2200 /* 2201 This variable is used to prevent counting of MatAssemblyBegin() that 2202 are called from within a MatAssemblyEnd(). 2203 */ 2204 static int MatAssemblyEnd_InUse = 0; 2205 #undef __FUNC__ 2206 #define __FUNC__ "MatAssemblyBegin" 2207 /*@ 2208 MatAssemblyBegin - Begins assembling the matrix. This routine should 2209 be called after completing all calls to MatSetValues(). 2210 2211 Collective on Mat 2212 2213 Input Parameters: 2214 + mat - the matrix 2215 - type - type of assembly, either MAT_FLUSH_ASSEMBLY or MAT_FINAL_ASSEMBLY 2216 2217 Notes: 2218 MatSetValues() generally caches the values. The matrix is ready to 2219 use only after MatAssemblyBegin() and MatAssemblyEnd() have been called. 2220 Use MAT_FLUSH_ASSEMBLY when switching between ADD_VALUES and INSERT_VALUES 2221 in MatSetValues(); use MAT_FINAL_ASSEMBLY for the final assembly before 2222 using the matrix. 2223 2224 Level: beginner 2225 2226 .keywords: matrix, assembly, assemble, begin 2227 2228 .seealso: MatAssemblyEnd(), MatSetValues() 2229 @*/ 2230 int MatAssemblyBegin(Mat mat,MatAssemblyType type) 2231 { 2232 int ierr; 2233 2234 PetscFunctionBegin; 2235 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2236 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix.\n did you forget to call MatSetUnfactored()?"); 2237 if (mat->assembled) { 2238 mat->was_assembled = PETSC_TRUE; 2239 mat->assembled = PETSC_FALSE; 2240 } 2241 if (!MatAssemblyEnd_InUse) { 2242 PLogEventBegin(MAT_AssemblyBegin,mat,0,0,0); 2243 if (mat->ops->assemblybegin){ierr = (*mat->ops->assemblybegin)(mat,type);CHKERRQ(ierr);} 2244 PLogEventEnd(MAT_AssemblyBegin,mat,0,0,0); 2245 } else { 2246 if (mat->ops->assemblybegin){ierr = (*mat->ops->assemblybegin)(mat,type);CHKERRQ(ierr);} 2247 } 2248 PetscFunctionReturn(0); 2249 } 2250 2251 2252 #undef __FUNC__ 2253 #define __FUNC__ "MatView_Private" 2254 /* 2255 Processes command line options to determine if/how a matrix 2256 is to be viewed. Called by MatAssemblyEnd() and MatLoad(). 2257 */ 2258 int MatView_Private(Mat mat) 2259 { 2260 int ierr,flg; 2261 2262 PetscFunctionBegin; 2263 ierr = OptionsHasName(PETSC_NULL,"-mat_view_info",&flg); CHKERRQ(ierr); 2264 if (flg) { 2265 ierr = ViewerPushFormat(VIEWER_STDOUT_(mat->comm),VIEWER_FORMAT_ASCII_INFO,0);CHKERRQ(ierr); 2266 ierr = MatView(mat,VIEWER_STDOUT_(mat->comm)); CHKERRQ(ierr); 2267 ierr = ViewerPopFormat(VIEWER_STDOUT_(mat->comm));CHKERRQ(ierr); 2268 } 2269 ierr = OptionsHasName(PETSC_NULL,"-mat_view_info_detailed",&flg);CHKERRQ(ierr); 2270 if (flg) { 2271 ierr = ViewerPushFormat(VIEWER_STDOUT_(mat->comm),VIEWER_FORMAT_ASCII_INFO_LONG,0);CHKERRQ(ierr); 2272 ierr = MatView(mat,VIEWER_STDOUT_(mat->comm)); CHKERRQ(ierr); 2273 ierr = ViewerPopFormat(VIEWER_STDOUT_(mat->comm));CHKERRQ(ierr); 2274 } 2275 ierr = OptionsHasName(PETSC_NULL,"-mat_view",&flg); CHKERRQ(ierr); 2276 if (flg) { 2277 ierr = MatView(mat,VIEWER_STDOUT_(mat->comm)); CHKERRQ(ierr); 2278 } 2279 ierr = OptionsHasName(PETSC_NULL,"-mat_view_matlab",&flg); CHKERRQ(ierr); 2280 if (flg) { 2281 ierr = ViewerPushFormat(VIEWER_STDOUT_(mat->comm),VIEWER_FORMAT_ASCII_MATLAB,"M");CHKERRQ(ierr); 2282 ierr = MatView(mat,VIEWER_STDOUT_(mat->comm)); CHKERRQ(ierr); 2283 ierr = ViewerPopFormat(VIEWER_STDOUT_(mat->comm));CHKERRQ(ierr); 2284 } 2285 ierr = OptionsHasName(PETSC_NULL,"-mat_view_draw",&flg); CHKERRQ(ierr); 2286 if (flg) { 2287 ierr = OptionsHasName(0,"-mat_view_contour",&flg); CHKERRQ(ierr); 2288 if (flg) { 2289 ViewerPushFormat(VIEWER_DRAW_(mat->comm),VIEWER_FORMAT_DRAW_CONTOUR,0);CHKERRQ(ierr); 2290 } 2291 ierr = MatView(mat,VIEWER_DRAW_(mat->comm)); CHKERRQ(ierr); 2292 ierr = ViewerFlush(VIEWER_DRAW_(mat->comm)); CHKERRQ(ierr); 2293 if (flg) { 2294 ViewerPopFormat(VIEWER_DRAW_(mat->comm));CHKERRQ(ierr); 2295 } 2296 } 2297 ierr = OptionsHasName(PETSC_NULL,"-mat_view_socket",&flg); CHKERRQ(ierr); 2298 if (flg) { 2299 ierr = MatView(mat,VIEWER_SOCKET_(mat->comm)); CHKERRQ(ierr); 2300 ierr = ViewerFlush(VIEWER_SOCKET_(mat->comm)); CHKERRQ(ierr); 2301 } 2302 PetscFunctionReturn(0); 2303 } 2304 2305 #undef __FUNC__ 2306 #define __FUNC__ "MatAssemblyEnd" 2307 /*@ 2308 MatAssemblyEnd - Completes assembling the matrix. This routine should 2309 be called after MatAssemblyBegin(). 2310 2311 Collective on Mat 2312 2313 Input Parameters: 2314 + mat - the matrix 2315 - type - type of assembly, either MAT_FLUSH_ASSEMBLY or MAT_FINAL_ASSEMBLY 2316 2317 Options Database Keys: 2318 + -mat_view_info - Prints info on matrix at conclusion of MatEndAssembly() 2319 . -mat_view_info_detailed - Prints more detailed info 2320 . -mat_view - Prints matrix in ASCII format 2321 . -mat_view_matlab - Prints matrix in Matlab format 2322 . -mat_view_draw - Draws nonzero structure of matrix, using MatView() and DrawOpenX(). 2323 . -display <name> - Sets display name (default is host) 2324 - -draw_pause <sec> - Sets number of seconds to pause after display 2325 2326 Notes: 2327 MatSetValues() generally caches the values. The matrix is ready to 2328 use only after MatAssemblyBegin() and MatAssemblyEnd() have been called. 2329 Use MAT_FLUSH_ASSEMBLY when switching between ADD_VALUES and INSERT_VALUES 2330 in MatSetValues(); use MAT_FINAL_ASSEMBLY for the final assembly before 2331 using the matrix. 2332 2333 Level: beginner 2334 2335 .keywords: matrix, assembly, assemble, end 2336 2337 .seealso: MatAssemblyBegin(), MatSetValues(), DrawOpenX(), MatView() 2338 @*/ 2339 int MatAssemblyEnd(Mat mat,MatAssemblyType type) 2340 { 2341 int ierr; 2342 static int inassm = 0; 2343 2344 PetscFunctionBegin; 2345 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2346 2347 inassm++; 2348 MatAssemblyEnd_InUse++; 2349 PLogEventBegin(MAT_AssemblyEnd,mat,0,0,0); 2350 if (mat->ops->assemblyend) { 2351 ierr = (*mat->ops->assemblyend)(mat,type); CHKERRQ(ierr); 2352 } 2353 2354 /* Flush assembly is not a true assembly */ 2355 if (type != MAT_FLUSH_ASSEMBLY) { 2356 mat->assembled = PETSC_TRUE; mat->num_ass++; 2357 } 2358 mat->insertmode = NOT_SET_VALUES; 2359 PLogEventEnd(MAT_AssemblyEnd,mat,0,0,0); 2360 MatAssemblyEnd_InUse--; 2361 2362 if (inassm == 1 && type != MAT_FLUSH_ASSEMBLY) { 2363 ierr = MatView_Private(mat); CHKERRQ(ierr); 2364 } 2365 inassm--; 2366 PetscFunctionReturn(0); 2367 } 2368 2369 2370 #undef __FUNC__ 2371 #define __FUNC__ "MatCompress" 2372 /*@ 2373 MatCompress - Tries to store the matrix in as little space as 2374 possible. May fail if memory is already fully used, since it 2375 tries to allocate new space. 2376 2377 Collective on Mat 2378 2379 Input Parameters: 2380 . mat - the matrix 2381 2382 .keywords: matrix, compress 2383 @*/ 2384 int MatCompress(Mat mat) 2385 { 2386 int ierr; 2387 2388 PetscFunctionBegin; 2389 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2390 if (mat->ops->compress) {ierr = (*mat->ops->compress)(mat);CHKERRQ(ierr);} 2391 PetscFunctionReturn(0); 2392 } 2393 2394 #undef __FUNC__ 2395 #define __FUNC__ "MatSetOption" 2396 /*@ 2397 MatSetOption - Sets a parameter option for a matrix. Some options 2398 may be specific to certain storage formats. Some options 2399 determine how values will be inserted (or added). Sorted, 2400 row-oriented input will generally assemble the fastest. The default 2401 is row-oriented, nonsorted input. 2402 2403 Collective on Mat 2404 2405 Input Parameters: 2406 + mat - the matrix 2407 - option - the option, one of those listed below (and possibly others), 2408 e.g., MAT_ROWS_SORTED, MAT_NEW_NONZERO_LOCATION_ERROR 2409 2410 Options Describing Matrix Structure: 2411 + MAT_SYMMETRIC - symmetric in terms of both structure and value 2412 - MAT_STRUCTURALLY_SYMMETRIC - symmetric nonzero structure 2413 2414 Options For Use with MatSetValues(): 2415 Insert a logically dense subblock, which can be 2416 + MAT_ROW_ORIENTED - row-oriented 2417 . MAT_COLUMN_ORIENTED - column-oriented 2418 . MAT_ROWS_SORTED - sorted by row 2419 . MAT_ROWS_UNSORTED - not sorted by row 2420 . MAT_COLUMNS_SORTED - sorted by column 2421 - MAT_COLUMNS_UNSORTED - not sorted by column 2422 2423 Not these options reflect the data you pass in with MatSetValues(); it has 2424 nothing to do with how the data is stored internally in the matrix 2425 data structure. 2426 2427 When (re)assembling a matrix, we can restrict the input for 2428 efficiency/debugging purposes. These options include 2429 + MAT_NO_NEW_NONZERO_LOCATIONS - additional insertions will not be 2430 allowed if they generate a new nonzero 2431 . MAT_YES_NEW_NONZERO_LOCATIONS - additional insertions will be allowed 2432 . MAT_NO_NEW_DIAGONALS - additional insertions will not be allowed if 2433 they generate a nonzero in a new diagonal (for block diagonal format only) 2434 . MAT_YES_NEW_DIAGONALS - new diagonals will be allowed (for block diagonal format only) 2435 . MAT_IGNORE_OFF_PROC_ENTRIES - drops off-processor entries 2436 . MAT_NEW_NONZERO_LOCATION_ERROR - generates an error for new matrix entry 2437 - MAT_USE_HASH_TABLE - uses a hash table to speed up matrix assembly 2438 2439 Notes: 2440 Some options are relevant only for particular matrix types and 2441 are thus ignored by others. Other options are not supported by 2442 certain matrix types and will generate an error message if set. 2443 2444 If using a Fortran 77 module to compute a matrix, one may need to 2445 use the column-oriented option (or convert to the row-oriented 2446 format). 2447 2448 MAT_NO_NEW_NONZERO_LOCATIONS indicates that any add or insertion 2449 that would generate a new entry in the nonzero structure is instead 2450 ignored. Thus, if memory has not alredy been allocated for this particular 2451 data, then the insertion is ignored. For dense matrices, in which 2452 the entire array is allocated, no entries are ever ignored. 2453 2454 MAT_NEW_NONZERO_LOCATION_ERROR indicates that any add or insertion 2455 that would generate a new entry in the nonzero structure instead produces 2456 an error. (Currently supported for AIJ and BAIJ formats only.) 2457 This is a useful flag when using SAME_NONZERO_PATTERN in calling 2458 SLESSetOperators() to ensure that the nonzero pattern truely does 2459 remain unchanged. 2460 2461 MAT_NEW_NONZERO_ALLOCATION_ERROR indicates that any add or insertion 2462 that would generate a new entry that has not been preallocated will 2463 instead produce an error. (Currently supported for AIJ and BAIJ formats 2464 only.) This is a useful flag when debugging matrix memory preallocation. 2465 2466 MAT_IGNORE_OFF_PROC_ENTRIES indicates entries destined for 2467 other processors should be dropped, rather than stashed. 2468 This is useful if you know that the "owning" processor is also 2469 always generating the correct matrix entries, so that PETSc need 2470 not transfer duplicate entries generated on another processor. 2471 2472 MAT_USE_HASH_TABLE indicates that a hash table be used to improve the 2473 searches during matrix assembly. When this flag is set, the hash table 2474 is created during the first Matrix Assembly. This hash table is 2475 used the next time through, during MatSetVaules()/MatSetVaulesBlocked() 2476 to improve the searching of indices. MAT_NO_NEW_NONZERO_LOCATIONS flag 2477 should be used with MAT_USE_HASH_TABLE flag. This option is currently 2478 supported by MATMPIBAIJ format only. 2479 2480 .keywords: matrix, option, row-oriented, column-oriented, sorted, nonzero 2481 @*/ 2482 int MatSetOption(Mat mat,MatOption op) 2483 { 2484 int ierr; 2485 2486 PetscFunctionBegin; 2487 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2488 if (mat->ops->setoption) {ierr = (*mat->ops->setoption)(mat,op);CHKERRQ(ierr);} 2489 PetscFunctionReturn(0); 2490 } 2491 2492 #undef __FUNC__ 2493 #define __FUNC__ "MatZeroEntries" 2494 /*@ 2495 MatZeroEntries - Zeros all entries of a matrix. For sparse matrices 2496 this routine retains the old nonzero structure. 2497 2498 Collective on Mat 2499 2500 Input Parameters: 2501 . mat - the matrix 2502 2503 .keywords: matrix, zero, entries 2504 2505 .seealso: MatZeroRows() 2506 @*/ 2507 int MatZeroEntries(Mat mat) 2508 { 2509 int ierr; 2510 2511 PetscFunctionBegin; 2512 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2513 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2514 if (!mat->ops->zeroentries) SETERRQ(PETSC_ERR_SUP,0,""); 2515 2516 PLogEventBegin(MAT_ZeroEntries,mat,0,0,0); 2517 ierr = (*mat->ops->zeroentries)(mat); CHKERRQ(ierr); 2518 PLogEventEnd(MAT_ZeroEntries,mat,0,0,0); 2519 PetscFunctionReturn(0); 2520 } 2521 2522 #undef __FUNC__ 2523 #define __FUNC__ "MatZeroRows" 2524 /*@ 2525 MatZeroRows - Zeros all entries (except possibly the main diagonal) 2526 of a set of rows of a matrix. 2527 2528 Collective on Mat 2529 2530 Input Parameters: 2531 + mat - the matrix 2532 . is - index set of rows to remove 2533 - diag - pointer to value put in all diagonals of eliminated rows. 2534 Note that diag is not a pointer to an array, but merely a 2535 pointer to a single value. 2536 2537 Notes: 2538 For the AIJ matrix formats this removes the old nonzero structure, 2539 but does not release memory. For the dense and block diagonal 2540 formats this does not alter the nonzero structure. 2541 2542 The user can set a value in the diagonal entry (or for the AIJ and 2543 row formats can optionally remove the main diagonal entry from the 2544 nonzero structure as well, by passing a null pointer as the final 2545 argument). 2546 2547 For the parallel case, all processes that share the matrix (i.e., 2548 those in the communicator used for matrix creation) MUST call this 2549 routine, regardless of whether any rows being zeroed are owned by 2550 them. 2551 2552 .keywords: matrix, zero, rows, boundary conditions 2553 2554 .seealso: MatZeroEntries(), 2555 @*/ 2556 int MatZeroRows(Mat mat,IS is, Scalar *diag) 2557 { 2558 int ierr; 2559 2560 PetscFunctionBegin; 2561 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2562 PetscValidHeaderSpecific(is,IS_COOKIE); 2563 if (diag) PetscValidScalarPointer(diag); 2564 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2565 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2566 if (!mat->ops->zerorows) SETERRQ(PETSC_ERR_SUP,0,""); 2567 2568 ierr = (*mat->ops->zerorows)(mat,is,diag); CHKERRQ(ierr); 2569 ierr = MatView_Private(mat); CHKERRQ(ierr); 2570 PetscFunctionReturn(0); 2571 } 2572 2573 #undef __FUNC__ 2574 #define __FUNC__ "MatZeroRowsLocal" 2575 /*@ 2576 MatZeroRowsLocal - Zeros all entries (except possibly the main diagonal) 2577 of a set of rows of a matrix; using local numbering of rows. 2578 2579 Collective on Mat 2580 2581 Input Parameters: 2582 + mat - the matrix 2583 . is - index set of rows to remove 2584 - diag - pointer to value put in all diagonals of eliminated rows. 2585 Note that diag is not a pointer to an array, but merely a 2586 pointer to a single value. 2587 2588 Notes: 2589 For the AIJ matrix formats this removes the old nonzero structure, 2590 but does not release memory. For the dense and block diagonal 2591 formats this does not alter the nonzero structure. 2592 2593 The user can set a value in the diagonal entry (or for the AIJ and 2594 row formats can optionally remove the main diagonal entry from the 2595 nonzero structure as well, by passing a null pointer as the final 2596 argument). 2597 2598 .keywords: matrix, zero, rows, boundary conditions 2599 2600 .seealso: MatZeroEntries(), 2601 @*/ 2602 int MatZeroRowsLocal(Mat mat,IS is, Scalar *diag) 2603 { 2604 int ierr; 2605 IS newis; 2606 2607 PetscFunctionBegin; 2608 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2609 PetscValidHeaderSpecific(is,IS_COOKIE); 2610 if (diag) PetscValidScalarPointer(diag); 2611 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2612 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2613 if (!mat->ops->zerorows) SETERRQ(PETSC_ERR_SUP,0,""); 2614 2615 ierr = ISLocalToGlobalMappingApplyIS(mat->mapping,is,&newis); CHKERRQ(ierr); 2616 ierr = (*mat->ops->zerorows)(mat,newis,diag); CHKERRQ(ierr); 2617 ierr = ISDestroy(newis); 2618 PetscFunctionReturn(0); 2619 } 2620 2621 #undef __FUNC__ 2622 #define __FUNC__ "MatGetSize" 2623 /*@ 2624 MatGetSize - Returns the numbers of rows and columns in a matrix. 2625 2626 Not Collective 2627 2628 Input Parameter: 2629 . mat - the matrix 2630 2631 Output Parameters: 2632 + m - the number of global rows 2633 - n - the number of global columns 2634 2635 .keywords: matrix, dimension, size, rows, columns, global, get 2636 2637 .seealso: MatGetLocalSize() 2638 @*/ 2639 int MatGetSize(Mat mat,int *m,int* n) 2640 { 2641 int ierr; 2642 2643 PetscFunctionBegin; 2644 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2645 ierr = (*mat->ops->getsize)(mat,m,n);CHKERRQ(ierr); 2646 PetscFunctionReturn(0); 2647 } 2648 2649 #undef __FUNC__ 2650 #define __FUNC__ "MatGetLocalSize" 2651 /*@ 2652 MatGetLocalSize - Returns the number of rows and columns in a matrix 2653 stored locally. This information may be implementation dependent, so 2654 use with care. 2655 2656 Not Collective 2657 2658 Input Parameters: 2659 . mat - the matrix 2660 2661 Output Parameters: 2662 + m - the number of local rows 2663 - n - the number of local columns 2664 2665 .keywords: matrix, dimension, size, local, rows, columns, get 2666 2667 .seealso: MatGetSize() 2668 @*/ 2669 int MatGetLocalSize(Mat mat,int *m,int* n) 2670 { 2671 int ierr; 2672 2673 PetscFunctionBegin; 2674 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2675 ierr = (*mat->ops->getlocalsize)(mat,m,n);CHKERRQ(ierr); 2676 PetscFunctionReturn(0); 2677 } 2678 2679 #undef __FUNC__ 2680 #define __FUNC__ "MatGetOwnershipRange" 2681 /*@ 2682 MatGetOwnershipRange - Returns the range of matrix rows owned by 2683 this processor, assuming that the matrix is laid out with the first 2684 n1 rows on the first processor, the next n2 rows on the second, etc. 2685 For certain parallel layouts this range may not be well defined. 2686 2687 Not Collective 2688 2689 Input Parameters: 2690 . mat - the matrix 2691 2692 Output Parameters: 2693 + m - the global index of the first local row 2694 - n - one more than the global index of the last local row 2695 2696 .keywords: matrix, get, range, ownership 2697 @*/ 2698 int MatGetOwnershipRange(Mat mat,int *m,int* n) 2699 { 2700 int ierr; 2701 2702 PetscFunctionBegin; 2703 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2704 PetscValidIntPointer(m); 2705 PetscValidIntPointer(n); 2706 if (!mat->ops->getownershiprange) SETERRQ(PETSC_ERR_SUP,0,""); 2707 ierr = (*mat->ops->getownershiprange)(mat,m,n);CHKERRQ(ierr); 2708 PetscFunctionReturn(0); 2709 } 2710 2711 #undef __FUNC__ 2712 #define __FUNC__ "MatILUFactorSymbolic" 2713 /*@ 2714 MatILUFactorSymbolic - Performs symbolic ILU factorization of a matrix. 2715 Uses levels of fill only, not drop tolerance. Use MatLUFactorNumeric() 2716 to complete the factorization. 2717 2718 Collective on Mat 2719 2720 Input Parameters: 2721 + mat - the matrix 2722 . row - row permutation 2723 . column - column permutation 2724 - info - structure containing 2725 $ levels - number of levels of fill. 2726 $ expected fill - as ratio of original fill. 2727 $ 1 or 0 - indicating force fill on diagonal (improves robustness for matrices 2728 missing diagonal entries) 2729 2730 Output Parameters: 2731 . fact - new matrix that has been symbolically factored 2732 2733 Notes: 2734 See the file ${PETSC_DIR}/Performace for additional information about 2735 choosing the fill factor for better efficiency. 2736 2737 .keywords: matrix, factor, incomplete, ILU, symbolic, fill 2738 2739 .seealso: MatLUFactorSymbolic(), MatLUFactorNumeric() 2740 @*/ 2741 int MatILUFactorSymbolic(Mat mat,IS row,IS col,MatILUInfo *info,Mat *fact) 2742 { 2743 int ierr; 2744 2745 PetscFunctionBegin; 2746 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2747 PetscValidPointer(fact); 2748 if (info && info->levels < 0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,0,"Levels of fill negative %d",info->levels); 2749 if (info && info->fill < 1.0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,0,"Expected fill less than 1.0 %g",info->fill); 2750 if (!mat->ops->ilufactorsymbolic) SETERRQ(PETSC_ERR_SUP,0,"Only MatCreateMPIRowbs() matrices support parallel ILU"); 2751 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2752 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2753 2754 PLogEventBegin(MAT_ILUFactorSymbolic,mat,row,col,0); 2755 ierr = (*mat->ops->ilufactorsymbolic)(mat,row,col,info,fact); CHKERRQ(ierr); 2756 PLogEventEnd(MAT_ILUFactorSymbolic,mat,row,col,0); 2757 PetscFunctionReturn(0); 2758 } 2759 2760 #undef __FUNC__ 2761 #define __FUNC__ "MatIncompleteCholeskyFactorSymbolic" 2762 /*@ 2763 MatIncompleteCholeskyFactorSymbolic - Performs symbolic incomplete 2764 Cholesky factorization for a symmetric matrix. Use 2765 MatCholeskyFactorNumeric() to complete the factorization. 2766 2767 Collective on Mat 2768 2769 Input Parameters: 2770 + mat - the matrix 2771 . perm - row and column permutation 2772 . fill - levels of fill 2773 - f - expected fill as ratio of original fill 2774 2775 Output Parameter: 2776 . fact - the factored matrix 2777 2778 Note: Currently only no-fill factorization is supported. 2779 2780 .keywords: matrix, factor, incomplete, ICC, Cholesky, symbolic, fill 2781 2782 .seealso: MatCholeskyFactorNumeric(), MatCholeskyFactor() 2783 @*/ 2784 int MatIncompleteCholeskyFactorSymbolic(Mat mat,IS perm,double f,int fill,Mat *fact) 2785 { 2786 int ierr; 2787 2788 PetscFunctionBegin; 2789 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2790 PetscValidPointer(fact); 2791 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2792 if (fill < 0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,0,"Fill negative %d",fill); 2793 if (f < 1.0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,0,"Expected fill less than 1.0 %g",f); 2794 if (!mat->ops->incompletecholeskyfactorsymbolic) SETERRQ(PETSC_ERR_SUP,0,"Currently only MatCreateMPIRowbs() matrices support ICC in parallel"); 2795 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2796 2797 PLogEventBegin(MAT_IncompleteCholeskyFactorSymbolic,mat,perm,0,0); 2798 ierr = (*mat->ops->incompletecholeskyfactorsymbolic)(mat,perm,f,fill,fact);CHKERRQ(ierr); 2799 PLogEventEnd(MAT_IncompleteCholeskyFactorSymbolic,mat,perm,0,0); 2800 PetscFunctionReturn(0); 2801 } 2802 2803 #undef __FUNC__ 2804 #define __FUNC__ "MatGetArray" 2805 /*@C 2806 MatGetArray - Returns a pointer to the element values in the matrix. 2807 The result of this routine is dependent on the underlying matrix data 2808 structure, and may not even work for certain matrix types. You MUST 2809 call MatRestoreArray() when you no longer need to access the array. 2810 2811 Not Collective 2812 2813 Input Parameter: 2814 . mat - the matrix 2815 2816 Output Parameter: 2817 . v - the location of the values 2818 2819 Currently returns an array only for the dense formats, giving access to 2820 the local portion of the matrix in the usual Fortran column-oriented format. 2821 2822 Fortran Note: 2823 This routine is used differently from Fortran 2824 $ Mat mat 2825 $ Scalar mat_array(1) 2826 $ PetscOffset i_mat 2827 $ int ierr 2828 $ call MatGetArray(mat,mat_array,i_mat,ierr) 2829 $ 2830 $ Access first local entry in matrix with 2831 $ value = mat_array(i_mat + 1) 2832 $ (note here array is treated as one dimensional) 2833 $ ...... other code 2834 $ call MatRestoreArray(mat,mat_array,i_mat,ierr) 2835 2836 See the Fortran chapter of the users manual and 2837 petsc/src/mat/examples/tests for details. 2838 2839 .keywords: matrix, array, elements, values 2840 2841 .seealso: MatRestoreArray() 2842 @*/ 2843 int MatGetArray(Mat mat,Scalar **v) 2844 { 2845 int ierr; 2846 2847 PetscFunctionBegin; 2848 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2849 PetscValidPointer(v); 2850 if (!mat->ops->getarray) SETERRQ(PETSC_ERR_SUP,0,""); 2851 ierr = (*mat->ops->getarray)(mat,v);CHKERRQ(ierr); 2852 PetscFunctionReturn(0); 2853 } 2854 2855 #undef __FUNC__ 2856 #define __FUNC__ "MatRestoreArray" 2857 /*@C 2858 MatRestoreArray - Restores the matrix after MatGetArray() has been called. 2859 2860 Not Collective 2861 2862 Input Parameter: 2863 + mat - the matrix 2864 - v - the location of the values 2865 2866 Fortran Note: 2867 This routine is used differently from Fortran 2868 $ Mat mat 2869 $ Scalar mat_array(1) 2870 $ PetscOffset i_mat 2871 $ int ierr 2872 $ call MatGetArray(mat,mat_array,i_mat,ierr) 2873 $ 2874 $ Access first local entry in matrix with 2875 $ value = mat_array(i_mat + 1) 2876 $ (note here the array is treated as one dimensional) 2877 $ 2878 $ ...... other code 2879 $ call MatRestoreArray(mat,mat_array,i_mat,ierr) 2880 2881 See the Fortran chapter of the users manual and 2882 petsc/src/mat/examples/tests for details 2883 2884 .keywords: matrix, array, elements, values, restore 2885 2886 .seealso: MatGetArray() 2887 @*/ 2888 int MatRestoreArray(Mat mat,Scalar **v) 2889 { 2890 int ierr; 2891 2892 PetscFunctionBegin; 2893 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2894 PetscValidPointer(v); 2895 if (!mat->ops->restorearray) SETERRQ(PETSC_ERR_SUP,0,""); 2896 ierr = (*mat->ops->restorearray)(mat,v);CHKERRQ(ierr); 2897 PetscFunctionReturn(0); 2898 } 2899 2900 #undef __FUNC__ 2901 #define __FUNC__ "MatGetSubMatrices" 2902 /*@C 2903 MatGetSubMatrices - Extracts several submatrices from a matrix. If submat 2904 points to an array of valid matrices, they may be reused to store the new 2905 submatrices. 2906 2907 Collective on Mat 2908 2909 Input Parameters: 2910 + mat - the matrix 2911 . n - the number of submatrixes to be extracted 2912 . irow, icol - index sets of rows and columns to extract 2913 - scall - either MAT_INITIAL_MATRIX or MAT_REUSE_MATRIX 2914 2915 Output Parameter: 2916 . submat - the array of submatrices 2917 2918 Notes: 2919 MatGetSubMatrices() can extract only sequential submatrices 2920 (from both sequential and parallel matrices). Use MatGetSubMatrix() 2921 to extract a parallel submatrix. 2922 2923 When extracting submatrices from a parallel matrix, each processor can 2924 form a different submatrix by setting the rows and columns of its 2925 individual index sets according to the local submatrix desired. 2926 2927 When finished using the submatrices, the user should destroy 2928 them with MatDestroySubMatrices(). 2929 2930 MAT_REUSE_MATRIX can only be used when the nonzero structure of the 2931 original matrix has not changed from that last call to MatGetSubMatrices() 2932 2933 Fortran Note: 2934 The Fortran interface is slightly different from that given below, it 2935 requires one to pass in as submat a Mat (integer) array of size at least m. 2936 2937 .keywords: matrix, get, submatrix, submatrices 2938 2939 .seealso: MatDestroyMatrices(), MatGetSubMatrix() 2940 @*/ 2941 int MatGetSubMatrices(Mat mat,int n,IS *irow,IS *icol,MatReuse scall,Mat **submat) 2942 { 2943 int ierr; 2944 2945 PetscFunctionBegin; 2946 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2947 if (!mat->ops->getsubmatrices) SETERRQ(PETSC_ERR_SUP,0,""); 2948 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2949 2950 PLogEventBegin(MAT_GetSubMatrices,mat,0,0,0); 2951 ierr = (*mat->ops->getsubmatrices)(mat,n,irow,icol,scall,submat); CHKERRQ(ierr); 2952 PLogEventEnd(MAT_GetSubMatrices,mat,0,0,0); 2953 2954 PetscFunctionReturn(0); 2955 } 2956 2957 #undef __FUNC__ 2958 #define __FUNC__ "MatDestroyMatrices" 2959 /*@C 2960 MatDestroyMatrices - Destroys a set of matrices obtained with MatGetSubMatrices(). 2961 2962 Collective on Mat 2963 2964 Input Parameters: 2965 + n - the number of local matrices 2966 - mat - the matrices 2967 2968 .keywords: matrix, destroy, submatrix, submatrices 2969 2970 .seealso: MatGetSubMatrices() 2971 @*/ 2972 int MatDestroyMatrices(int n,Mat **mat) 2973 { 2974 int ierr,i; 2975 2976 PetscFunctionBegin; 2977 if (n < 0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,1,"Trying to destroy negative number of matrices %d",n); 2978 PetscValidPointer(mat); 2979 for ( i=0; i<n; i++ ) { 2980 ierr = MatDestroy((*mat)[i]); CHKERRQ(ierr); 2981 } 2982 if (n) PetscFree(*mat); 2983 PetscFunctionReturn(0); 2984 } 2985 2986 #undef __FUNC__ 2987 #define __FUNC__ "MatIncreaseOverlap" 2988 /*@ 2989 MatIncreaseOverlap - Given a set of submatrices indicated by index sets, 2990 replaces the index sets by larger ones that represent submatrices with 2991 additional overlap. 2992 2993 Collective on Mat 2994 2995 Input Parameters: 2996 + mat - the matrix 2997 . n - the number of index sets 2998 . is - the array of pointers to index sets 2999 - ov - the additional overlap requested 3000 3001 .keywords: matrix, overlap, Schwarz 3002 3003 .seealso: MatGetSubMatrices() 3004 @*/ 3005 int MatIncreaseOverlap(Mat mat,int n, IS *is,int ov) 3006 { 3007 int ierr; 3008 3009 PetscFunctionBegin; 3010 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3011 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 3012 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 3013 3014 if (ov == 0) PetscFunctionReturn(0); 3015 if (!mat->ops->increaseoverlap) SETERRQ(PETSC_ERR_SUP,0,""); 3016 PLogEventBegin(MAT_IncreaseOverlap,mat,0,0,0); 3017 ierr = (*mat->ops->increaseoverlap)(mat,n,is,ov); CHKERRQ(ierr); 3018 PLogEventEnd(MAT_IncreaseOverlap,mat,0,0,0); 3019 PetscFunctionReturn(0); 3020 } 3021 3022 #undef __FUNC__ 3023 #define __FUNC__ "MatPrintHelp" 3024 /*@ 3025 MatPrintHelp - Prints all the options for the matrix. 3026 3027 Collective on Mat 3028 3029 Input Parameter: 3030 . mat - the matrix 3031 3032 Options Database Keys: 3033 + -help - Prints matrix options 3034 - -h - Prints matrix options 3035 3036 .keywords: mat, help 3037 3038 .seealso: MatCreate(), MatCreateXXX() 3039 @*/ 3040 int MatPrintHelp(Mat mat) 3041 { 3042 static int called = 0; 3043 int ierr; 3044 MPI_Comm comm; 3045 3046 PetscFunctionBegin; 3047 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3048 3049 comm = mat->comm; 3050 if (!called) { 3051 (*PetscHelpPrintf)(comm,"General matrix options:\n"); 3052 (*PetscHelpPrintf)(comm," -mat_view_info: view basic matrix info during MatAssemblyEnd()\n"); 3053 (*PetscHelpPrintf)(comm," -mat_view_info_detailed: view detailed matrix info during MatAssemblyEnd()\n"); 3054 (*PetscHelpPrintf)(comm," -mat_view_draw: draw nonzero matrix structure during MatAssemblyEnd()\n"); 3055 (*PetscHelpPrintf)(comm," -draw_pause <sec>: set seconds of display pause\n"); 3056 (*PetscHelpPrintf)(comm," -display <name>: set alternate display\n"); 3057 called = 1; 3058 } 3059 if (mat->ops->printhelp) { 3060 ierr = (*mat->ops->printhelp)(mat); CHKERRQ(ierr); 3061 } 3062 PetscFunctionReturn(0); 3063 } 3064 3065 #undef __FUNC__ 3066 #define __FUNC__ "MatGetBlockSize" 3067 /*@ 3068 MatGetBlockSize - Returns the matrix block size; useful especially for the 3069 block row and block diagonal formats. 3070 3071 Not Collective 3072 3073 Input Parameter: 3074 . mat - the matrix 3075 3076 Output Parameter: 3077 . bs - block size 3078 3079 Notes: 3080 Block diagonal formats are MATSEQBDIAG, MATMPIBDIAG. 3081 Block row formats are MATSEQBAIJ, MATMPIBAIJ 3082 3083 .keywords: matrix, get, block, size 3084 3085 .seealso: MatCreateSeqBAIJ(), MatCreateMPIBAIJ(), MatCreateSeqBDiag(), MatCreateMPIBDiag() 3086 @*/ 3087 int MatGetBlockSize(Mat mat,int *bs) 3088 { 3089 int ierr; 3090 3091 PetscFunctionBegin; 3092 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3093 PetscValidIntPointer(bs); 3094 if (!mat->ops->getblocksize) SETERRQ(PETSC_ERR_SUP,0,""); 3095 ierr = (*mat->ops->getblocksize)(mat,bs);CHKERRQ(ierr); 3096 PetscFunctionReturn(0); 3097 } 3098 3099 #undef __FUNC__ 3100 #define __FUNC__ "MatGetRowIJ" 3101 /*@C 3102 MatGetRowIJ - Returns the compressed row storage i and j indices for sequential matrices. 3103 EXPERTS ONLY. 3104 3105 Collective on Mat 3106 3107 Input Parameters: 3108 + mat - the matrix 3109 . shift - 0 or 1 indicating we want the indices starting at 0 or 1 3110 - symmetric - PETSC_TRUE or PETSC_FALSE indicating the matrix data structure should be 3111 symmetrized 3112 3113 Output Parameters: 3114 + n - number of rows in the (possibly compressed) matrix 3115 . ia - the row pointers 3116 . ja - the column indices 3117 - done - PETSC_TRUE or PETSC_FALSE, indicating whether the values have been returned 3118 3119 .seealso: MatGetColumnIJ(), MatRestoreRowIJ() 3120 @*/ 3121 int MatGetRowIJ(Mat mat,int shift,PetscTruth symmetric,int *n,int **ia,int** ja,PetscTruth *done) 3122 { 3123 int ierr; 3124 3125 PetscFunctionBegin; 3126 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3127 if (ia) PetscValidIntPointer(ia); 3128 if (ja) PetscValidIntPointer(ja); 3129 PetscValidIntPointer(done); 3130 if (!mat->ops->getrowij) *done = PETSC_FALSE; 3131 else { 3132 *done = PETSC_TRUE; 3133 ierr = (*mat->ops->getrowij)(mat,shift,symmetric,n,ia,ja,done); CHKERRQ(ierr); 3134 } 3135 PetscFunctionReturn(0); 3136 } 3137 3138 #undef __FUNC__ 3139 #define __FUNC__ "MatGetColumnIJ" 3140 /*@C 3141 MatGetColumnIJ - Returns the compressed column storage i and j indices for sequential matrices. 3142 EXPERTS ONLY. 3143 3144 Collective on Mat 3145 3146 Input Parameters: 3147 + mat - the matrix 3148 . shift - 1 or zero indicating we want the indices starting at 0 or 1 3149 - symmetric - PETSC_TRUE or PETSC_FALSE indicating the matrix data structure should be 3150 symmetrized 3151 3152 Output Parameters: 3153 + n - number of columns in the (possibly compressed) matrix 3154 . ia - the column pointers 3155 . ja - the row indices 3156 - done - PETSC_TRUE or PETSC_FALSE, indicating whether the values have been returned 3157 3158 .seealso: MatGetRowIJ(), MatRestoreColumnIJ() 3159 @*/ 3160 int MatGetColumnIJ(Mat mat,int shift,PetscTruth symmetric,int *n,int **ia,int** ja,PetscTruth *done) 3161 { 3162 int ierr; 3163 3164 PetscFunctionBegin; 3165 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3166 if (ia) PetscValidIntPointer(ia); 3167 if (ja) PetscValidIntPointer(ja); 3168 PetscValidIntPointer(done); 3169 3170 if (!mat->ops->getcolumnij) *done = PETSC_FALSE; 3171 else { 3172 *done = PETSC_TRUE; 3173 ierr = (*mat->ops->getcolumnij)(mat,shift,symmetric,n,ia,ja,done); CHKERRQ(ierr); 3174 } 3175 PetscFunctionReturn(0); 3176 } 3177 3178 #undef __FUNC__ 3179 #define __FUNC__ "MatRestoreRowIJ" 3180 /*@C 3181 MatRestoreRowIJ - Call after you are completed with the ia,ja indices obtained with 3182 MatGetRowIJ(). EXPERTS ONLY. 3183 3184 Collective on Mat 3185 3186 Input Parameters: 3187 + mat - the matrix 3188 . shift - 1 or zero indicating we want the indices starting at 0 or 1 3189 - symmetric - PETSC_TRUE or PETSC_FALSE indicating the matrix data structure should be 3190 symmetrized 3191 3192 Output Parameters: 3193 + n - size of (possibly compressed) matrix 3194 . ia - the row pointers 3195 . ja - the column indices 3196 - done - PETSC_TRUE or PETSC_FALSE indicated that the values have been returned 3197 3198 .seealso: MatGetRowIJ(), MatRestoreColumnIJ() 3199 @*/ 3200 int MatRestoreRowIJ(Mat mat,int shift,PetscTruth symmetric,int *n,int **ia,int** ja,PetscTruth *done) 3201 { 3202 int ierr; 3203 3204 PetscFunctionBegin; 3205 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3206 if (ia) PetscValidIntPointer(ia); 3207 if (ja) PetscValidIntPointer(ja); 3208 PetscValidIntPointer(done); 3209 3210 if (!mat->ops->restorerowij) *done = PETSC_FALSE; 3211 else { 3212 *done = PETSC_TRUE; 3213 ierr = (*mat->ops->restorerowij)(mat,shift,symmetric,n,ia,ja,done); CHKERRQ(ierr); 3214 } 3215 PetscFunctionReturn(0); 3216 } 3217 3218 #undef __FUNC__ 3219 #define __FUNC__ "MatRestoreColumnIJ" 3220 /*@C 3221 MatRestoreColumnIJ - Call after you are completed with the ia,ja indices obtained with 3222 MatGetColumnIJ(). EXPERTS ONLY. 3223 3224 Collective on Mat 3225 3226 Input Parameters: 3227 + mat - the matrix 3228 . shift - 1 or zero indicating we want the indices starting at 0 or 1 3229 - symmetric - PETSC_TRUE or PETSC_FALSE indicating the matrix data structure should be 3230 symmetrized 3231 3232 Output Parameters: 3233 + n - size of (possibly compressed) matrix 3234 . ia - the column pointers 3235 . ja - the row indices 3236 - done - PETSC_TRUE or PETSC_FALSE indicated that the values have been returned 3237 3238 .seealso: MatGetColumnIJ(), MatRestoreRowIJ() 3239 @*/ 3240 int MatRestoreColumnIJ(Mat mat,int shift,PetscTruth symmetric,int *n,int **ia,int** ja,PetscTruth *done) 3241 { 3242 int ierr; 3243 3244 PetscFunctionBegin; 3245 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3246 if (ia) PetscValidIntPointer(ia); 3247 if (ja) PetscValidIntPointer(ja); 3248 PetscValidIntPointer(done); 3249 3250 if (!mat->ops->restorecolumnij) *done = PETSC_FALSE; 3251 else { 3252 *done = PETSC_TRUE; 3253 ierr = (*mat->ops->restorecolumnij)(mat,shift,symmetric,n,ia,ja,done); CHKERRQ(ierr); 3254 } 3255 PetscFunctionReturn(0); 3256 } 3257 3258 #undef __FUNC__ 3259 #define __FUNC__ "MatColoringPatch" 3260 /*@C 3261 MatColoringPatch - EXPERTS ONLY, used inside matrix coloring routines that 3262 use MatGetRowIJ() and/or MatGetColumnIJ(). 3263 3264 Collective on Mat 3265 3266 Input Parameters: 3267 + mat - the matrix 3268 . n - number of colors 3269 - colorarray - array indicating color for each column 3270 3271 Output Parameters: 3272 . iscoloring - coloring generated using colorarray information 3273 3274 .seealso: MatGetRowIJ(), MatGetColumnIJ() 3275 3276 .keywords: mat, coloring, patch 3277 @*/ 3278 int MatColoringPatch(Mat mat,int n,int *colorarray,ISColoring *iscoloring) 3279 { 3280 int ierr; 3281 3282 PetscFunctionBegin; 3283 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3284 PetscValidIntPointer(colorarray); 3285 3286 if (!mat->ops->coloringpatch) {SETERRQ(PETSC_ERR_SUP,0,"");} 3287 else { 3288 ierr = (*mat->ops->coloringpatch)(mat,n,colorarray,iscoloring); CHKERRQ(ierr); 3289 } 3290 PetscFunctionReturn(0); 3291 } 3292 3293 3294 #undef __FUNC__ 3295 #define __FUNC__ "MatSetUnfactored" 3296 /*@ 3297 MatSetUnfactored - Resets a factored matrix to be treated as unfactored. 3298 3299 Collective on Mat 3300 3301 Input Parameter: 3302 . mat - the factored matrix to be reset 3303 3304 Notes: 3305 This routine should be used only with factored matrices formed by in-place 3306 factorization via ILU(0) (or by in-place LU factorization for the MATSEQDENSE 3307 format). This option can save memory, for example, when solving nonlinear 3308 systems with a matrix-free Newton-Krylov method and a matrix-based, in-place 3309 ILU(0) preconditioner. 3310 3311 Note that one can specify in-place ILU(0) factorization by calling 3312 $ PCType(pc,PCILU); 3313 $ PCILUSeUseInPlace(pc); 3314 or by using the options -pc_type ilu -pc_ilu_in_place 3315 3316 In-place factorization ILU(0) can also be used as a local 3317 solver for the blocks within the block Jacobi or additive Schwarz 3318 methods (runtime option: -sub_pc_ilu_in_place). See the discussion 3319 of these preconditioners in the users manual for details on setting 3320 local solver options. 3321 3322 .seealso: PCILUSetUseInPlace(), PCLUSetUseInPlace() 3323 3324 .keywords: matrix-free, in-place ILU, in-place LU 3325 @*/ 3326 int MatSetUnfactored(Mat mat) 3327 { 3328 int ierr; 3329 3330 PetscFunctionBegin; 3331 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3332 mat->factor = 0; 3333 if (!mat->ops->setunfactored) PetscFunctionReturn(0); 3334 ierr = (*mat->ops->setunfactored)(mat); CHKERRQ(ierr); 3335 PetscFunctionReturn(0); 3336 } 3337 3338 #undef __FUNC__ 3339 #define __FUNC__ "MatGetType" 3340 /*@C 3341 MatGetType - Gets the matrix type and name (as a string) from the matrix. 3342 3343 Not Collective 3344 3345 Input Parameter: 3346 . mat - the matrix 3347 3348 Output Parameter: 3349 + type - the matrix type (or use PETSC_NULL) 3350 - name - name of matrix type (or use PETSC_NULL) 3351 3352 .keywords: matrix, get, type, name 3353 @*/ 3354 int MatGetType(Mat mat,MatType *type,char **name) 3355 { 3356 int itype = (int)mat->type; 3357 char *matname[10]; 3358 3359 PetscFunctionBegin; 3360 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3361 3362 if (type) *type = (MatType) mat->type; 3363 if (name) { 3364 /* Note: Be sure that this list corresponds to the enum in mat.h */ 3365 matname[0] = "MATSEQDENSE"; 3366 matname[1] = "MATSEQAIJ"; 3367 matname[2] = "MATMPIAIJ"; 3368 matname[3] = "MATSHELL"; 3369 matname[4] = "MATMPIROWBS"; 3370 matname[5] = "MATSEQBDIAG"; 3371 matname[6] = "MATMPIBDIAG"; 3372 matname[7] = "MATMPIDENSE"; 3373 matname[8] = "MATSEQBAIJ"; 3374 matname[9] = "MATMPIBAIJ"; 3375 3376 if (itype < 0 || itype > 9) *name = "Unknown matrix type"; 3377 else *name = matname[itype]; 3378 } 3379 PetscFunctionReturn(0); 3380 } 3381 3382 /*MC 3383 MatGetArrayF90 - Accesses a matrix array from Fortran90. 3384 3385 Synopsis: 3386 MatGetArrayF90(Mat x,{Scalar, pointer :: xx_v(:)},integer ierr) 3387 3388 Not collective 3389 3390 Input Parameter: 3391 . x - matrix 3392 3393 Output Parameters: 3394 + xx_v - the Fortran90 pointer to the array 3395 - ierr - error code 3396 3397 Example of Usage: 3398 .vb 3399 Scalar, pointer xx_v(:) 3400 .... 3401 call MatGetArrayF90(x,xx_v,ierr) 3402 a = xx_v(3) 3403 call MatRestoreArrayF90(x,xx_v,ierr) 3404 .ve 3405 3406 Notes: 3407 Not yet supported for all F90 compilers 3408 3409 .seealso: MatRestoreArrayF90(), MatGetArray(), MatRestoreArray() 3410 3411 .keywords: matrix, array, f90 3412 M*/ 3413 3414 /*MC 3415 MatRestoreArrayF90 - Restores a matrix array that has been 3416 accessed with MatGetArrayF90(). 3417 3418 Synopsis: 3419 MatRestoreArrayF90(Mat x,{Scalar, pointer :: xx_v(:)},integer ierr) 3420 3421 Not collective 3422 3423 Input Parameters: 3424 + x - matrix 3425 - xx_v - the Fortran90 pointer to the array 3426 3427 Output Parameter: 3428 . ierr - error code 3429 3430 Example of Usage: 3431 .vb 3432 Scalar, pointer xx_v(:) 3433 .... 3434 call MatGetArrayF90(x,xx_v,ierr) 3435 a = xx_v(3) 3436 call MatRestoreArrayF90(x,xx_v,ierr) 3437 .ve 3438 3439 Notes: 3440 Not yet supported for all F90 compilers 3441 3442 .seealso: MatGetArrayF90(), MatGetArray(), MatRestoreArray() 3443 3444 .keywords: matrix, array, f90 3445 M*/ 3446 3447 3448 #undef __FUNC__ 3449 #define __FUNC__ "MatGetSubMatrix" 3450 /*@ 3451 MatGetSubMatrix - Gets a single submatrix on the same number of processors 3452 as the original matrix. 3453 3454 Collective on Mat 3455 3456 Input Parameters: 3457 + mat - the original matrix 3458 . isrow - rows this processor should obtain 3459 . iscol - columns for all processors you wish to keep 3460 . csize - number of columns "local" to this processor (does nothing for sequential 3461 matrices). This should match the result from VecGetLocalSize(x,...) if you 3462 plan to use the matrix in a A*x or use PETSC_DECIDE 3463 - cll - either MAT_INITIAL_MATRIX or MAT_REUSE_MATRIX 3464 3465 Output Parameter: 3466 . newmat - the new submatrix, of the same type as the old 3467 3468 .seealso: MatGetSubMatrices() 3469 3470 @*/ 3471 int MatGetSubMatrix(Mat mat,IS isrow,IS iscol,int csize,MatReuse cll,Mat *newmat) 3472 { 3473 int ierr, size; 3474 Mat *local; 3475 3476 PetscFunctionBegin; 3477 MPI_Comm_size(mat->comm,&size); 3478 3479 /* if original matrix is on just one processor then use submatrix generated */ 3480 if (!mat->ops->getsubmatrix && size == 1 && cll == MAT_REUSE_MATRIX) { 3481 ierr = MatGetSubMatrices(mat,1,&isrow,&iscol,MAT_REUSE_MATRIX,&newmat);CHKERRQ(ierr); 3482 PetscFunctionReturn(0); 3483 } else if (!mat->ops->getsubmatrix && size == 1) { 3484 ierr = MatGetSubMatrices(mat,1,&isrow,&iscol,MAT_INITIAL_MATRIX,&local);CHKERRQ(ierr); 3485 *newmat = *local; 3486 PetscFree(local); 3487 PetscFunctionReturn(0); 3488 } 3489 3490 if (!mat->ops->getsubmatrix) SETERRQ(PETSC_ERR_SUP,0,"Not currently implemented"); 3491 ierr = (*mat->ops->getsubmatrix)(mat,isrow,iscol,csize,cll,newmat);CHKERRQ(ierr); 3492 PetscFunctionReturn(0); 3493 } 3494 3495 #undef __FUNC__ 3496 #define __FUNC__ "MatGetMaps" 3497 /*@C 3498 MatGetMaps - Returns the maps associated with the matrix. 3499 3500 Not Collective 3501 3502 Input Parameter: 3503 . mat - the matrix 3504 3505 Output Parameters: 3506 + rmap - the row (right) map 3507 - cmap - the column (left) map 3508 3509 .keywords: matrix, get, map 3510 @*/ 3511 int MatGetMaps(Mat mat,Map *rmap,Map *cmap) 3512 { 3513 int ierr; 3514 3515 PetscFunctionBegin; 3516 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3517 ierr = (*mat->ops->getmaps)(mat,rmap,cmap);CHKERRQ(ierr); 3518 PetscFunctionReturn(0); 3519 } 3520 3521 /* 3522 Version that works for all PETSc matrices 3523 */ 3524 #undef __FUNC__ 3525 #define __FUNC__ "MatGetMaps_Petsc" 3526 int MatGetMaps_Petsc(Mat mat,Map *rmap,Map *cmap) 3527 { 3528 PetscFunctionBegin; 3529 if (rmap) *rmap = mat->rmap; 3530 if (cmap) *cmap = mat->cmap; 3531 PetscFunctionReturn(0); 3532 } 3533