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