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