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