1 /*$Id: matrix.c,v 1.353 1999/11/10 03:19:04 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 Level: intermediate 2742 2743 .keywords: matrix, option, row-oriented, column-oriented, sorted, nonzero 2744 @*/ 2745 int MatSetOption(Mat mat,MatOption op) 2746 { 2747 int ierr; 2748 2749 PetscFunctionBegin; 2750 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2751 if (mat->ops->setoption) {ierr = (*mat->ops->setoption)(mat,op);CHKERRQ(ierr);} 2752 PetscFunctionReturn(0); 2753 } 2754 2755 #undef __FUNC__ 2756 #define __FUNC__ "MatZeroEntries" 2757 /*@ 2758 MatZeroEntries - Zeros all entries of a matrix. For sparse matrices 2759 this routine retains the old nonzero structure. 2760 2761 Collective on Mat 2762 2763 Input Parameters: 2764 . mat - the matrix 2765 2766 Level: intermediate 2767 2768 .keywords: matrix, zero, entries 2769 2770 .seealso: MatZeroRows() 2771 @*/ 2772 int MatZeroEntries(Mat mat) 2773 { 2774 int ierr; 2775 2776 PetscFunctionBegin; 2777 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2778 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2779 if (!mat->ops->zeroentries) SETERRQ(PETSC_ERR_SUP,0,""); 2780 2781 PLogEventBegin(MAT_ZeroEntries,mat,0,0,0); 2782 ierr = (*mat->ops->zeroentries)(mat);CHKERRQ(ierr); 2783 PLogEventEnd(MAT_ZeroEntries,mat,0,0,0); 2784 PetscFunctionReturn(0); 2785 } 2786 2787 #undef __FUNC__ 2788 #define __FUNC__ "MatZeroRows" 2789 /*@C 2790 MatZeroRows - Zeros all entries (except possibly the main diagonal) 2791 of a set of rows of a matrix. 2792 2793 Collective on Mat 2794 2795 Input Parameters: 2796 + mat - the matrix 2797 . is - index set of rows to remove 2798 - diag - pointer to value put in all diagonals of eliminated rows. 2799 Note that diag is not a pointer to an array, but merely a 2800 pointer to a single value. 2801 2802 Notes: 2803 For the AIJ and BAIJ matrix formats this removes the old nonzero structure, 2804 but does not release memory. For the dense and block diagonal 2805 formats this does not alter the nonzero structure. 2806 2807 If the option MatSetOption(mat,MAT_KEEP_ZEROED_ROWS) the nonzero structure 2808 of the matrix is not changed (even for AIJ and BAIJ matrices) the values are 2809 merely zeroed. 2810 2811 The user can set a value in the diagonal entry (or for the AIJ and 2812 row formats can optionally remove the main diagonal entry from the 2813 nonzero structure as well, by passing a null pointer (PETSC_NULL 2814 in C or PETSC_NULL_SCALAR in Fortran) as the final argument). 2815 2816 For the parallel case, all processes that share the matrix (i.e., 2817 those in the communicator used for matrix creation) MUST call this 2818 routine, regardless of whether any rows being zeroed are owned by 2819 them. 2820 2821 2822 Level: intermediate 2823 2824 .keywords: matrix, zero, rows, boundary conditions 2825 2826 .seealso: MatZeroEntries(), MatZeroRowsLocal(), MatSetOption() 2827 @*/ 2828 int MatZeroRows(Mat mat,IS is, Scalar *diag) 2829 { 2830 int ierr; 2831 2832 PetscFunctionBegin; 2833 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2834 PetscValidHeaderSpecific(is,IS_COOKIE); 2835 if (diag) PetscValidScalarPointer(diag); 2836 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2837 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2838 if (!mat->ops->zerorows) SETERRQ(PETSC_ERR_SUP,0,""); 2839 2840 ierr = (*mat->ops->zerorows)(mat,is,diag);CHKERRQ(ierr); 2841 ierr = MatView_Private(mat);CHKERRQ(ierr); 2842 PetscFunctionReturn(0); 2843 } 2844 2845 #undef __FUNC__ 2846 #define __FUNC__ "MatZeroRowsLocal" 2847 /*@C 2848 MatZeroRowsLocal - Zeros all entries (except possibly the main diagonal) 2849 of a set of rows of a matrix; using local numbering of rows. 2850 2851 Collective on Mat 2852 2853 Input Parameters: 2854 + mat - the matrix 2855 . is - index set of rows to remove 2856 - diag - pointer to value put in all diagonals of eliminated rows. 2857 Note that diag is not a pointer to an array, but merely a 2858 pointer to a single value. 2859 2860 Notes: 2861 For the AIJ matrix formats this removes the old nonzero structure, 2862 but does not release memory. For the dense and block diagonal 2863 formats this does not alter the nonzero structure. 2864 2865 The user can set a value in the diagonal entry (or for the AIJ and 2866 row formats can optionally remove the main diagonal entry from the 2867 nonzero structure as well, by passing a null pointer (PETSC_NULL 2868 in C or PETSC_NULL_SCALAR in Fortran) as the final argument). 2869 2870 Level: intermediate 2871 2872 .keywords: matrix, zero, rows, boundary conditions 2873 2874 .seealso: MatZeroEntries(), MatZeroRows() 2875 @*/ 2876 int MatZeroRowsLocal(Mat mat,IS is, Scalar *diag) 2877 { 2878 int ierr; 2879 IS newis; 2880 2881 PetscFunctionBegin; 2882 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2883 PetscValidHeaderSpecific(is,IS_COOKIE); 2884 if (diag) PetscValidScalarPointer(diag); 2885 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2886 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2887 if (!mat->ops->zerorows) SETERRQ(PETSC_ERR_SUP,0,""); 2888 2889 ierr = ISLocalToGlobalMappingApplyIS(mat->mapping,is,&newis);CHKERRQ(ierr); 2890 ierr = (*mat->ops->zerorows)(mat,newis,diag);CHKERRQ(ierr); 2891 ierr = ISDestroy(newis);CHKERRQ(ierr); 2892 PetscFunctionReturn(0); 2893 } 2894 2895 #undef __FUNC__ 2896 #define __FUNC__ "MatGetSize" 2897 /*@ 2898 MatGetSize - Returns the numbers of rows and columns in a matrix. 2899 2900 Not Collective 2901 2902 Input Parameter: 2903 . mat - the matrix 2904 2905 Output Parameters: 2906 + m - the number of global rows 2907 - n - the number of global columns 2908 2909 Level: beginner 2910 2911 .keywords: matrix, dimension, size, rows, columns, global, get 2912 2913 .seealso: MatGetLocalSize() 2914 @*/ 2915 int MatGetSize(Mat mat,int *m,int* n) 2916 { 2917 int ierr; 2918 2919 PetscFunctionBegin; 2920 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2921 ierr = (*mat->ops->getsize)(mat,m,n);CHKERRQ(ierr); 2922 PetscFunctionReturn(0); 2923 } 2924 2925 #undef __FUNC__ 2926 #define __FUNC__ "MatGetLocalSize" 2927 /*@ 2928 MatGetLocalSize - Returns the number of rows and columns in a matrix 2929 stored locally. This information may be implementation dependent, so 2930 use with care. 2931 2932 Not Collective 2933 2934 Input Parameters: 2935 . mat - the matrix 2936 2937 Output Parameters: 2938 + m - the number of local rows 2939 - n - the number of local columns 2940 2941 Level: beginner 2942 2943 .keywords: matrix, dimension, size, local, rows, columns, get 2944 2945 .seealso: MatGetSize() 2946 @*/ 2947 int MatGetLocalSize(Mat mat,int *m,int* n) 2948 { 2949 int ierr; 2950 2951 PetscFunctionBegin; 2952 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2953 ierr = (*mat->ops->getlocalsize)(mat,m,n);CHKERRQ(ierr); 2954 PetscFunctionReturn(0); 2955 } 2956 2957 #undef __FUNC__ 2958 #define __FUNC__ "MatGetOwnershipRange" 2959 /*@ 2960 MatGetOwnershipRange - Returns the range of matrix rows owned by 2961 this processor, assuming that the matrix is laid out with the first 2962 n1 rows on the first processor, the next n2 rows on the second, etc. 2963 For certain parallel layouts this range may not be well defined. 2964 2965 Not Collective 2966 2967 Input Parameters: 2968 . mat - the matrix 2969 2970 Output Parameters: 2971 + m - the global index of the first local row 2972 - n - one more than the global index of the last local row 2973 2974 Level: beginner 2975 2976 .keywords: matrix, get, range, ownership 2977 @*/ 2978 int MatGetOwnershipRange(Mat mat,int *m,int* n) 2979 { 2980 int ierr; 2981 2982 PetscFunctionBegin; 2983 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2984 PetscValidIntPointer(m); 2985 PetscValidIntPointer(n); 2986 if (!mat->ops->getownershiprange) SETERRQ(PETSC_ERR_SUP,0,""); 2987 ierr = (*mat->ops->getownershiprange)(mat,m,n);CHKERRQ(ierr); 2988 PetscFunctionReturn(0); 2989 } 2990 2991 #undef __FUNC__ 2992 #define __FUNC__ "MatILUFactorSymbolic" 2993 /*@ 2994 MatILUFactorSymbolic - Performs symbolic ILU factorization of a matrix. 2995 Uses levels of fill only, not drop tolerance. Use MatLUFactorNumeric() 2996 to complete the factorization. 2997 2998 Collective on Mat 2999 3000 Input Parameters: 3001 + mat - the matrix 3002 . row - row permutation 3003 . column - column permutation 3004 - info - structure containing 3005 $ levels - number of levels of fill. 3006 $ expected fill - as ratio of original fill. 3007 $ 1 or 0 - indicating force fill on diagonal (improves robustness for matrices 3008 missing diagonal entries) 3009 3010 Output Parameters: 3011 . fact - new matrix that has been symbolically factored 3012 3013 Notes: 3014 See the users manual for additional information about 3015 choosing the fill factor for better efficiency. 3016 3017 Most users should employ the simplified SLES interface for linear solvers 3018 instead of working directly with matrix algebra routines such as this. 3019 See, e.g., SLESCreate(). 3020 3021 Level: developer 3022 3023 .keywords: matrix, factor, incomplete, ILU, symbolic, fill 3024 3025 .seealso: MatLUFactorSymbolic(), MatLUFactorNumeric(), MatCholeskyFactor() 3026 MatGetOrdering() 3027 3028 @*/ 3029 int MatILUFactorSymbolic(Mat mat,IS row,IS col,MatILUInfo *info,Mat *fact) 3030 { 3031 int ierr; 3032 3033 PetscFunctionBegin; 3034 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3035 PetscValidPointer(fact); 3036 if (info && info->levels < 0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,0,"Levels of fill negative %d",info->levels); 3037 if (info && info->fill < 1.0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,0,"Expected fill less than 1.0 %g",info->fill); 3038 if (!mat->ops->ilufactorsymbolic) SETERRQ(PETSC_ERR_SUP,0,"Only MatCreateMPIRowbs() matrices support parallel ILU"); 3039 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 3040 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 3041 3042 PLogEventBegin(MAT_ILUFactorSymbolic,mat,row,col,0); 3043 ierr = (*mat->ops->ilufactorsymbolic)(mat,row,col,info,fact);CHKERRQ(ierr); 3044 PLogEventEnd(MAT_ILUFactorSymbolic,mat,row,col,0); 3045 PetscFunctionReturn(0); 3046 } 3047 3048 #undef __FUNC__ 3049 #define __FUNC__ "MatIncompleteCholeskyFactorSymbolic" 3050 /*@ 3051 MatIncompleteCholeskyFactorSymbolic - Performs symbolic incomplete 3052 Cholesky factorization for a symmetric matrix. Use 3053 MatCholeskyFactorNumeric() to complete the factorization. 3054 3055 Collective on Mat 3056 3057 Input Parameters: 3058 + mat - the matrix 3059 . perm - row and column permutation 3060 . fill - levels of fill 3061 - f - expected fill as ratio of original fill 3062 3063 Output Parameter: 3064 . fact - the factored matrix 3065 3066 Notes: 3067 Currently only no-fill factorization is supported. 3068 3069 Most users should employ the simplified SLES interface for linear solvers 3070 instead of working directly with matrix algebra routines such as this. 3071 See, e.g., SLESCreate(). 3072 3073 Level: developer 3074 3075 .keywords: matrix, factor, incomplete, ICC, Cholesky, symbolic, fill 3076 3077 .seealso: MatCholeskyFactorNumeric(), MatCholeskyFactor() 3078 @*/ 3079 int MatIncompleteCholeskyFactorSymbolic(Mat mat,IS perm,double f,int fill,Mat *fact) 3080 { 3081 int ierr; 3082 3083 PetscFunctionBegin; 3084 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3085 PetscValidPointer(fact); 3086 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 3087 if (fill < 0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,0,"Fill negative %d",fill); 3088 if (f < 1.0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,0,"Expected fill less than 1.0 %g",f); 3089 if (!mat->ops->incompletecholeskyfactorsymbolic) SETERRQ(PETSC_ERR_SUP,0,"Currently only MatCreateMPIRowbs() matrices support ICC in parallel"); 3090 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 3091 3092 PLogEventBegin(MAT_IncompleteCholeskyFactorSymbolic,mat,perm,0,0); 3093 ierr = (*mat->ops->incompletecholeskyfactorsymbolic)(mat,perm,f,fill,fact);CHKERRQ(ierr); 3094 PLogEventEnd(MAT_IncompleteCholeskyFactorSymbolic,mat,perm,0,0); 3095 PetscFunctionReturn(0); 3096 } 3097 3098 #undef __FUNC__ 3099 #define __FUNC__ "MatGetArray" 3100 /*@C 3101 MatGetArray - Returns a pointer to the element values in the matrix. 3102 The result of this routine is dependent on the underlying matrix data 3103 structure, and may not even work for certain matrix types. You MUST 3104 call MatRestoreArray() when you no longer need to access the array. 3105 3106 Not Collective 3107 3108 Input Parameter: 3109 . mat - the matrix 3110 3111 Output Parameter: 3112 . v - the location of the values 3113 3114 Currently returns an array only for the dense formats, giving access to 3115 the local portion of the matrix in the usual Fortran column-oriented format. 3116 3117 Fortran Note: 3118 This routine is used differently from Fortran, e.g., 3119 .vb 3120 Mat mat 3121 Scalar mat_array(1) 3122 PetscOffset i_mat 3123 int ierr 3124 call MatGetArray(mat,mat_array,i_mat,ierr) 3125 3126 C Access first local entry in matrix; note that array is 3127 C treated as one dimensional 3128 value = mat_array(i_mat + 1) 3129 3130 [... other code ...] 3131 call MatRestoreArray(mat,mat_array,i_mat,ierr) 3132 .ve 3133 3134 See the Fortran chapter of the users manual and 3135 petsc/src/mat/examples/tests for details. 3136 3137 Level: advanced 3138 3139 .keywords: matrix, array, elements, values 3140 3141 .seealso: MatRestoreArray(), MatGetArrayF90() 3142 @*/ 3143 int MatGetArray(Mat mat,Scalar **v) 3144 { 3145 int ierr; 3146 3147 PetscFunctionBegin; 3148 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3149 PetscValidPointer(v); 3150 if (!mat->ops->getarray) SETERRQ(PETSC_ERR_SUP,0,""); 3151 ierr = (*mat->ops->getarray)(mat,v);CHKERRQ(ierr); 3152 PetscFunctionReturn(0); 3153 } 3154 3155 #undef __FUNC__ 3156 #define __FUNC__ "MatRestoreArray" 3157 /*@C 3158 MatRestoreArray - Restores the matrix after MatGetArray() has been called. 3159 3160 Not Collective 3161 3162 Input Parameter: 3163 + mat - the matrix 3164 - v - the location of the values 3165 3166 Fortran Note: 3167 This routine is used differently from Fortran, e.g., 3168 .vb 3169 Mat mat 3170 Scalar mat_array(1) 3171 PetscOffset i_mat 3172 int ierr 3173 call MatGetArray(mat,mat_array,i_mat,ierr) 3174 3175 C Access first local entry in matrix; note that array is 3176 C treated as one dimensional 3177 value = mat_array(i_mat + 1) 3178 3179 [... other code ...] 3180 call MatRestoreArray(mat,mat_array,i_mat,ierr) 3181 .ve 3182 3183 See the Fortran chapter of the users manual and 3184 petsc/src/mat/examples/tests for details 3185 3186 Level: advanced 3187 3188 .keywords: matrix, array, elements, values, restore 3189 3190 .seealso: MatGetArray(), MatRestoreArrayF90() 3191 @*/ 3192 int MatRestoreArray(Mat mat,Scalar **v) 3193 { 3194 int ierr; 3195 3196 PetscFunctionBegin; 3197 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3198 PetscValidPointer(v); 3199 if (!mat->ops->restorearray) SETERRQ(PETSC_ERR_SUP,0,""); 3200 ierr = (*mat->ops->restorearray)(mat,v);CHKERRQ(ierr); 3201 PetscFunctionReturn(0); 3202 } 3203 3204 #undef __FUNC__ 3205 #define __FUNC__ "MatGetSubMatrices" 3206 /*@C 3207 MatGetSubMatrices - Extracts several submatrices from a matrix. If submat 3208 points to an array of valid matrices, they may be reused to store the new 3209 submatrices. 3210 3211 Collective on Mat 3212 3213 Input Parameters: 3214 + mat - the matrix 3215 . n - the number of submatrixes to be extracted (on this processor, may be zero) 3216 . irow, icol - index sets of rows and columns to extract 3217 - scall - either MAT_INITIAL_MATRIX or MAT_REUSE_MATRIX 3218 3219 Output Parameter: 3220 . submat - the array of submatrices 3221 3222 Notes: 3223 MatGetSubMatrices() can extract only sequential submatrices 3224 (from both sequential and parallel matrices). Use MatGetSubMatrix() 3225 to extract a parallel submatrix. 3226 3227 When extracting submatrices from a parallel matrix, each processor can 3228 form a different submatrix by setting the rows and columns of its 3229 individual index sets according to the local submatrix desired. 3230 3231 When finished using the submatrices, the user should destroy 3232 them with MatDestroySubMatrices(). 3233 3234 MAT_REUSE_MATRIX can only be used when the nonzero structure of the 3235 original matrix has not changed from that last call to MatGetSubMatrices(). 3236 3237 Fortran Note: 3238 The Fortran interface is slightly different from that given below; it 3239 requires one to pass in as submat a Mat (integer) array of size at least m. 3240 3241 Level: advanced 3242 3243 .keywords: matrix, get, submatrix, submatrices 3244 3245 .seealso: MatDestroyMatrices(), MatGetSubMatrix(), MatGetRow(), MatGetDiagonal() 3246 @*/ 3247 int MatGetSubMatrices(Mat mat,int n,IS *irow,IS *icol,MatReuse scall,Mat **submat) 3248 { 3249 int ierr; 3250 3251 PetscFunctionBegin; 3252 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3253 if (!mat->ops->getsubmatrices) SETERRQ(PETSC_ERR_SUP,0,""); 3254 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 3255 3256 PLogEventBegin(MAT_GetSubMatrices,mat,0,0,0); 3257 ierr = (*mat->ops->getsubmatrices)(mat,n,irow,icol,scall,submat);CHKERRQ(ierr); 3258 PLogEventEnd(MAT_GetSubMatrices,mat,0,0,0); 3259 3260 PetscFunctionReturn(0); 3261 } 3262 3263 #undef __FUNC__ 3264 #define __FUNC__ "MatDestroyMatrices" 3265 /*@C 3266 MatDestroyMatrices - Destroys a set of matrices obtained with MatGetSubMatrices(). 3267 3268 Collective on Mat 3269 3270 Input Parameters: 3271 + n - the number of local matrices 3272 - mat - the matrices 3273 3274 Level: advanced 3275 3276 .keywords: matrix, destroy, submatrix, submatrices 3277 3278 .seealso: MatGetSubMatrices() 3279 @*/ 3280 int MatDestroyMatrices(int n,Mat **mat) 3281 { 3282 int ierr,i; 3283 3284 PetscFunctionBegin; 3285 if (n < 0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,1,"Trying to destroy negative number of matrices %d",n); 3286 PetscValidPointer(mat); 3287 for ( i=0; i<n; i++ ) { 3288 ierr = MatDestroy((*mat)[i]);CHKERRQ(ierr); 3289 } 3290 if (n) {ierr = PetscFree(*mat);CHKERRQ(ierr);} 3291 PetscFunctionReturn(0); 3292 } 3293 3294 #undef __FUNC__ 3295 #define __FUNC__ "MatIncreaseOverlap" 3296 /*@ 3297 MatIncreaseOverlap - Given a set of submatrices indicated by index sets, 3298 replaces the index sets by larger ones that represent submatrices with 3299 additional overlap. 3300 3301 Collective on Mat 3302 3303 Input Parameters: 3304 + mat - the matrix 3305 . n - the number of index sets 3306 . is - the array of pointers to index sets 3307 - ov - the additional overlap requested 3308 3309 Level: developer 3310 3311 .keywords: matrix, overlap, Schwarz 3312 3313 .seealso: MatGetSubMatrices() 3314 @*/ 3315 int MatIncreaseOverlap(Mat mat,int n, IS *is,int ov) 3316 { 3317 int ierr; 3318 3319 PetscFunctionBegin; 3320 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3321 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 3322 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 3323 3324 if (ov == 0) PetscFunctionReturn(0); 3325 if (!mat->ops->increaseoverlap) SETERRQ(PETSC_ERR_SUP,0,""); 3326 PLogEventBegin(MAT_IncreaseOverlap,mat,0,0,0); 3327 ierr = (*mat->ops->increaseoverlap)(mat,n,is,ov);CHKERRQ(ierr); 3328 PLogEventEnd(MAT_IncreaseOverlap,mat,0,0,0); 3329 PetscFunctionReturn(0); 3330 } 3331 3332 #undef __FUNC__ 3333 #define __FUNC__ "MatPrintHelp" 3334 /*@ 3335 MatPrintHelp - Prints all the options for the matrix. 3336 3337 Collective on Mat 3338 3339 Input Parameter: 3340 . mat - the matrix 3341 3342 Options Database Keys: 3343 + -help - Prints matrix options 3344 - -h - Prints matrix options 3345 3346 Level: developer 3347 3348 .keywords: mat, help 3349 3350 .seealso: MatCreate(), MatCreateXXX() 3351 @*/ 3352 int MatPrintHelp(Mat mat) 3353 { 3354 static PetscTruth called = PETSC_FALSE; 3355 int ierr; 3356 MPI_Comm comm; 3357 3358 PetscFunctionBegin; 3359 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3360 3361 comm = mat->comm; 3362 if (!called) { 3363 ierr = (*PetscHelpPrintf)(comm,"General matrix options:\n");CHKERRQ(ierr); 3364 ierr = (*PetscHelpPrintf)(comm," -mat_view_info: view basic matrix info during MatAssemblyEnd()\n");CHKERRQ(ierr); 3365 ierr = (*PetscHelpPrintf)(comm," -mat_view_info_detailed: view detailed matrix info during MatAssemblyEnd()\n");CHKERRQ(ierr); 3366 ierr = (*PetscHelpPrintf)(comm," -mat_view_draw: draw nonzero matrix structure during MatAssemblyEnd()\n");CHKERRQ(ierr); 3367 ierr = (*PetscHelpPrintf)(comm," -draw_pause <sec>: set seconds of display pause\n");CHKERRQ(ierr); 3368 ierr = (*PetscHelpPrintf)(comm," -display <name>: set alternate display\n");CHKERRQ(ierr); 3369 called = PETSC_TRUE; 3370 } 3371 if (mat->ops->printhelp) { 3372 ierr = (*mat->ops->printhelp)(mat);CHKERRQ(ierr); 3373 } 3374 PetscFunctionReturn(0); 3375 } 3376 3377 #undef __FUNC__ 3378 #define __FUNC__ "MatGetBlockSize" 3379 /*@ 3380 MatGetBlockSize - Returns the matrix block size; useful especially for the 3381 block row and block diagonal formats. 3382 3383 Not Collective 3384 3385 Input Parameter: 3386 . mat - the matrix 3387 3388 Output Parameter: 3389 . bs - block size 3390 3391 Notes: 3392 Block diagonal formats are MATSEQBDIAG, MATMPIBDIAG. 3393 Block row formats are MATSEQBAIJ, MATMPIBAIJ 3394 3395 Level: intermediate 3396 3397 .keywords: matrix, get, block, size 3398 3399 .seealso: MatCreateSeqBAIJ(), MatCreateMPIBAIJ(), MatCreateSeqBDiag(), MatCreateMPIBDiag() 3400 @*/ 3401 int MatGetBlockSize(Mat mat,int *bs) 3402 { 3403 int ierr; 3404 3405 PetscFunctionBegin; 3406 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3407 PetscValidIntPointer(bs); 3408 if (!mat->ops->getblocksize) SETERRQ(PETSC_ERR_SUP,0,""); 3409 ierr = (*mat->ops->getblocksize)(mat,bs);CHKERRQ(ierr); 3410 PetscFunctionReturn(0); 3411 } 3412 3413 #undef __FUNC__ 3414 #define __FUNC__ "MatGetRowIJ" 3415 /*@C 3416 MatGetRowIJ - Returns the compressed row storage i and j indices for sequential matrices. 3417 3418 Collective on Mat 3419 3420 Input Parameters: 3421 + mat - the matrix 3422 . shift - 0 or 1 indicating we want the indices starting at 0 or 1 3423 - symmetric - PETSC_TRUE or PETSC_FALSE indicating the matrix data structure should be 3424 symmetrized 3425 3426 Output Parameters: 3427 + n - number of rows in the (possibly compressed) matrix 3428 . ia - the row pointers 3429 . ja - the column indices 3430 - done - PETSC_TRUE or PETSC_FALSE, indicating whether the values have been returned 3431 3432 Level: developer 3433 3434 .seealso: MatGetColumnIJ(), MatRestoreRowIJ() 3435 @*/ 3436 int MatGetRowIJ(Mat mat,int shift,PetscTruth symmetric,int *n,int **ia,int** ja,PetscTruth *done) 3437 { 3438 int ierr; 3439 3440 PetscFunctionBegin; 3441 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3442 if (ia) PetscValidIntPointer(ia); 3443 if (ja) PetscValidIntPointer(ja); 3444 PetscValidIntPointer(done); 3445 if (!mat->ops->getrowij) *done = PETSC_FALSE; 3446 else { 3447 *done = PETSC_TRUE; 3448 ierr = (*mat->ops->getrowij)(mat,shift,symmetric,n,ia,ja,done);CHKERRQ(ierr); 3449 } 3450 PetscFunctionReturn(0); 3451 } 3452 3453 #undef __FUNC__ 3454 #define __FUNC__ "MatGetColumnIJ" 3455 /*@C 3456 MatGetColumnIJ - Returns the compressed column storage i and j indices for sequential matrices. 3457 3458 Collective on Mat 3459 3460 Input Parameters: 3461 + mat - the matrix 3462 . shift - 1 or zero indicating we want the indices starting at 0 or 1 3463 - symmetric - PETSC_TRUE or PETSC_FALSE indicating the matrix data structure should be 3464 symmetrized 3465 3466 Output Parameters: 3467 + n - number of columns in the (possibly compressed) matrix 3468 . ia - the column pointers 3469 . ja - the row indices 3470 - done - PETSC_TRUE or PETSC_FALSE, indicating whether the values have been returned 3471 3472 Level: developer 3473 3474 .seealso: MatGetRowIJ(), MatRestoreColumnIJ() 3475 @*/ 3476 int MatGetColumnIJ(Mat mat,int shift,PetscTruth symmetric,int *n,int **ia,int** ja,PetscTruth *done) 3477 { 3478 int ierr; 3479 3480 PetscFunctionBegin; 3481 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3482 if (ia) PetscValidIntPointer(ia); 3483 if (ja) PetscValidIntPointer(ja); 3484 PetscValidIntPointer(done); 3485 3486 if (!mat->ops->getcolumnij) *done = PETSC_FALSE; 3487 else { 3488 *done = PETSC_TRUE; 3489 ierr = (*mat->ops->getcolumnij)(mat,shift,symmetric,n,ia,ja,done);CHKERRQ(ierr); 3490 } 3491 PetscFunctionReturn(0); 3492 } 3493 3494 #undef __FUNC__ 3495 #define __FUNC__ "MatRestoreRowIJ" 3496 /*@C 3497 MatRestoreRowIJ - Call after you are completed with the ia,ja indices obtained with 3498 MatGetRowIJ(). 3499 3500 Collective on Mat 3501 3502 Input Parameters: 3503 + mat - the matrix 3504 . shift - 1 or zero indicating we want the indices starting at 0 or 1 3505 - symmetric - PETSC_TRUE or PETSC_FALSE indicating the matrix data structure should be 3506 symmetrized 3507 3508 Output Parameters: 3509 + n - size of (possibly compressed) matrix 3510 . ia - the row pointers 3511 . ja - the column indices 3512 - done - PETSC_TRUE or PETSC_FALSE indicated that the values have been returned 3513 3514 Level: developer 3515 3516 .seealso: MatGetRowIJ(), MatRestoreColumnIJ() 3517 @*/ 3518 int MatRestoreRowIJ(Mat mat,int shift,PetscTruth symmetric,int *n,int **ia,int** ja,PetscTruth *done) 3519 { 3520 int ierr; 3521 3522 PetscFunctionBegin; 3523 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3524 if (ia) PetscValidIntPointer(ia); 3525 if (ja) PetscValidIntPointer(ja); 3526 PetscValidIntPointer(done); 3527 3528 if (!mat->ops->restorerowij) *done = PETSC_FALSE; 3529 else { 3530 *done = PETSC_TRUE; 3531 ierr = (*mat->ops->restorerowij)(mat,shift,symmetric,n,ia,ja,done);CHKERRQ(ierr); 3532 } 3533 PetscFunctionReturn(0); 3534 } 3535 3536 #undef __FUNC__ 3537 #define __FUNC__ "MatRestoreColumnIJ" 3538 /*@C 3539 MatRestoreColumnIJ - Call after you are completed with the ia,ja indices obtained with 3540 MatGetColumnIJ(). 3541 3542 Collective on Mat 3543 3544 Input Parameters: 3545 + mat - the matrix 3546 . shift - 1 or zero indicating we want the indices starting at 0 or 1 3547 - symmetric - PETSC_TRUE or PETSC_FALSE indicating the matrix data structure should be 3548 symmetrized 3549 3550 Output Parameters: 3551 + n - size of (possibly compressed) matrix 3552 . ia - the column pointers 3553 . ja - the row indices 3554 - done - PETSC_TRUE or PETSC_FALSE indicated that the values have been returned 3555 3556 Level: developer 3557 3558 .seealso: MatGetColumnIJ(), MatRestoreRowIJ() 3559 @*/ 3560 int MatRestoreColumnIJ(Mat mat,int shift,PetscTruth symmetric,int *n,int **ia,int** ja,PetscTruth *done) 3561 { 3562 int ierr; 3563 3564 PetscFunctionBegin; 3565 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3566 if (ia) PetscValidIntPointer(ia); 3567 if (ja) PetscValidIntPointer(ja); 3568 PetscValidIntPointer(done); 3569 3570 if (!mat->ops->restorecolumnij) *done = PETSC_FALSE; 3571 else { 3572 *done = PETSC_TRUE; 3573 ierr = (*mat->ops->restorecolumnij)(mat,shift,symmetric,n,ia,ja,done);CHKERRQ(ierr); 3574 } 3575 PetscFunctionReturn(0); 3576 } 3577 3578 #undef __FUNC__ 3579 #define __FUNC__ "MatColoringPatch" 3580 /*@C 3581 MatColoringPatch -Used inside matrix coloring routines that 3582 use MatGetRowIJ() and/or MatGetColumnIJ(). 3583 3584 Collective on Mat 3585 3586 Input Parameters: 3587 + mat - the matrix 3588 . n - number of colors 3589 - colorarray - array indicating color for each column 3590 3591 Output Parameters: 3592 . iscoloring - coloring generated using colorarray information 3593 3594 Level: developer 3595 3596 .seealso: MatGetRowIJ(), MatGetColumnIJ() 3597 3598 .keywords: mat, coloring, patch 3599 @*/ 3600 int MatColoringPatch(Mat mat,int n,int *colorarray,ISColoring *iscoloring) 3601 { 3602 int ierr; 3603 3604 PetscFunctionBegin; 3605 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3606 PetscValidIntPointer(colorarray); 3607 3608 if (!mat->ops->coloringpatch) {SETERRQ(PETSC_ERR_SUP,0,"");} 3609 else { 3610 ierr = (*mat->ops->coloringpatch)(mat,n,colorarray,iscoloring);CHKERRQ(ierr); 3611 } 3612 PetscFunctionReturn(0); 3613 } 3614 3615 3616 #undef __FUNC__ 3617 #define __FUNC__ "MatSetUnfactored" 3618 /*@ 3619 MatSetUnfactored - Resets a factored matrix to be treated as unfactored. 3620 3621 Collective on Mat 3622 3623 Input Parameter: 3624 . mat - the factored matrix to be reset 3625 3626 Notes: 3627 This routine should be used only with factored matrices formed by in-place 3628 factorization via ILU(0) (or by in-place LU factorization for the MATSEQDENSE 3629 format). This option can save memory, for example, when solving nonlinear 3630 systems with a matrix-free Newton-Krylov method and a matrix-based, in-place 3631 ILU(0) preconditioner. 3632 3633 Note that one can specify in-place ILU(0) factorization by calling 3634 .vb 3635 PCType(pc,PCILU); 3636 PCILUSeUseInPlace(pc); 3637 .ve 3638 or by using the options -pc_type ilu -pc_ilu_in_place 3639 3640 In-place factorization ILU(0) can also be used as a local 3641 solver for the blocks within the block Jacobi or additive Schwarz 3642 methods (runtime option: -sub_pc_ilu_in_place). See the discussion 3643 of these preconditioners in the users manual for details on setting 3644 local solver options. 3645 3646 Most users should employ the simplified SLES interface for linear solvers 3647 instead of working directly with matrix algebra routines such as this. 3648 See, e.g., SLESCreate(). 3649 3650 Level: developer 3651 3652 .seealso: PCILUSetUseInPlace(), PCLUSetUseInPlace() 3653 3654 .keywords: matrix-free, in-place ILU, in-place LU 3655 @*/ 3656 int MatSetUnfactored(Mat mat) 3657 { 3658 int ierr; 3659 3660 PetscFunctionBegin; 3661 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3662 mat->factor = 0; 3663 if (!mat->ops->setunfactored) PetscFunctionReturn(0); 3664 ierr = (*mat->ops->setunfactored)(mat);CHKERRQ(ierr); 3665 PetscFunctionReturn(0); 3666 } 3667 3668 #undef __FUNC__ 3669 #define __FUNC__ "MatGetType" 3670 /*@C 3671 MatGetType - Gets the matrix type and name (as a string) from the matrix. 3672 3673 Not Collective 3674 3675 Input Parameter: 3676 . mat - the matrix 3677 3678 Output Parameter: 3679 + type - the matrix type (or use PETSC_NULL) 3680 - name - name of matrix type (or use PETSC_NULL) 3681 3682 Level: intermediate 3683 3684 .keywords: matrix, get, type, name 3685 @*/ 3686 int MatGetType(Mat mat,MatType *type,char **name) 3687 { 3688 int itype = (int)mat->type; 3689 char *matname[10]; 3690 3691 PetscFunctionBegin; 3692 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3693 3694 if (type) *type = (MatType) mat->type; 3695 if (name) { 3696 /* Note: Be sure that this list corresponds to the enum in mat.h */ 3697 matname[0] = "MATSEQDENSE"; 3698 matname[1] = "MATSEQAIJ"; 3699 matname[2] = "MATMPIAIJ"; 3700 matname[3] = "MATSHELL"; 3701 matname[4] = "MATMPIROWBS"; 3702 matname[5] = "MATSEQBDIAG"; 3703 matname[6] = "MATMPIBDIAG"; 3704 matname[7] = "MATMPIDENSE"; 3705 matname[8] = "MATSEQBAIJ"; 3706 matname[9] = "MATMPIBAIJ"; 3707 3708 if (itype < 0 || itype > 9) *name = "Unknown matrix type"; 3709 else *name = matname[itype]; 3710 } 3711 PetscFunctionReturn(0); 3712 } 3713 3714 /*MC 3715 MatGetArrayF90 - Accesses a matrix array from Fortran90. 3716 3717 Synopsis: 3718 MatGetArrayF90(Mat x,{Scalar, pointer :: xx_v(:)},integer ierr) 3719 3720 Not collective 3721 3722 Input Parameter: 3723 . x - matrix 3724 3725 Output Parameters: 3726 + xx_v - the Fortran90 pointer to the array 3727 - ierr - error code 3728 3729 Example of Usage: 3730 .vb 3731 Scalar, pointer xx_v(:) 3732 .... 3733 call MatGetArrayF90(x,xx_v,ierr) 3734 a = xx_v(3) 3735 call MatRestoreArrayF90(x,xx_v,ierr) 3736 .ve 3737 3738 Notes: 3739 Not yet supported for all F90 compilers 3740 3741 Level: advanced 3742 3743 .seealso: MatRestoreArrayF90(), MatGetArray(), MatRestoreArray() 3744 3745 .keywords: matrix, array, f90 3746 M*/ 3747 3748 /*MC 3749 MatRestoreArrayF90 - Restores a matrix array that has been 3750 accessed with MatGetArrayF90(). 3751 3752 Synopsis: 3753 MatRestoreArrayF90(Mat x,{Scalar, pointer :: xx_v(:)},integer ierr) 3754 3755 Not collective 3756 3757 Input Parameters: 3758 + x - matrix 3759 - xx_v - the Fortran90 pointer to the array 3760 3761 Output Parameter: 3762 . ierr - error code 3763 3764 Example of Usage: 3765 .vb 3766 Scalar, pointer xx_v(:) 3767 .... 3768 call MatGetArrayF90(x,xx_v,ierr) 3769 a = xx_v(3) 3770 call MatRestoreArrayF90(x,xx_v,ierr) 3771 .ve 3772 3773 Notes: 3774 Not yet supported for all F90 compilers 3775 3776 Level: advanced 3777 3778 .seealso: MatGetArrayF90(), MatGetArray(), MatRestoreArray() 3779 3780 .keywords: matrix, array, f90 3781 M*/ 3782 3783 3784 #undef __FUNC__ 3785 #define __FUNC__ "MatGetSubMatrix" 3786 /*@ 3787 MatGetSubMatrix - Gets a single submatrix on the same number of processors 3788 as the original matrix. 3789 3790 Collective on Mat 3791 3792 Input Parameters: 3793 + mat - the original matrix 3794 . isrow - rows this processor should obtain 3795 . iscol - columns for all processors you wish to keep 3796 . csize - number of columns "local" to this processor (does nothing for sequential 3797 matrices). This should match the result from VecGetLocalSize(x,...) if you 3798 plan to use the matrix in a A*x; alternatively, you can use PETSC_DECIDE 3799 - cll - either MAT_INITIAL_MATRIX or MAT_REUSE_MATRIX 3800 3801 Output Parameter: 3802 . newmat - the new submatrix, of the same type as the old 3803 3804 Level: advanced 3805 3806 .keywords: matrix, get, submatrix, submatrices 3807 3808 .seealso: MatGetSubMatrices() 3809 @*/ 3810 int MatGetSubMatrix(Mat mat,IS isrow,IS iscol,int csize,MatReuse cll,Mat *newmat) 3811 { 3812 int ierr, size; 3813 Mat *local; 3814 3815 PetscFunctionBegin; 3816 ierr = MPI_Comm_size(mat->comm,&size);CHKERRQ(ierr); 3817 3818 /* if original matrix is on just one processor then use submatrix generated */ 3819 if (!mat->ops->getsubmatrix && size == 1 && cll == MAT_REUSE_MATRIX) { 3820 ierr = MatGetSubMatrices(mat,1,&isrow,&iscol,MAT_REUSE_MATRIX,&newmat);CHKERRQ(ierr); 3821 PetscFunctionReturn(0); 3822 } else if (!mat->ops->getsubmatrix && size == 1) { 3823 ierr = MatGetSubMatrices(mat,1,&isrow,&iscol,MAT_INITIAL_MATRIX,&local);CHKERRQ(ierr); 3824 *newmat = *local; 3825 ierr = PetscFree(local);CHKERRQ(ierr); 3826 PetscFunctionReturn(0); 3827 } 3828 3829 if (!mat->ops->getsubmatrix) SETERRQ(PETSC_ERR_SUP,0,"Not currently implemented"); 3830 ierr = (*mat->ops->getsubmatrix)(mat,isrow,iscol,csize,cll,newmat);CHKERRQ(ierr); 3831 PetscFunctionReturn(0); 3832 } 3833 3834 #undef __FUNC__ 3835 #define __FUNC__ "MatGetMaps" 3836 /*@C 3837 MatGetMaps - Returns the maps associated with the matrix. 3838 3839 Not Collective 3840 3841 Input Parameter: 3842 . mat - the matrix 3843 3844 Output Parameters: 3845 + rmap - the row (right) map 3846 - cmap - the column (left) map 3847 3848 Level: developer 3849 3850 .keywords: matrix, get, map 3851 @*/ 3852 int MatGetMaps(Mat mat,Map *rmap,Map *cmap) 3853 { 3854 int ierr; 3855 3856 PetscFunctionBegin; 3857 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3858 ierr = (*mat->ops->getmaps)(mat,rmap,cmap);CHKERRQ(ierr); 3859 PetscFunctionReturn(0); 3860 } 3861 3862 /* 3863 Version that works for all PETSc matrices 3864 */ 3865 #undef __FUNC__ 3866 #define __FUNC__ "MatGetMaps_Petsc" 3867 int MatGetMaps_Petsc(Mat mat,Map *rmap,Map *cmap) 3868 { 3869 PetscFunctionBegin; 3870 if (rmap) *rmap = mat->rmap; 3871 if (cmap) *cmap = mat->cmap; 3872 PetscFunctionReturn(0); 3873 } 3874 3875 #undef __FUNC__ 3876 #define __FUNC__ "MatSetStashInitialSize" 3877 /*@ 3878 MatSetStashInitialSize - sets the sizes of the matrix stash, that is 3879 used during the assembly process to store values that belong to 3880 other processors. 3881 3882 Not Collective 3883 3884 Input Parameters: 3885 + mat - the matrix 3886 . size - the initial size of the stash. 3887 - bsize - the initial size of the block-stash(if used). 3888 3889 Options Database Keys: 3890 + -matstash_initial_size <size> or <size0,size1,...sizep-1> 3891 - -matstash_block_initial_size <bsize> or <bsize0,bsize1,...bsizep-1> 3892 3893 Level: intermediate 3894 3895 Notes: 3896 The block-stash is used for values set with VecSetValuesBlocked() while 3897 the stash is used for values set with VecSetValues() 3898 3899 Run with the option -log_info and look for output of the form 3900 MatAssemblyBegin_MPIXXX:Stash has MM entries, uses nn mallocs. 3901 to determine the appropriate value, MM, to use for size and 3902 MatAssemblyBegin_MPIXXX:Block-Stash has BMM entries, uses nn mallocs. 3903 to determine the value, BMM to use for bsize 3904 3905 .keywords: matrix, stash, assembly 3906 @*/ 3907 int MatSetStashInitialSize(Mat mat,int size, int bsize) 3908 { 3909 int ierr; 3910 3911 PetscFunctionBegin; 3912 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3913 ierr = MatStashSetInitialSize_Private(&mat->stash,size);CHKERRQ(ierr); 3914 ierr = MatStashSetInitialSize_Private(&mat->bstash,bsize);CHKERRQ(ierr); 3915 PetscFunctionReturn(0); 3916 } 3917 3918 #undef __FUNC__ 3919 #define __FUNC__ "MatInterpolateAdd" 3920 /*@ 3921 MatInterpolateAdd - w = y + A*x or A'*x depending on the shape of 3922 the matrix 3923 3924 Collective on Mat 3925 3926 Input Parameters: 3927 + mat - the matrix 3928 . x,y - the vectors 3929 - w - where the result is stored 3930 3931 Level: intermediate 3932 3933 Notes: 3934 w may be the same vector as y. 3935 3936 This allows one to use either the restriction or interpolation (its transpose) 3937 matrix to do the interpolation 3938 3939 .keywords: interpolate, 3940 3941 .seealso: MatMultAdd(), MatMultTransposeAdd(), MatRestrict() 3942 3943 @*/ 3944 int MatInterpolateAdd(Mat A,Vec x,Vec y,Vec w) 3945 { 3946 int M,N,ierr; 3947 3948 PetscFunctionBegin; 3949 ierr = MatGetSize(A,&M,&N);CHKERRQ(ierr); 3950 if (N > M) { 3951 ierr = MatMultTransposeAdd(A,x,y,w);CHKERRQ(ierr); 3952 } else { 3953 ierr = MatMultAdd(A,x,y,w);CHKERRQ(ierr); 3954 } 3955 PetscFunctionReturn(0); 3956 } 3957 3958 #undef __FUNC__ 3959 #define __FUNC__ "MatInterpolate" 3960 /*@ 3961 MatInterpolate - y = A*x or A'*x depending on the shape of 3962 the matrix 3963 3964 Collective on Mat 3965 3966 Input Parameters: 3967 + mat - the matrix 3968 - x,y - the vectors 3969 3970 Level: intermediate 3971 3972 Notes: 3973 This allows one to use either the restriction or interpolation (its transpose) 3974 matrix to do the interpolation 3975 3976 .keywords: interpolate, 3977 3978 .seealso: MatMultAdd(), MatMultTransposeAdd(), MatRestrict() 3979 3980 @*/ 3981 int MatInterpolate(Mat A,Vec x,Vec y) 3982 { 3983 int M,N,ierr; 3984 3985 PetscFunctionBegin; 3986 ierr = MatGetSize(A,&M,&N);CHKERRQ(ierr); 3987 if (N > M) { 3988 ierr = MatMultTranspose(A,x,y);CHKERRQ(ierr); 3989 } else { 3990 ierr = MatMult(A,x,y);CHKERRQ(ierr); 3991 } 3992 PetscFunctionReturn(0); 3993 } 3994 3995 #undef __FUNC__ 3996 #define __FUNC__ "MatRestrict" 3997 /*@ 3998 MatRestrict - y = A*x or A'*x 3999 4000 Collective on Mat 4001 4002 Input Parameters: 4003 + mat - the matrix 4004 - x,y - the vectors 4005 4006 Level: intermediate 4007 4008 Notes: 4009 This allows one to use either the restriction or interpolation (its transpose) 4010 matrix to do the restriction 4011 4012 .keywords: interpolate, 4013 4014 .seealso: MatMultAdd(), MatMultTransposeAdd(), MatInterpolate() 4015 4016 @*/ 4017 int MatRestrict(Mat A,Vec x,Vec y) 4018 { 4019 int M,N,ierr; 4020 4021 PetscFunctionBegin; 4022 ierr = MatGetSize(A,&M,&N);CHKERRQ(ierr); 4023 if (N > M) { 4024 ierr = MatMult(A,x,y);CHKERRQ(ierr); 4025 } else { 4026 ierr = MatMultTranspose(A,x,y);CHKERRQ(ierr); 4027 } 4028 PetscFunctionReturn(0); 4029 } 4030