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