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