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