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