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