1 #ifdef PETSC_RCS_HEADER 2 static char vcid[] = "$Id: matrix.c,v 1.311 1999/01/12 20:44:41 curfman 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 - ViewerSocketOpen() - Outputs matrix to Socket viewer. 159 Currently only the sequential dense and AIJ 160 matrix types support the Socket 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 ViewerSocketOpen(), 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 SETERRQ2(PETSC_ERR_SUP,0,"Number column/row indices must be <= 2048: are %d %d",nrow,ncol); 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 SETERRQ2(PETSC_ERR_SUP,0,"Number column/row indices must be <= 2048: are %d %d",nrow,ncol); 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) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim %d %d",mat->N,x->N); 725 if (mat->M != y->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec y: global dim %d %d",mat->M,y->N); 726 if (mat->m != y->n) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec y: local dim %d %d",mat->m,y->n); 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) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim %d %d",mat->M,x->N); 768 if (mat->N != y->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec y: global dim %d %d",mat->N,y->N); 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) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec v1: global dim %d %d",mat->N,v1->N); 807 if (mat->M != v2->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec v2: global dim %d %d",mat->M,v2->N); 808 if (mat->M != v3->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec v3: global dim %d %d",mat->M,v3->N); 809 if (mat->m != v3->n) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec v3: local dim %d %d",mat->m,v3->n); 810 if (mat->m != v2->n) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec v2: local dim %d %d",mat->m,v2->n); 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) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec v1: global dim %d %d",mat->M,v1->N); 853 if (mat->N != v2->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec v2: global dim %d %d",mat->N,v2->N); 854 if (mat->N != v3->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec v3: global dim %d %d",mat->N,v3->N); 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 SETERRQ4(PETSC_ERR_ARG_SIZ,0,"Mat mat,Mat *fact: global dimensions are different %d should = %d %d should = %d", 1125 mat->M,(*fact)->M,mat->N,(*fact)->N); 1126 } 1127 if (!(*fact)->ops->lufactornumeric) SETERRQ(PETSC_ERR_SUP,0,""); 1128 1129 PLogEventBegin(MAT_LUFactorNumeric,mat,*fact,0,0); 1130 ierr = (*(*fact)->ops->lufactornumeric)(mat,fact); CHKERRQ(ierr); 1131 PLogEventEnd(MAT_LUFactorNumeric,mat,*fact,0,0); 1132 ierr = OptionsHasName(PETSC_NULL,"-mat_view_draw",&flg); CHKERRQ(ierr); 1133 if (flg) { 1134 ierr = OptionsHasName(0,"-mat_view_contour",&flg); CHKERRQ(ierr); 1135 if (flg) { 1136 ViewerPushFormat(VIEWER_DRAWX_(mat->comm),VIEWER_FORMAT_DRAW_CONTOUR,0);CHKERRQ(ierr); 1137 } 1138 ierr = MatView(*fact,VIEWER_DRAWX_(mat->comm)); CHKERRQ(ierr); 1139 ierr = ViewerFlush(VIEWER_DRAWX_(mat->comm)); CHKERRQ(ierr); 1140 if (flg) { 1141 ViewerPopFormat(VIEWER_DRAWX_(mat->comm));CHKERRQ(ierr); 1142 } 1143 } 1144 PetscFunctionReturn(0); 1145 } 1146 1147 #undef __FUNC__ 1148 #define __FUNC__ "MatCholeskyFactor" 1149 /*@ 1150 MatCholeskyFactor - Performs in-place Cholesky factorization of a 1151 symmetric matrix. 1152 1153 Collective on Mat 1154 1155 Input Parameters: 1156 + mat - the matrix 1157 . perm - row and column permutations 1158 - f - expected fill as ratio of original fill 1159 1160 Notes: 1161 See MatLUFactor() for the nonsymmetric case. See also 1162 MatCholeskyFactorSymbolic(), and MatCholeskyFactorNumeric(). 1163 1164 .keywords: matrix, factor, in-place, Cholesky 1165 1166 .seealso: MatLUFactor(), MatCholeskyFactorSymbolic(), MatCholeskyFactorNumeric() 1167 @*/ 1168 int MatCholeskyFactor(Mat mat,IS perm,double f) 1169 { 1170 int ierr; 1171 1172 PetscFunctionBegin; 1173 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1174 if (mat->M != mat->N) SETERRQ(PETSC_ERR_ARG_WRONG,0,"Matrix must be square"); 1175 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 1176 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 1177 if (!mat->ops->choleskyfactor) SETERRQ(PETSC_ERR_SUP,0,""); 1178 1179 PLogEventBegin(MAT_CholeskyFactor,mat,perm,0,0); 1180 ierr = (*mat->ops->choleskyfactor)(mat,perm,f); CHKERRQ(ierr); 1181 PLogEventEnd(MAT_CholeskyFactor,mat,perm,0,0); 1182 PetscFunctionReturn(0); 1183 } 1184 1185 #undef __FUNC__ 1186 #define __FUNC__ "MatCholeskyFactorSymbolic" 1187 /*@ 1188 MatCholeskyFactorSymbolic - Performs symbolic Cholesky factorization 1189 of a symmetric matrix. 1190 1191 Collective on Mat 1192 1193 Input Parameters: 1194 + mat - the matrix 1195 . perm - row and column permutations 1196 - f - expected fill as ratio of original 1197 1198 Output Parameter: 1199 . fact - the factored matrix 1200 1201 Notes: 1202 See MatLUFactorSymbolic() for the nonsymmetric case. See also 1203 MatCholeskyFactor() and MatCholeskyFactorNumeric(). 1204 1205 .keywords: matrix, factor, factorization, symbolic, Cholesky 1206 1207 .seealso: MatLUFactorSymbolic(), MatCholeskyFactor(), MatCholeskyFactorNumeric() 1208 @*/ 1209 int MatCholeskyFactorSymbolic(Mat mat,IS perm,double f,Mat *fact) 1210 { 1211 int ierr; 1212 1213 PetscFunctionBegin; 1214 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1215 PetscValidPointer(fact); 1216 if (mat->M != mat->N) SETERRQ(PETSC_ERR_ARG_WRONG,0,"Matrix must be square"); 1217 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 1218 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 1219 if (!mat->ops->choleskyfactorsymbolic) SETERRQ(PETSC_ERR_SUP,0,""); 1220 1221 PLogEventBegin(MAT_CholeskyFactorSymbolic,mat,perm,0,0); 1222 ierr = (*mat->ops->choleskyfactorsymbolic)(mat,perm,f,fact); CHKERRQ(ierr); 1223 PLogEventEnd(MAT_CholeskyFactorSymbolic,mat,perm,0,0); 1224 PetscFunctionReturn(0); 1225 } 1226 1227 #undef __FUNC__ 1228 #define __FUNC__ "MatCholeskyFactorNumeric" 1229 /*@ 1230 MatCholeskyFactorNumeric - Performs numeric Cholesky factorization 1231 of a symmetric matrix. Call this routine after first calling 1232 MatCholeskyFactorSymbolic(). 1233 1234 Collective on Mat 1235 1236 Input Parameter: 1237 . mat - the initial matrix 1238 1239 Output Parameter: 1240 . fact - the factored matrix 1241 1242 .keywords: matrix, factor, numeric, Cholesky 1243 1244 .seealso: MatCholeskyFactorSymbolic(), MatCholeskyFactor(), MatLUFactorNumeric() 1245 @*/ 1246 int MatCholeskyFactorNumeric(Mat mat,Mat *fact) 1247 { 1248 int ierr; 1249 1250 PetscFunctionBegin; 1251 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1252 PetscValidPointer(fact); 1253 if (!mat->ops->choleskyfactornumeric) SETERRQ(PETSC_ERR_SUP,0,""); 1254 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 1255 if (mat->M != (*fact)->M || mat->N != (*fact)->N) { 1256 SETERRQ4(PETSC_ERR_ARG_SIZ,0,"Mat mat,Mat *fact: global dim %d should = %d %d should = %d", 1257 mat->M,(*fact)->M,mat->N,(*fact)->N); 1258 } 1259 1260 PLogEventBegin(MAT_CholeskyFactorNumeric,mat,*fact,0,0); 1261 ierr = (*mat->ops->choleskyfactornumeric)(mat,fact); CHKERRQ(ierr); 1262 PLogEventEnd(MAT_CholeskyFactorNumeric,mat,*fact,0,0); 1263 PetscFunctionReturn(0); 1264 } 1265 1266 /* ----------------------------------------------------------------*/ 1267 #undef __FUNC__ 1268 #define __FUNC__ "MatSolve" 1269 /*@ 1270 MatSolve - Solves A x = b, given a factored matrix. 1271 1272 Collective on Mat and Vec 1273 1274 Input Parameters: 1275 + mat - the factored matrix 1276 - b - the right-hand-side vector 1277 1278 Output Parameter: 1279 . x - the result vector 1280 1281 Notes: 1282 The vectors b and x cannot be the same. I.e., one cannot 1283 call MatSolve(A,x,x). 1284 1285 .keywords: matrix, linear system, solve, LU, Cholesky, triangular solve 1286 1287 .seealso: MatSolveAdd(), MatSolveTrans(), MatSolveTransAdd() 1288 @*/ 1289 int MatSolve(Mat mat,Vec b,Vec x) 1290 { 1291 int ierr; 1292 1293 PetscFunctionBegin; 1294 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1295 PetscValidHeaderSpecific(b,VEC_COOKIE); PetscValidHeaderSpecific(x,VEC_COOKIE); 1296 if (x == b) SETERRQ(PETSC_ERR_ARG_IDN,0,"x and b must be different vectors"); 1297 if (!mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Unfactored matrix"); 1298 if (mat->N != x->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim %d %d",mat->N,x->N); 1299 if (mat->M != b->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: global dim %d %d",mat->M,b->N); 1300 if (mat->m != b->n) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: local dim %d %d",mat->m,b->n); 1301 if (mat->M == 0 && mat->N == 0) PetscFunctionReturn(0); 1302 1303 if (!mat->ops->solve) SETERRQ(PETSC_ERR_SUP,0,""); 1304 PLogEventBegin(MAT_Solve,mat,b,x,0); 1305 ierr = (*mat->ops->solve)(mat,b,x); CHKERRQ(ierr); 1306 PLogEventEnd(MAT_Solve,mat,b,x,0); 1307 PetscFunctionReturn(0); 1308 } 1309 1310 #undef __FUNC__ 1311 #define __FUNC__ "MatForwardSolve" 1312 /* @ 1313 MatForwardSolve - Solves L x = b, given a factored matrix, A = LU. 1314 1315 Collective on Mat and Vec 1316 1317 Input Parameters: 1318 + mat - the factored matrix 1319 - b - the right-hand-side vector 1320 1321 Output Parameter: 1322 . x - the result vector 1323 1324 Notes: 1325 MatSolve() should be used for most applications, as it performs 1326 a forward solve followed by a backward solve. 1327 1328 The vectors b and x cannot be the same. I.e., one cannot 1329 call MatForwardSolve(A,x,x). 1330 1331 .keywords: matrix, forward, LU, Cholesky, triangular solve 1332 1333 .seealso: MatSolve(), MatBackwardSolve() 1334 @ */ 1335 int MatForwardSolve(Mat mat,Vec b,Vec x) 1336 { 1337 int ierr; 1338 1339 PetscFunctionBegin; 1340 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1341 PetscValidHeaderSpecific(b,VEC_COOKIE); PetscValidHeaderSpecific(x,VEC_COOKIE); 1342 if (x == b) SETERRQ(PETSC_ERR_ARG_IDN,0,"x and b must be different vectors"); 1343 if (!mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Unfactored matrix"); 1344 if (!mat->ops->forwardsolve) SETERRQ(PETSC_ERR_SUP,0,""); 1345 if (mat->N != x->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim %d %d",mat->N,x->N); 1346 if (mat->M != b->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: global dim %d %d",mat->M,b->N); 1347 if (mat->m != b->n) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: local dim %d %d",mat->m,b->n); 1348 1349 PLogEventBegin(MAT_ForwardSolve,mat,b,x,0); 1350 ierr = (*mat->ops->forwardsolve)(mat,b,x); CHKERRQ(ierr); 1351 PLogEventEnd(MAT_ForwardSolve,mat,b,x,0); 1352 PetscFunctionReturn(0); 1353 } 1354 1355 #undef __FUNC__ 1356 #define __FUNC__ "MatBackwardSolve" 1357 /* @ 1358 MatBackwardSolve - Solves U x = b, given a factored matrix, A = LU. 1359 1360 Collective on Mat and Vec 1361 1362 Input Parameters: 1363 + mat - the factored matrix 1364 - b - the right-hand-side vector 1365 1366 Output Parameter: 1367 . x - the result vector 1368 1369 Notes: 1370 MatSolve() should be used for most applications, as it performs 1371 a forward solve followed by a backward solve. 1372 1373 The vectors b and x cannot be the same. I.e., one cannot 1374 call MatBackwardSolve(A,x,x). 1375 1376 .keywords: matrix, backward, LU, Cholesky, triangular solve 1377 1378 .seealso: MatSolve(), MatForwardSolve() 1379 @ */ 1380 int MatBackwardSolve(Mat mat,Vec b,Vec x) 1381 { 1382 int ierr; 1383 1384 PetscFunctionBegin; 1385 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1386 PetscValidHeaderSpecific(b,VEC_COOKIE); PetscValidHeaderSpecific(x,VEC_COOKIE); 1387 if (x == b) SETERRQ(PETSC_ERR_ARG_IDN,0,"x and b must be different vectors"); 1388 if (!mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Unfactored matrix"); 1389 if (!mat->ops->backwardsolve) SETERRQ(PETSC_ERR_SUP,0,""); 1390 if (mat->N != x->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim %d %d",mat->N,x->N); 1391 if (mat->M != b->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: global dim %d %d",mat->M,b->N); 1392 if (mat->m != b->n) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: local dim %d %d",mat->m,b->n); 1393 1394 PLogEventBegin(MAT_BackwardSolve,mat,b,x,0); 1395 ierr = (*mat->ops->backwardsolve)(mat,b,x); CHKERRQ(ierr); 1396 PLogEventEnd(MAT_BackwardSolve,mat,b,x,0); 1397 PetscFunctionReturn(0); 1398 } 1399 1400 #undef __FUNC__ 1401 #define __FUNC__ "MatSolveAdd" 1402 /*@ 1403 MatSolveAdd - Computes x = y + inv(A)*b, given a factored matrix. 1404 1405 Collective on Mat and Vec 1406 1407 Input Parameters: 1408 + mat - the factored matrix 1409 . b - the right-hand-side vector 1410 - y - the vector to be added to 1411 1412 Output Parameter: 1413 . x - the result vector 1414 1415 Notes: 1416 The vectors b and x cannot be the same. I.e., one cannot 1417 call MatSolveAdd(A,x,y,x). 1418 1419 .keywords: matrix, linear system, solve, LU, Cholesky, add 1420 1421 .seealso: MatSolve(), MatSolveTrans(), MatSolveTransAdd() 1422 @*/ 1423 int MatSolveAdd(Mat mat,Vec b,Vec y,Vec x) 1424 { 1425 Scalar one = 1.0; 1426 Vec tmp; 1427 int ierr; 1428 1429 PetscFunctionBegin; 1430 PetscValidHeaderSpecific(mat,MAT_COOKIE);PetscValidHeaderSpecific(y,VEC_COOKIE); 1431 PetscValidHeaderSpecific(b,VEC_COOKIE); PetscValidHeaderSpecific(x,VEC_COOKIE); 1432 if (x == b) SETERRQ(PETSC_ERR_ARG_IDN,0,"x and b must be different vectors"); 1433 if (!mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Unfactored matrix"); 1434 if (mat->N != x->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim %d %d",mat->N,x->N); 1435 if (mat->M != b->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: global dim %d %d",mat->M,b->N); 1436 if (mat->M != y->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec y: global dim %d %d",mat->M,y->N); 1437 if (mat->m != b->n) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: local dim %d %d",mat->m,b->n); 1438 if (x->n != y->n) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Vec x,Vec y: local dim %d %d",x->n,y->n); 1439 1440 PLogEventBegin(MAT_SolveAdd,mat,b,x,y); 1441 if (mat->ops->solveadd) { 1442 ierr = (*mat->ops->solveadd)(mat,b,y,x); CHKERRQ(ierr); 1443 } else { 1444 /* do the solve then the add manually */ 1445 if (x != y) { 1446 ierr = MatSolve(mat,b,x); CHKERRQ(ierr); 1447 ierr = VecAXPY(&one,y,x); CHKERRQ(ierr); 1448 } else { 1449 ierr = VecDuplicate(x,&tmp); CHKERRQ(ierr); 1450 PLogObjectParent(mat,tmp); 1451 ierr = VecCopy(x,tmp); CHKERRQ(ierr); 1452 ierr = MatSolve(mat,b,x); CHKERRQ(ierr); 1453 ierr = VecAXPY(&one,tmp,x); CHKERRQ(ierr); 1454 ierr = VecDestroy(tmp); CHKERRQ(ierr); 1455 } 1456 } 1457 PLogEventEnd(MAT_SolveAdd,mat,b,x,y); 1458 PetscFunctionReturn(0); 1459 } 1460 1461 #undef __FUNC__ 1462 #define __FUNC__ "MatSolveTrans" 1463 /*@ 1464 MatSolveTrans - Solves A' x = b, given a factored matrix. 1465 1466 Collective on Mat and Vec 1467 1468 Input Parameters: 1469 + mat - the factored matrix 1470 - b - the right-hand-side vector 1471 1472 Output Parameter: 1473 . x - the result vector 1474 1475 Notes: 1476 The vectors b and x cannot be the same. I.e., one cannot 1477 call MatSolveTrans(A,x,x). 1478 1479 .keywords: matrix, linear system, solve, LU, Cholesky, transpose 1480 1481 .seealso: MatSolve(), MatSolveAdd(), MatSolveTransAdd() 1482 @*/ 1483 int MatSolveTrans(Mat mat,Vec b,Vec x) 1484 { 1485 int ierr; 1486 1487 PetscFunctionBegin; 1488 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1489 PetscValidHeaderSpecific(b,VEC_COOKIE); PetscValidHeaderSpecific(x,VEC_COOKIE); 1490 if (!mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Unfactored matrix"); 1491 if (x == b) SETERRQ(PETSC_ERR_ARG_IDN,0,"x and b must be different vectors"); 1492 if (!mat->ops->solvetrans) SETERRQ(PETSC_ERR_SUP,0,""); 1493 if (mat->M != x->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim %d %d",mat->M,x->N); 1494 if (mat->N != b->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: global dim %d %d",mat->N,b->N); 1495 1496 PLogEventBegin(MAT_SolveTrans,mat,b,x,0); 1497 ierr = (*mat->ops->solvetrans)(mat,b,x); CHKERRQ(ierr); 1498 PLogEventEnd(MAT_SolveTrans,mat,b,x,0); 1499 PetscFunctionReturn(0); 1500 } 1501 1502 #undef __FUNC__ 1503 #define __FUNC__ "MatSolveTransAdd" 1504 /*@ 1505 MatSolveTransAdd - Computes x = y + inv(trans(A)) b, given a 1506 factored matrix. 1507 1508 Collective on Mat and Vec 1509 1510 Input Parameters: 1511 + mat - the factored matrix 1512 . b - the right-hand-side vector 1513 - y - the vector to be added to 1514 1515 Output Parameter: 1516 . x - the result vector 1517 1518 Notes: 1519 The vectors b and x cannot be the same. I.e., one cannot 1520 call MatSolveTransAdd(A,x,y,x). 1521 1522 .keywords: matrix, linear system, solve, LU, Cholesky, transpose, add 1523 1524 .seealso: MatSolve(), MatSolveAdd(), MatSolveTrans() 1525 @*/ 1526 int MatSolveTransAdd(Mat mat,Vec b,Vec y,Vec x) 1527 { 1528 Scalar one = 1.0; 1529 int ierr; 1530 Vec tmp; 1531 1532 PetscFunctionBegin; 1533 PetscValidHeaderSpecific(mat,MAT_COOKIE);PetscValidHeaderSpecific(y,VEC_COOKIE); 1534 PetscValidHeaderSpecific(b,VEC_COOKIE); PetscValidHeaderSpecific(x,VEC_COOKIE); 1535 if (x == b) SETERRQ(PETSC_ERR_ARG_IDN,0,"x and b must be different vectors"); 1536 if (!mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Unfactored matrix"); 1537 if (mat->M != x->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim %d %d",mat->M,x->N); 1538 if (mat->N != b->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: global dim %d %d",mat->N,b->N); 1539 if (mat->N != y->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec y: global dim %d %d",mat->N,y->N); 1540 if (x->n != y->n) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Vec x,Vec y: local dim %d %d",x->n,y->n); 1541 1542 PLogEventBegin(MAT_SolveTransAdd,mat,b,x,y); 1543 if (mat->ops->solvetransadd) { 1544 ierr = (*mat->ops->solvetransadd)(mat,b,y,x); CHKERRQ(ierr); 1545 } else { 1546 /* do the solve then the add manually */ 1547 if (x != y) { 1548 ierr = MatSolveTrans(mat,b,x); CHKERRQ(ierr); 1549 ierr = VecAXPY(&one,y,x); CHKERRQ(ierr); 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) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec x: global dim %d %d",mat->N,x->N); 1621 if (mat->M != b->N) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: global dim %d %d",mat->M,b->N); 1622 if (mat->m != b->n) SETERRQ2(PETSC_ERR_ARG_SIZ,0,"Mat mat,Vec b: local dim %d %d",mat->m,b->n); 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) SETERRQ4(PETSC_ERR_ARG_SIZ,0,"Mat A,Mat B: global dim %d %d",A->M,B->M, 1685 A->N,B->N); 1686 1687 PLogEventBegin(MAT_Copy,A,B,0,0); 1688 if (A->ops->copy) { 1689 ierr = (*A->ops->copy)(A,B,str); CHKERRQ(ierr); 1690 } else { /* generic conversion */ 1691 ierr = MatCopy_Basic(A,B,str); CHKERRQ(ierr); 1692 } 1693 PLogEventEnd(MAT_Copy,A,B,0,0); 1694 PetscFunctionReturn(0); 1695 } 1696 1697 static int MatConvertersSet = 0; 1698 static int (*MatConverters[MAX_MATRIX_TYPES][MAX_MATRIX_TYPES])(Mat,MatType,Mat*) = 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 {0,0,0,0,0,0,0,0,0,0,0,0},{0,0,0,0,0,0,0,0,0,0,0,0}}; 1705 1706 #undef __FUNC__ 1707 #define __FUNC__ "MatConvertRegister" 1708 /*@C 1709 MatConvertRegister - Allows one to register a routine that converts between 1710 two matrix types. 1711 1712 Not Collective 1713 1714 Input Parameters: 1715 . intype - the type of matrix (defined in include/mat.h), for example, MATSEQAIJ. 1716 . outtype - new matrix type, or MATSAME 1717 1718 .seealso: MatConvertRegisterAll() 1719 @*/ 1720 int MatConvertRegister(MatType intype,MatType outtype,int (*converter)(Mat,MatType,Mat*)) 1721 { 1722 PetscFunctionBegin; 1723 MatConverters[intype][outtype] = converter; 1724 MatConvertersSet = 1; 1725 PetscFunctionReturn(0); 1726 } 1727 1728 #undef __FUNC__ 1729 #define __FUNC__ "MatConvert" 1730 /*@C 1731 MatConvert - Converts a matrix to another matrix, either of the same 1732 or different type. 1733 1734 Collective on Mat 1735 1736 Input Parameters: 1737 + mat - the matrix 1738 - newtype - new matrix type. Use MATSAME to create a new matrix of the 1739 same type as the original matrix. 1740 1741 Output Parameter: 1742 . M - pointer to place new matrix 1743 1744 Notes: 1745 MatConvert() first creates a new matrix and then copies the data from 1746 the first matrix. A related routine is MatCopy(), which copies the matrix 1747 entries of one matrix to another already existing matrix context. 1748 1749 .keywords: matrix, copy, convert 1750 1751 .seealso: MatCopy(), MatDuplicate() 1752 @*/ 1753 int MatConvert(Mat mat,MatType newtype,Mat *M) 1754 { 1755 int ierr; 1756 1757 PetscFunctionBegin; 1758 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1759 PetscValidPointer(M); 1760 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 1761 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 1762 1763 if (newtype > MAX_MATRIX_TYPES || newtype < -1) { 1764 SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,1,"Not a valid matrix type"); 1765 } 1766 *M = 0; 1767 1768 if (!MatConvertersSet) { 1769 ierr = MatLoadRegisterAll(); CHKERRQ(ierr); 1770 } 1771 1772 PLogEventBegin(MAT_Convert,mat,0,0,0); 1773 if ((newtype == mat->type || newtype == MATSAME) && mat->ops->duplicate) { 1774 ierr = (*mat->ops->duplicate)(mat,MAT_COPY_VALUES,M); CHKERRQ(ierr); 1775 } else { 1776 if (!MatConvertersSet) { 1777 ierr = MatConvertRegisterAll(); CHKERRQ(ierr); 1778 } 1779 if (!MatConverters[mat->type][newtype]) { 1780 SETERRQ(PETSC_ERR_ARG_WRONG,1,"Invalid matrix type, or matrix converter not registered"); 1781 } 1782 ierr = (*MatConverters[mat->type][newtype])(mat,newtype,M); CHKERRQ(ierr); 1783 } 1784 PLogEventEnd(MAT_Convert,mat,0,0,0); 1785 PetscFunctionReturn(0); 1786 } 1787 1788 #undef __FUNC__ 1789 #define __FUNC__ "MatDuplicate" 1790 /*@C 1791 MatDuplicate - Duplicates a matrix including the non-zero structure. 1792 1793 Collective on Mat 1794 1795 Input Parameters: 1796 + mat - the matrix 1797 - op - either MAT_DO_NO_COPY_VALUES or MAT_COPY_VALUES, cause it to copy nonzero 1798 values as well or not 1799 1800 Output Parameter: 1801 . M - pointer to place new matrix 1802 1803 .keywords: matrix, copy, convert, duplicate 1804 1805 .seealso: MatCopy(), MatConvert() 1806 @*/ 1807 int MatDuplicate(Mat mat,MatDuplicateOption op,Mat *M) 1808 { 1809 int ierr; 1810 1811 PetscFunctionBegin; 1812 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1813 PetscValidPointer(M); 1814 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 1815 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 1816 1817 *M = 0; 1818 PLogEventBegin(MAT_Convert,mat,0,0,0); 1819 if (!mat->ops->duplicate) { 1820 SETERRQ(PETSC_ERR_SUP,1,"Not written for this matrix type"); 1821 } 1822 ierr = (*mat->ops->duplicate)(mat,op,M); CHKERRQ(ierr); 1823 PLogEventEnd(MAT_Convert,mat,0,0,0); 1824 PetscFunctionReturn(0); 1825 } 1826 1827 #undef __FUNC__ 1828 #define __FUNC__ "MatGetDiagonal" 1829 /*@ 1830 MatGetDiagonal - Gets the diagonal of a matrix. 1831 1832 Collective on Mat and Vec 1833 1834 Input Parameters: 1835 + mat - the matrix 1836 - v - the vector for storing the diagonal 1837 1838 Output Parameter: 1839 . v - the diagonal of the matrix 1840 1841 Notes: 1842 For the SeqAIJ matrix format, this routine may also be called 1843 on a LU factored matrix; in that case it routines the reciprocal of 1844 the diagonal entries in U. It returns the entries permuted by the 1845 row and column permutation used during the symbolic factorization. 1846 1847 .keywords: matrix, get, diagonal 1848 @*/ 1849 int MatGetDiagonal(Mat mat,Vec v) 1850 { 1851 int ierr; 1852 1853 PetscFunctionBegin; 1854 PetscValidHeaderSpecific(mat,MAT_COOKIE);PetscValidHeaderSpecific(v,VEC_COOKIE); 1855 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 1856 if (!mat->ops->getdiagonal) SETERRQ(PETSC_ERR_SUP,0,""); 1857 ierr = (*mat->ops->getdiagonal)(mat,v);CHKERRQ(ierr); 1858 PetscFunctionReturn(0); 1859 } 1860 1861 #undef __FUNC__ 1862 #define __FUNC__ "MatTranspose" 1863 /*@C 1864 MatTranspose - Computes an in-place or out-of-place transpose of a matrix. 1865 1866 Collective on Mat 1867 1868 Input Parameter: 1869 . mat - the matrix to transpose 1870 1871 Output Parameters: 1872 . B - the transpose (or pass in PETSC_NULL for an in-place transpose) 1873 1874 .keywords: matrix, transpose 1875 1876 .seealso: MatMultTrans(), MatMultTransAdd() 1877 @*/ 1878 int MatTranspose(Mat mat,Mat *B) 1879 { 1880 int ierr; 1881 1882 PetscFunctionBegin; 1883 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1884 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 1885 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 1886 if (!mat->ops->transpose) SETERRQ(PETSC_ERR_SUP,0,""); 1887 ierr = (*mat->ops->transpose)(mat,B);CHKERRQ(ierr); 1888 PetscFunctionReturn(0); 1889 } 1890 1891 #undef __FUNC__ 1892 #define __FUNC__ "MatPermute" 1893 /*@C 1894 MatPermute - Creates a new matrix with rows and columns permuted from the 1895 original. 1896 1897 Collective on Mat 1898 1899 Input Parameters: 1900 + mat - the matrix to permute 1901 . row - row permutation 1902 - col - column permutation 1903 1904 Output Parameters: 1905 . B - the permuted matrix 1906 1907 .keywords: matrix, transpose 1908 1909 .seealso: MatGetReordering() 1910 @*/ 1911 int MatPermute(Mat mat,IS row,IS col,Mat *B) 1912 { 1913 int ierr; 1914 1915 PetscFunctionBegin; 1916 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1917 PetscValidHeaderSpecific(row,IS_COOKIE); 1918 PetscValidHeaderSpecific(col,IS_COOKIE); 1919 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 1920 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 1921 if (!mat->ops->permute) SETERRQ(PETSC_ERR_SUP,0,""); 1922 ierr = (*mat->ops->permute)(mat,row,col,B);CHKERRQ(ierr); 1923 PetscFunctionReturn(0); 1924 } 1925 1926 #undef __FUNC__ 1927 #define __FUNC__ "MatEqual" 1928 /*@ 1929 MatEqual - Compares two matrices. 1930 1931 Collective on Mat 1932 1933 Input Parameters: 1934 + A - the first matrix 1935 - B - the second matrix 1936 1937 Output Parameter: 1938 . flg - PETSC_TRUE if the matrices are equal; PETSC_FALSE otherwise. 1939 1940 .keywords: matrix, equal, equivalent 1941 @*/ 1942 int MatEqual(Mat A,Mat B,PetscTruth *flg) 1943 { 1944 int ierr; 1945 1946 PetscFunctionBegin; 1947 PetscValidHeaderSpecific(A,MAT_COOKIE); PetscValidHeaderSpecific(B,MAT_COOKIE); 1948 PetscValidIntPointer(flg); 1949 if (!A->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 1950 if (!B->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 1951 if (A->M != B->M || A->N != B->N) SETERRQ4(PETSC_ERR_ARG_SIZ,0,"Mat A,Mat B: global dim %d %d %d %d", 1952 A->M,B->M,A->N,B->N); 1953 if (!A->ops->equal) SETERRQ(PETSC_ERR_SUP,0,""); 1954 ierr = (*A->ops->equal)(A,B,flg);CHKERRQ(ierr); 1955 PetscFunctionReturn(0); 1956 } 1957 1958 #undef __FUNC__ 1959 #define __FUNC__ "MatDiagonalScale" 1960 /*@ 1961 MatDiagonalScale - Scales a matrix on the left and right by diagonal 1962 matrices that are stored as vectors. Either of the two scaling 1963 matrices can be PETSC_NULL. 1964 1965 Collective on Mat 1966 1967 Input Parameters: 1968 + mat - the matrix to be scaled 1969 . l - the left scaling vector (or PETSC_NULL) 1970 - r - the right scaling vector (or PETSC_NULL) 1971 1972 Notes: 1973 MatDiagonalScale() computes A = LAR, where 1974 L = a diagonal matrix, R = a diagonal matrix 1975 1976 .keywords: matrix, diagonal, scale 1977 1978 .seealso: MatDiagonalScale() 1979 @*/ 1980 int MatDiagonalScale(Mat mat,Vec l,Vec r) 1981 { 1982 int ierr; 1983 1984 PetscFunctionBegin; 1985 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1986 if (!mat->ops->diagonalscale) SETERRQ(PETSC_ERR_SUP,0,""); 1987 if (l) PetscValidHeaderSpecific(l,VEC_COOKIE); 1988 if (r) PetscValidHeaderSpecific(r,VEC_COOKIE); 1989 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 1990 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 1991 1992 PLogEventBegin(MAT_Scale,mat,0,0,0); 1993 ierr = (*mat->ops->diagonalscale)(mat,l,r); CHKERRQ(ierr); 1994 PLogEventEnd(MAT_Scale,mat,0,0,0); 1995 PetscFunctionReturn(0); 1996 } 1997 1998 #undef __FUNC__ 1999 #define __FUNC__ "MatScale" 2000 /*@ 2001 MatScale - Scales all elements of a matrix by a given number. 2002 2003 Collective on Mat 2004 2005 Input Parameters: 2006 + mat - the matrix to be scaled 2007 - a - the scaling value 2008 2009 Output Parameter: 2010 . mat - the scaled matrix 2011 2012 .keywords: matrix, scale 2013 2014 .seealso: MatDiagonalScale() 2015 @*/ 2016 int MatScale(Scalar *a,Mat mat) 2017 { 2018 int ierr; 2019 2020 PetscFunctionBegin; 2021 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2022 PetscValidScalarPointer(a); 2023 if (!mat->ops->scale) SETERRQ(PETSC_ERR_SUP,0,""); 2024 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2025 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2026 2027 PLogEventBegin(MAT_Scale,mat,0,0,0); 2028 ierr = (*mat->ops->scale)(a,mat); CHKERRQ(ierr); 2029 PLogEventEnd(MAT_Scale,mat,0,0,0); 2030 PetscFunctionReturn(0); 2031 } 2032 2033 #undef __FUNC__ 2034 #define __FUNC__ "MatNorm" 2035 /*@ 2036 MatNorm - Calculates various norms of a matrix. 2037 2038 Collective on Mat 2039 2040 Input Parameters: 2041 + mat - the matrix 2042 - type - the type of norm, NORM_1, NORM_2, NORM_FROBENIUS, NORM_INFINITY 2043 2044 Output Parameters: 2045 . norm - the resulting norm 2046 2047 .keywords: matrix, norm, Frobenius 2048 @*/ 2049 int MatNorm(Mat mat,NormType type,double *norm) 2050 { 2051 int ierr; 2052 2053 PetscFunctionBegin; 2054 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2055 PetscValidScalarPointer(norm); 2056 2057 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2058 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2059 if (!mat->ops->norm) SETERRQ(PETSC_ERR_SUP,0,"Not for this matrix type"); 2060 ierr = (*mat->ops->norm)(mat,type,norm);CHKERRQ(ierr); 2061 PetscFunctionReturn(0); 2062 } 2063 2064 /* 2065 This variable is used to prevent counting of MatAssemblyBegin() that 2066 are called from within a MatAssemblyEnd(). 2067 */ 2068 static int MatAssemblyEnd_InUse = 0; 2069 #undef __FUNC__ 2070 #define __FUNC__ "MatAssemblyBegin" 2071 /*@ 2072 MatAssemblyBegin - Begins assembling the matrix. This routine should 2073 be called after completing all calls to MatSetValues(). 2074 2075 Collective on Mat 2076 2077 Input Parameters: 2078 + mat - the matrix 2079 - type - type of assembly, either MAT_FLUSH_ASSEMBLY or MAT_FINAL_ASSEMBLY 2080 2081 Notes: 2082 MatSetValues() generally caches the values. The matrix is ready to 2083 use only after MatAssemblyBegin() and MatAssemblyEnd() have been called. 2084 Use MAT_FLUSH_ASSEMBLY when switching between ADD_VALUES and INSERT_VALUES 2085 in MatSetValues(); use MAT_FINAL_ASSEMBLY for the final assembly before 2086 using the matrix. 2087 2088 .keywords: matrix, assembly, assemble, begin 2089 2090 .seealso: MatAssemblyEnd(), MatSetValues() 2091 @*/ 2092 int MatAssemblyBegin(Mat mat,MatAssemblyType type) 2093 { 2094 int ierr; 2095 2096 PetscFunctionBegin; 2097 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2098 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix.\n did you forget to call MatSetUnfactored()?"); 2099 if (mat->assembled) { 2100 mat->was_assembled = PETSC_TRUE; 2101 mat->assembled = PETSC_FALSE; 2102 } 2103 if (!MatAssemblyEnd_InUse) { 2104 PLogEventBegin(MAT_AssemblyBegin,mat,0,0,0); 2105 if (mat->ops->assemblybegin){ierr = (*mat->ops->assemblybegin)(mat,type);CHKERRQ(ierr);} 2106 PLogEventEnd(MAT_AssemblyBegin,mat,0,0,0); 2107 } else { 2108 if (mat->ops->assemblybegin){ierr = (*mat->ops->assemblybegin)(mat,type);CHKERRQ(ierr);} 2109 } 2110 PetscFunctionReturn(0); 2111 } 2112 2113 2114 #undef __FUNC__ 2115 #define __FUNC__ "MatView_Private" 2116 /* 2117 Processes command line options to determine if/how a matrix 2118 is to be viewed. Called by MatAssemblyEnd() and MatLoad(). 2119 */ 2120 int MatView_Private(Mat mat) 2121 { 2122 int ierr,flg; 2123 2124 PetscFunctionBegin; 2125 ierr = OptionsHasName(PETSC_NULL,"-mat_view_info",&flg); CHKERRQ(ierr); 2126 if (flg) { 2127 ierr = ViewerPushFormat(VIEWER_STDOUT_(mat->comm),VIEWER_FORMAT_ASCII_INFO,0);CHKERRQ(ierr); 2128 ierr = MatView(mat,VIEWER_STDOUT_(mat->comm)); CHKERRQ(ierr); 2129 ierr = ViewerPopFormat(VIEWER_STDOUT_(mat->comm));CHKERRQ(ierr); 2130 } 2131 ierr = OptionsHasName(PETSC_NULL,"-mat_view_info_detailed",&flg);CHKERRQ(ierr); 2132 if (flg) { 2133 ierr = ViewerPushFormat(VIEWER_STDOUT_(mat->comm),VIEWER_FORMAT_ASCII_INFO_LONG,0);CHKERRQ(ierr); 2134 ierr = MatView(mat,VIEWER_STDOUT_(mat->comm)); CHKERRQ(ierr); 2135 ierr = ViewerPopFormat(VIEWER_STDOUT_(mat->comm));CHKERRQ(ierr); 2136 } 2137 ierr = OptionsHasName(PETSC_NULL,"-mat_view",&flg); CHKERRQ(ierr); 2138 if (flg) { 2139 ierr = MatView(mat,VIEWER_STDOUT_(mat->comm)); CHKERRQ(ierr); 2140 } 2141 ierr = OptionsHasName(PETSC_NULL,"-mat_view_matlab",&flg); CHKERRQ(ierr); 2142 if (flg) { 2143 ierr = ViewerPushFormat(VIEWER_STDOUT_(mat->comm),VIEWER_FORMAT_ASCII_MATLAB,"M");CHKERRQ(ierr); 2144 ierr = MatView(mat,VIEWER_STDOUT_(mat->comm)); CHKERRQ(ierr); 2145 ierr = ViewerPopFormat(VIEWER_STDOUT_(mat->comm));CHKERRQ(ierr); 2146 } 2147 ierr = OptionsHasName(PETSC_NULL,"-mat_view_draw",&flg); CHKERRQ(ierr); 2148 if (flg) { 2149 ierr = OptionsHasName(0,"-mat_view_contour",&flg); CHKERRQ(ierr); 2150 if (flg) { 2151 ViewerPushFormat(VIEWER_DRAWX_(mat->comm),VIEWER_FORMAT_DRAW_CONTOUR,0);CHKERRQ(ierr); 2152 } 2153 ierr = MatView(mat,VIEWER_DRAWX_(mat->comm)); CHKERRQ(ierr); 2154 ierr = ViewerFlush(VIEWER_DRAWX_(mat->comm)); CHKERRQ(ierr); 2155 if (flg) { 2156 ViewerPopFormat(VIEWER_DRAWX_(mat->comm));CHKERRQ(ierr); 2157 } 2158 } 2159 PetscFunctionReturn(0); 2160 } 2161 2162 #undef __FUNC__ 2163 #define __FUNC__ "MatAssemblyEnd" 2164 /*@ 2165 MatAssemblyEnd - Completes assembling the matrix. This routine should 2166 be called after MatAssemblyBegin(). 2167 2168 Collective on Mat 2169 2170 Input Parameters: 2171 + mat - the matrix 2172 - type - type of assembly, either MAT_FLUSH_ASSEMBLY or MAT_FINAL_ASSEMBLY 2173 2174 Options Database Keys: 2175 + -mat_view_info - Prints info on matrix at conclusion of MatEndAssembly() 2176 . -mat_view_info_detailed - Prints more detailed info 2177 . -mat_view - Prints matrix in ASCII format 2178 . -mat_view_matlab - Prints matrix in Matlab format 2179 . -mat_view_draw - Draws nonzero structure of matrix, using MatView() and DrawOpenX(). 2180 . -display <name> - Sets display name (default is host) 2181 - -draw_pause <sec> - Sets number of seconds to pause after display 2182 2183 Notes: 2184 MatSetValues() generally caches the values. The matrix is ready to 2185 use only after MatAssemblyBegin() and MatAssemblyEnd() have been called. 2186 Use MAT_FLUSH_ASSEMBLY when switching between ADD_VALUES and INSERT_VALUES 2187 in MatSetValues(); use MAT_FINAL_ASSEMBLY for the final assembly before 2188 using the matrix. 2189 2190 .keywords: matrix, assembly, assemble, end 2191 2192 .seealso: MatAssemblyBegin(), MatSetValues(), DrawOpenX(), MatView() 2193 @*/ 2194 int MatAssemblyEnd(Mat mat,MatAssemblyType type) 2195 { 2196 int ierr; 2197 static int inassm = 0; 2198 2199 PetscFunctionBegin; 2200 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2201 2202 inassm++; 2203 MatAssemblyEnd_InUse++; 2204 PLogEventBegin(MAT_AssemblyEnd,mat,0,0,0); 2205 if (mat->ops->assemblyend) { 2206 ierr = (*mat->ops->assemblyend)(mat,type); CHKERRQ(ierr); 2207 } 2208 2209 /* Flush assembly is not a true assembly */ 2210 if (type != MAT_FLUSH_ASSEMBLY) { 2211 mat->assembled = PETSC_TRUE; mat->num_ass++; 2212 } 2213 mat->insertmode = NOT_SET_VALUES; 2214 PLogEventEnd(MAT_AssemblyEnd,mat,0,0,0); 2215 MatAssemblyEnd_InUse--; 2216 2217 if (inassm == 1 && type != MAT_FLUSH_ASSEMBLY) { 2218 ierr = MatView_Private(mat); CHKERRQ(ierr); 2219 } 2220 inassm--; 2221 PetscFunctionReturn(0); 2222 } 2223 2224 2225 #undef __FUNC__ 2226 #define __FUNC__ "MatCompress" 2227 /*@ 2228 MatCompress - Tries to store the matrix in as little space as 2229 possible. May fail if memory is already fully used, since it 2230 tries to allocate new space. 2231 2232 Collective on Mat 2233 2234 Input Parameters: 2235 . mat - the matrix 2236 2237 .keywords: matrix, compress 2238 @*/ 2239 int MatCompress(Mat mat) 2240 { 2241 int ierr; 2242 2243 PetscFunctionBegin; 2244 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2245 if (mat->ops->compress) {ierr = (*mat->ops->compress)(mat);CHKERRQ(ierr);} 2246 PetscFunctionReturn(0); 2247 } 2248 2249 #undef __FUNC__ 2250 #define __FUNC__ "MatSetOption" 2251 /*@ 2252 MatSetOption - Sets a parameter option for a matrix. Some options 2253 may be specific to certain storage formats. Some options 2254 determine how values will be inserted (or added). Sorted, 2255 row-oriented input will generally assemble the fastest. The default 2256 is row-oriented, nonsorted input. 2257 2258 Collective on Mat 2259 2260 Input Parameters: 2261 + mat - the matrix 2262 - option - the option, one of those listed below (and possibly others), 2263 e.g., MAT_ROWS_SORTED, MAT_NEW_NONZERO_LOCATION_ERROR 2264 2265 Options Describing Matrix Structure: 2266 + MAT_SYMMETRIC - symmetric in terms of both structure and value 2267 - MAT_STRUCTURALLY_SYMMETRIC - symmetric nonzero structure 2268 2269 Options For Use with MatSetValues(): 2270 Insert a logically dense subblock, which can be 2271 + MAT_ROW_ORIENTED - row-oriented 2272 . MAT_COLUMN_ORIENTED - column-oriented 2273 . MAT_ROWS_SORTED - sorted by row 2274 . MAT_ROWS_UNSORTED - not sorted by row 2275 . MAT_COLUMNS_SORTED - sorted by column 2276 - MAT_COLUMNS_UNSORTED - not sorted by column 2277 2278 Not these options reflect the data you pass in with MatSetValues(); it has 2279 nothing to do with how the data is stored internally in the matrix 2280 data structure. 2281 2282 When (re)assembling a matrix, we can restrict the input for 2283 efficiency/debugging purposes. These options include 2284 + MAT_NO_NEW_NONZERO_LOCATIONS - additional insertions will not be 2285 allowed if they generate a new nonzero 2286 . MAT_YES_NEW_NONZERO_LOCATIONS - additional insertions will be allowed 2287 . MAT_NO_NEW_DIAGONALS - additional insertions will not be allowed if 2288 they generate a nonzero in a new diagonal (for block diagonal format only) 2289 . MAT_YES_NEW_DIAGONALS - new diagonals will be allowed (for block diagonal format only) 2290 . MAT_IGNORE_OFF_PROC_ENTRIES - drops off-processor entries 2291 . MAT_NEW_NONZERO_LOCATION_ERROR - generates an error for new matrix entry 2292 - MAT_USE_HASH_TABLE - uses a hash table to speed up matrix assembly 2293 2294 Notes: 2295 Some options are relevant only for particular matrix types and 2296 are thus ignored by others. Other options are not supported by 2297 certain matrix types and will generate an error message if set. 2298 2299 If using a Fortran 77 module to compute a matrix, one may need to 2300 use the column-oriented option (or convert to the row-oriented 2301 format). 2302 2303 MAT_NO_NEW_NONZERO_LOCATIONS indicates that any add or insertion 2304 that would generate a new entry in the nonzero structure is instead 2305 ignored. Thus, if memory has not alredy been allocated for this particular 2306 data, then the insertion is ignored. For dense matrices, in which 2307 the entire array is allocated, no entries are ever ignored. 2308 2309 MAT_NEW_NONZERO_LOCATION_ERROR indicates that any add or insertion 2310 that would generate a new entry in the nonzero structure instead produces 2311 an error. (Currently supported for AIJ and BAIJ formats only.) 2312 This is a useful flag when using SAME_NONZERO_PATTERN in calling 2313 SLESSetOperators() to ensure that the nonzero pattern truely does 2314 remain unchanged. 2315 2316 MAT_NEW_NONZERO_ALLOCATION_ERROR indicates that any add or insertion 2317 that would generate a new entry that has not been preallocated will 2318 instead produce an error. (Currently supported for AIJ and BAIJ formats 2319 only.) This is a useful flag when debugging matrix memory preallocation. 2320 2321 MAT_IGNORE_OFF_PROC_ENTRIES indicates entries destined for 2322 other processors should be dropped, rather than stashed. 2323 This is useful if you know that the "owning" processor is also 2324 always generating the correct matrix entries, so that PETSc need 2325 not transfer duplicate entries generated on another processor. 2326 2327 MAT_USE_HASH_TABLE indicates that a hash table be used to improve the 2328 searches during matrix assembly. When this flag is set, the hash table 2329 is created during the first Matrix Assembly. This hash table is 2330 used the next time through, during MatSetVaules()/MatSetVaulesBlocked() 2331 to improve the searching of indices. MAT_NO_NEW_NONZERO_LOCATIONS flag 2332 should be used with MAT_USE_HASH_TABLE flag. This option is currently 2333 supported by MATMPIBAIJ format only. 2334 2335 .keywords: matrix, option, row-oriented, column-oriented, sorted, nonzero 2336 @*/ 2337 int MatSetOption(Mat mat,MatOption op) 2338 { 2339 int ierr; 2340 2341 PetscFunctionBegin; 2342 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2343 if (mat->ops->setoption) {ierr = (*mat->ops->setoption)(mat,op);CHKERRQ(ierr);} 2344 PetscFunctionReturn(0); 2345 } 2346 2347 #undef __FUNC__ 2348 #define __FUNC__ "MatZeroEntries" 2349 /*@ 2350 MatZeroEntries - Zeros all entries of a matrix. For sparse matrices 2351 this routine retains the old nonzero structure. 2352 2353 Collective on Mat 2354 2355 Input Parameters: 2356 . mat - the matrix 2357 2358 .keywords: matrix, zero, entries 2359 2360 .seealso: MatZeroRows() 2361 @*/ 2362 int MatZeroEntries(Mat mat) 2363 { 2364 int ierr; 2365 2366 PetscFunctionBegin; 2367 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2368 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2369 if (!mat->ops->zeroentries) SETERRQ(PETSC_ERR_SUP,0,""); 2370 2371 PLogEventBegin(MAT_ZeroEntries,mat,0,0,0); 2372 ierr = (*mat->ops->zeroentries)(mat); CHKERRQ(ierr); 2373 PLogEventEnd(MAT_ZeroEntries,mat,0,0,0); 2374 PetscFunctionReturn(0); 2375 } 2376 2377 #undef __FUNC__ 2378 #define __FUNC__ "MatZeroRows" 2379 /*@ 2380 MatZeroRows - Zeros all entries (except possibly the main diagonal) 2381 of a set of rows of a matrix. 2382 2383 Collective on Mat 2384 2385 Input Parameters: 2386 + mat - the matrix 2387 . is - index set of rows to remove 2388 - diag - pointer to value put in all diagonals of eliminated rows. 2389 Note that diag is not a pointer to an array, but merely a 2390 pointer to a single value. 2391 2392 Notes: 2393 For the AIJ matrix formats this removes the old nonzero structure, 2394 but does not release memory. For the dense and block diagonal 2395 formats this does not alter the nonzero structure. 2396 2397 The user can set a value in the diagonal entry (or for the AIJ and 2398 row formats can optionally remove the main diagonal entry from the 2399 nonzero structure as well, by passing a null pointer as the final 2400 argument). 2401 2402 For the parallel case, all processes that share the matrix (i.e., 2403 those in the communicator used for matrix creation) MUST call this 2404 routine, regardless of whether any rows being zeroed are owned by 2405 them. 2406 2407 .keywords: matrix, zero, rows, boundary conditions 2408 2409 .seealso: MatZeroEntries(), 2410 @*/ 2411 int MatZeroRows(Mat mat,IS is, Scalar *diag) 2412 { 2413 int ierr; 2414 2415 PetscFunctionBegin; 2416 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2417 PetscValidHeaderSpecific(is,IS_COOKIE); 2418 if (diag) PetscValidScalarPointer(diag); 2419 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2420 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2421 if (!mat->ops->zerorows) SETERRQ(PETSC_ERR_SUP,0,""); 2422 2423 ierr = (*mat->ops->zerorows)(mat,is,diag); CHKERRQ(ierr); 2424 ierr = MatView_Private(mat); CHKERRQ(ierr); 2425 PetscFunctionReturn(0); 2426 } 2427 2428 #undef __FUNC__ 2429 #define __FUNC__ "MatZeroRowsLocal" 2430 /*@ 2431 MatZeroRowsLocal - Zeros all entries (except possibly the main diagonal) 2432 of a set of rows of a matrix; using local numbering of rows. 2433 2434 Collective on Mat 2435 2436 Input Parameters: 2437 + mat - the matrix 2438 . is - index set of rows to remove 2439 - diag - pointer to value put in all diagonals of eliminated rows. 2440 Note that diag is not a pointer to an array, but merely a 2441 pointer to a single value. 2442 2443 Notes: 2444 For the AIJ matrix formats this removes the old nonzero structure, 2445 but does not release memory. For the dense and block diagonal 2446 formats this does not alter the nonzero structure. 2447 2448 The user can set a value in the diagonal entry (or for the AIJ and 2449 row formats can optionally remove the main diagonal entry from the 2450 nonzero structure as well, by passing a null pointer as the final 2451 argument). 2452 2453 .keywords: matrix, zero, rows, boundary conditions 2454 2455 .seealso: MatZeroEntries(), 2456 @*/ 2457 int MatZeroRowsLocal(Mat mat,IS is, Scalar *diag) 2458 { 2459 int ierr; 2460 IS newis; 2461 2462 PetscFunctionBegin; 2463 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2464 PetscValidHeaderSpecific(is,IS_COOKIE); 2465 if (diag) PetscValidScalarPointer(diag); 2466 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2467 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2468 if (!mat->ops->zerorows) SETERRQ(PETSC_ERR_SUP,0,""); 2469 2470 ierr = ISLocalToGlobalMappingApplyIS(mat->mapping,is,&newis); CHKERRQ(ierr); 2471 ierr = (*mat->ops->zerorows)(mat,newis,diag); CHKERRQ(ierr); 2472 ierr = ISDestroy(newis); 2473 PetscFunctionReturn(0); 2474 } 2475 2476 #undef __FUNC__ 2477 #define __FUNC__ "MatGetSize" 2478 /*@ 2479 MatGetSize - Returns the numbers of rows and columns in a matrix. 2480 2481 Not Collective 2482 2483 Input Parameter: 2484 . mat - the matrix 2485 2486 Output Parameters: 2487 + m - the number of global rows 2488 - n - the number of global columns 2489 2490 .keywords: matrix, dimension, size, rows, columns, global, get 2491 2492 .seealso: MatGetLocalSize() 2493 @*/ 2494 int MatGetSize(Mat mat,int *m,int* n) 2495 { 2496 int ierr; 2497 2498 PetscFunctionBegin; 2499 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2500 ierr = (*mat->ops->getsize)(mat,m,n);CHKERRQ(ierr); 2501 PetscFunctionReturn(0); 2502 } 2503 2504 #undef __FUNC__ 2505 #define __FUNC__ "MatGetLocalSize" 2506 /*@ 2507 MatGetLocalSize - Returns the number of rows and columns in a matrix 2508 stored locally. This information may be implementation dependent, so 2509 use with care. 2510 2511 Not Collective 2512 2513 Input Parameters: 2514 . mat - the matrix 2515 2516 Output Parameters: 2517 + m - the number of local rows 2518 - n - the number of local columns 2519 2520 .keywords: matrix, dimension, size, local, rows, columns, get 2521 2522 .seealso: MatGetSize() 2523 @*/ 2524 int MatGetLocalSize(Mat mat,int *m,int* n) 2525 { 2526 int ierr; 2527 2528 PetscFunctionBegin; 2529 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2530 ierr = (*mat->ops->getlocalsize)(mat,m,n);CHKERRQ(ierr); 2531 PetscFunctionReturn(0); 2532 } 2533 2534 #undef __FUNC__ 2535 #define __FUNC__ "MatGetOwnershipRange" 2536 /*@ 2537 MatGetOwnershipRange - Returns the range of matrix rows owned by 2538 this processor, assuming that the matrix is laid out with the first 2539 n1 rows on the first processor, the next n2 rows on the second, etc. 2540 For certain parallel layouts this range may not be well defined. 2541 2542 Not Collective 2543 2544 Input Parameters: 2545 . mat - the matrix 2546 2547 Output Parameters: 2548 + m - the global index of the first local row 2549 - n - one more than the global index of the last local row 2550 2551 .keywords: matrix, get, range, ownership 2552 @*/ 2553 int MatGetOwnershipRange(Mat mat,int *m,int* n) 2554 { 2555 int ierr; 2556 2557 PetscFunctionBegin; 2558 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2559 PetscValidIntPointer(m); 2560 PetscValidIntPointer(n); 2561 if (!mat->ops->getownershiprange) SETERRQ(PETSC_ERR_SUP,0,""); 2562 ierr = (*mat->ops->getownershiprange)(mat,m,n);CHKERRQ(ierr); 2563 PetscFunctionReturn(0); 2564 } 2565 2566 #undef __FUNC__ 2567 #define __FUNC__ "MatILUFactorSymbolic" 2568 /*@ 2569 MatILUFactorSymbolic - Performs symbolic ILU factorization of a matrix. 2570 Uses levels of fill only, not drop tolerance. Use MatLUFactorNumeric() 2571 to complete the factorization. 2572 2573 Collective on Mat 2574 2575 Input Parameters: 2576 + mat - the matrix 2577 . row - row permutation 2578 . column - column permutation 2579 . fill - number of levels of fill 2580 - f - expected fill as ratio of the original number of nonzeros, 2581 for example 3.0; choosing this parameter well can result in 2582 more efficient use of time and space. Run your code with -log_info 2583 to determine an optimal value to use. 2584 2585 Output Parameters: 2586 . fact - new matrix that has been symbolically factored 2587 2588 Notes: 2589 See the file ${PETSC_DIR}/Performace for additional information about 2590 choosing the fill factor for better efficiency. 2591 2592 .keywords: matrix, factor, incomplete, ILU, symbolic, fill 2593 2594 .seealso: MatLUFactorSymbolic(), MatLUFactorNumeric() 2595 @*/ 2596 int MatILUFactorSymbolic(Mat mat,IS row,IS col,double f,int fill,Mat *fact) 2597 { 2598 int ierr; 2599 2600 PetscFunctionBegin; 2601 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2602 PetscValidPointer(fact); 2603 if (fill < 0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,0,"Levels of fill negative %d",fill); 2604 if (f < 1.0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,0,"Expected fill less than 1.0 %g",f); 2605 if (!mat->ops->ilufactorsymbolic) SETERRQ(PETSC_ERR_SUP,0,"Only MatCreateMPIRowbs() matrices support parallel ILU"); 2606 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2607 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2608 2609 PLogEventBegin(MAT_ILUFactorSymbolic,mat,row,col,0); 2610 ierr = (*mat->ops->ilufactorsymbolic)(mat,row,col,f,fill,fact); CHKERRQ(ierr); 2611 PLogEventEnd(MAT_ILUFactorSymbolic,mat,row,col,0); 2612 PetscFunctionReturn(0); 2613 } 2614 2615 #undef __FUNC__ 2616 #define __FUNC__ "MatIncompleteCholeskyFactorSymbolic" 2617 /*@ 2618 MatIncompleteCholeskyFactorSymbolic - Performs symbolic incomplete 2619 Cholesky factorization for a symmetric matrix. Use 2620 MatCholeskyFactorNumeric() to complete the factorization. 2621 2622 Collective on Mat 2623 2624 Input Parameters: 2625 + mat - the matrix 2626 . perm - row and column permutation 2627 . fill - levels of fill 2628 - f - expected fill as ratio of original fill 2629 2630 Output Parameter: 2631 . fact - the factored matrix 2632 2633 Note: Currently only no-fill factorization is supported. 2634 2635 .keywords: matrix, factor, incomplete, ICC, Cholesky, symbolic, fill 2636 2637 .seealso: MatCholeskyFactorNumeric(), MatCholeskyFactor() 2638 @*/ 2639 int MatIncompleteCholeskyFactorSymbolic(Mat mat,IS perm,double f,int fill,Mat *fact) 2640 { 2641 int ierr; 2642 2643 PetscFunctionBegin; 2644 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2645 PetscValidPointer(fact); 2646 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2647 if (fill < 0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,0,"Fill negative %d",fill); 2648 if (f < 1.0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,0,"Expected fill less than 1.0 %g",f); 2649 if (!mat->ops->incompletecholeskyfactorsymbolic) SETERRQ(PETSC_ERR_SUP,0,"Currently only MatCreateMPIRowbs() matrices support ICC in parallel"); 2650 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2651 2652 PLogEventBegin(MAT_IncompleteCholeskyFactorSymbolic,mat,perm,0,0); 2653 ierr = (*mat->ops->incompletecholeskyfactorsymbolic)(mat,perm,f,fill,fact);CHKERRQ(ierr); 2654 PLogEventEnd(MAT_IncompleteCholeskyFactorSymbolic,mat,perm,0,0); 2655 PetscFunctionReturn(0); 2656 } 2657 2658 #undef __FUNC__ 2659 #define __FUNC__ "MatGetArray" 2660 /*@C 2661 MatGetArray - Returns a pointer to the element values in the matrix. 2662 The result of this routine is dependent on the underlying matrix data 2663 structure, and may not even work for certain matrix types. You MUST 2664 call MatRestoreArray() when you no longer need to access the array. 2665 2666 Not Collective 2667 2668 Input Parameter: 2669 . mat - the matrix 2670 2671 Output Parameter: 2672 . v - the location of the values 2673 2674 Currently returns an array only for the dense formats, giving access to 2675 the local portion of the matrix in the usual Fortran column-oriented format. 2676 2677 Fortran Note: 2678 This routine is used differently from Fortran 2679 $ Mat mat 2680 $ Scalar mat_array(1) 2681 $ PetscOffset i_mat 2682 $ int ierr 2683 $ call MatGetArray(mat,mat_array,i_mat,ierr) 2684 $ 2685 $ Access first local entry in matrix with 2686 $ value = mat_array(i_mat + 1) 2687 $ (note here array is treated as one dimensional) 2688 $ ...... other code 2689 $ call MatRestoreArray(mat,mat_array,i_mat,ierr) 2690 2691 See the Fortran chapter of the users manual and 2692 petsc/src/mat/examples/tests for details. 2693 2694 .keywords: matrix, array, elements, values 2695 2696 .seealso: MatRestoreArray() 2697 @*/ 2698 int MatGetArray(Mat mat,Scalar **v) 2699 { 2700 int ierr; 2701 2702 PetscFunctionBegin; 2703 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2704 PetscValidPointer(v); 2705 if (!mat->ops->getarray) SETERRQ(PETSC_ERR_SUP,0,""); 2706 ierr = (*mat->ops->getarray)(mat,v);CHKERRQ(ierr); 2707 PetscFunctionReturn(0); 2708 } 2709 2710 #undef __FUNC__ 2711 #define __FUNC__ "MatRestoreArray" 2712 /*@C 2713 MatRestoreArray - Restores the matrix after MatGetArray() has been called. 2714 2715 Not Collective 2716 2717 Input Parameter: 2718 + mat - the matrix 2719 - v - the location of the values 2720 2721 Fortran Note: 2722 This routine is used differently from Fortran 2723 $ Mat mat 2724 $ Scalar mat_array(1) 2725 $ PetscOffset i_mat 2726 $ int ierr 2727 $ call MatGetArray(mat,mat_array,i_mat,ierr) 2728 $ 2729 $ Access first local entry in matrix with 2730 $ value = mat_array(i_mat + 1) 2731 $ (note here the array is treated as one dimensional) 2732 $ 2733 $ ...... other code 2734 $ call MatRestoreArray(mat,mat_array,i_mat,ierr) 2735 2736 See the Fortran chapter of the users manual and 2737 petsc/src/mat/examples/tests for details 2738 2739 .keywords: matrix, array, elements, values, restore 2740 2741 .seealso: MatGetArray() 2742 @*/ 2743 int MatRestoreArray(Mat mat,Scalar **v) 2744 { 2745 int ierr; 2746 2747 PetscFunctionBegin; 2748 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2749 PetscValidPointer(v); 2750 if (!mat->ops->restorearray) SETERRQ(PETSC_ERR_SUP,0,""); 2751 ierr = (*mat->ops->restorearray)(mat,v);CHKERRQ(ierr); 2752 PetscFunctionReturn(0); 2753 } 2754 2755 #undef __FUNC__ 2756 #define __FUNC__ "MatGetSubMatrices" 2757 /*@C 2758 MatGetSubMatrices - Extracts several submatrices from a matrix. If submat 2759 points to an array of valid matrices, they may be reused to store the new 2760 submatrices. 2761 2762 Collective on Mat 2763 2764 Input Parameters: 2765 + mat - the matrix 2766 . n - the number of submatrixes to be extracted 2767 . irow, icol - index sets of rows and columns to extract 2768 - scall - either MAT_INITIAL_MATRIX or MAT_REUSE_MATRIX 2769 2770 Output Parameter: 2771 . submat - the array of submatrices 2772 2773 Notes: 2774 MatGetSubMatrices() can extract only sequential submatrices 2775 (from both sequential and parallel matrices). Use MatGetSubMatrix() 2776 to extract a parallel submatrix. 2777 2778 When extracting submatrices from a parallel matrix, each processor can 2779 form a different submatrix by setting the rows and columns of its 2780 individual index sets according to the local submatrix desired. 2781 2782 When finished using the submatrices, the user should destroy 2783 them with MatDestroySubMatrices(). 2784 2785 MAT_REUSE_MATRIX can only be used when the nonzero structure of the 2786 original matrix has not changed from that last call to MatGetSubMatrices() 2787 2788 Fortran Note: 2789 The Fortran interface is slightly different from that given below, it 2790 requires one to pass in as submat a Mat (integer) array of size at least m. 2791 2792 .keywords: matrix, get, submatrix, submatrices 2793 2794 .seealso: MatDestroyMatrices(), MatGetSubMatrix() 2795 @*/ 2796 int MatGetSubMatrices(Mat mat,int n,IS *irow,IS *icol,MatReuse scall,Mat **submat) 2797 { 2798 int ierr; 2799 2800 PetscFunctionBegin; 2801 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2802 if (!mat->ops->getsubmatrices) SETERRQ(PETSC_ERR_SUP,0,""); 2803 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2804 2805 PLogEventBegin(MAT_GetSubMatrices,mat,0,0,0); 2806 ierr = (*mat->ops->getsubmatrices)(mat,n,irow,icol,scall,submat); CHKERRQ(ierr); 2807 PLogEventEnd(MAT_GetSubMatrices,mat,0,0,0); 2808 2809 PetscFunctionReturn(0); 2810 } 2811 2812 #undef __FUNC__ 2813 #define __FUNC__ "MatDestroyMatrices" 2814 /*@C 2815 MatDestroyMatrices - Destroys a set of matrices obtained with MatGetSubMatrices(). 2816 2817 Collective on Mat 2818 2819 Input Parameters: 2820 + n - the number of local matrices 2821 - mat - the matrices 2822 2823 .keywords: matrix, destroy, submatrix, submatrices 2824 2825 .seealso: MatGetSubMatrices() 2826 @*/ 2827 int MatDestroyMatrices(int n,Mat **mat) 2828 { 2829 int ierr,i; 2830 2831 PetscFunctionBegin; 2832 if (n < 0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,1,"Trying to destroy negative number of matrices %d",n); 2833 PetscValidPointer(mat); 2834 for ( i=0; i<n; i++ ) { 2835 ierr = MatDestroy((*mat)[i]); CHKERRQ(ierr); 2836 } 2837 if (n) PetscFree(*mat); 2838 PetscFunctionReturn(0); 2839 } 2840 2841 #undef __FUNC__ 2842 #define __FUNC__ "MatIncreaseOverlap" 2843 /*@ 2844 MatIncreaseOverlap - Given a set of submatrices indicated by index sets, 2845 replaces the index sets by larger ones that represent submatrices with 2846 additional overlap. 2847 2848 Collective on Mat 2849 2850 Input Parameters: 2851 + mat - the matrix 2852 . n - the number of index sets 2853 . is - the array of pointers to index sets 2854 - ov - the additional overlap requested 2855 2856 .keywords: matrix, overlap, Schwarz 2857 2858 .seealso: MatGetSubMatrices() 2859 @*/ 2860 int MatIncreaseOverlap(Mat mat,int n, IS *is,int ov) 2861 { 2862 int ierr; 2863 2864 PetscFunctionBegin; 2865 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2866 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2867 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2868 2869 if (ov == 0) PetscFunctionReturn(0); 2870 if (!mat->ops->increaseoverlap) SETERRQ(PETSC_ERR_SUP,0,""); 2871 PLogEventBegin(MAT_IncreaseOverlap,mat,0,0,0); 2872 ierr = (*mat->ops->increaseoverlap)(mat,n,is,ov); CHKERRQ(ierr); 2873 PLogEventEnd(MAT_IncreaseOverlap,mat,0,0,0); 2874 PetscFunctionReturn(0); 2875 } 2876 2877 #undef __FUNC__ 2878 #define __FUNC__ "MatPrintHelp" 2879 /*@ 2880 MatPrintHelp - Prints all the options for the matrix. 2881 2882 Collective on Mat 2883 2884 Input Parameter: 2885 . mat - the matrix 2886 2887 Options Database Keys: 2888 + -help - Prints matrix options 2889 - -h - Prints matrix options 2890 2891 .keywords: mat, help 2892 2893 .seealso: MatCreate(), MatCreateXXX() 2894 @*/ 2895 int MatPrintHelp(Mat mat) 2896 { 2897 static int called = 0; 2898 int ierr; 2899 MPI_Comm comm; 2900 2901 PetscFunctionBegin; 2902 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2903 2904 comm = mat->comm; 2905 if (!called) { 2906 (*PetscHelpPrintf)(comm,"General matrix options:\n"); 2907 (*PetscHelpPrintf)(comm," -mat_view_info: view basic matrix info during MatAssemblyEnd()\n"); 2908 (*PetscHelpPrintf)(comm," -mat_view_info_detailed: view detailed matrix info during MatAssemblyEnd()\n"); 2909 (*PetscHelpPrintf)(comm," -mat_view_draw: draw nonzero matrix structure during MatAssemblyEnd()\n"); 2910 (*PetscHelpPrintf)(comm," -draw_pause <sec>: set seconds of display pause\n"); 2911 (*PetscHelpPrintf)(comm," -display <name>: set alternate display\n"); 2912 called = 1; 2913 } 2914 if (mat->ops->printhelp) { 2915 ierr = (*mat->ops->printhelp)(mat); CHKERRQ(ierr); 2916 } 2917 PetscFunctionReturn(0); 2918 } 2919 2920 #undef __FUNC__ 2921 #define __FUNC__ "MatGetBlockSize" 2922 /*@ 2923 MatGetBlockSize - Returns the matrix block size; useful especially for the 2924 block row and block diagonal formats. 2925 2926 Not Collective 2927 2928 Input Parameter: 2929 . mat - the matrix 2930 2931 Output Parameter: 2932 . bs - block size 2933 2934 Notes: 2935 Block diagonal formats are MATSEQBDIAG, MATMPIBDIAG. 2936 Block row formats are MATSEQBAIJ, MATMPIBAIJ 2937 2938 .keywords: matrix, get, block, size 2939 2940 .seealso: MatCreateSeqBAIJ(), MatCreateMPIBAIJ(), MatCreateSeqBDiag(), MatCreateMPIBDiag() 2941 @*/ 2942 int MatGetBlockSize(Mat mat,int *bs) 2943 { 2944 int ierr; 2945 2946 PetscFunctionBegin; 2947 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2948 PetscValidIntPointer(bs); 2949 if (!mat->ops->getblocksize) SETERRQ(PETSC_ERR_SUP,0,""); 2950 ierr = (*mat->ops->getblocksize)(mat,bs);CHKERRQ(ierr); 2951 PetscFunctionReturn(0); 2952 } 2953 2954 #undef __FUNC__ 2955 #define __FUNC__ "MatGetRowIJ" 2956 /*@C 2957 MatGetRowIJ - Returns the compressed row storage i and j indices for sequential matrices. 2958 EXPERTS ONLY. 2959 2960 Collective on Mat 2961 2962 Input Parameters: 2963 + mat - the matrix 2964 . shift - 0 or 1 indicating we want the indices starting at 0 or 1 2965 - symmetric - PETSC_TRUE or PETSC_FALSE indicating the matrix data structure should be 2966 symmetrized 2967 2968 Output Parameters: 2969 + n - number of rows in the (possibly compressed) matrix 2970 . ia - the row pointers 2971 . ja - the column indices 2972 - done - PETSC_TRUE or PETSC_FALSE, indicating whether the values have been returned 2973 2974 .seealso: MatGetColumnIJ(), MatRestoreRowIJ() 2975 @*/ 2976 int MatGetRowIJ(Mat mat,int shift,PetscTruth symmetric,int *n,int **ia,int** ja,PetscTruth *done) 2977 { 2978 int ierr; 2979 2980 PetscFunctionBegin; 2981 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2982 if (ia) PetscValidIntPointer(ia); 2983 if (ja) PetscValidIntPointer(ja); 2984 PetscValidIntPointer(done); 2985 if (!mat->ops->getrowij) *done = PETSC_FALSE; 2986 else { 2987 *done = PETSC_TRUE; 2988 ierr = (*mat->ops->getrowij)(mat,shift,symmetric,n,ia,ja,done); CHKERRQ(ierr); 2989 } 2990 PetscFunctionReturn(0); 2991 } 2992 2993 #undef __FUNC__ 2994 #define __FUNC__ "MatGetColumnIJ" 2995 /*@C 2996 MatGetColumnIJ - Returns the compressed column storage i and j indices for sequential matrices. 2997 EXPERTS ONLY. 2998 2999 Collective on Mat 3000 3001 Input Parameters: 3002 + mat - the matrix 3003 . shift - 1 or zero indicating we want the indices starting at 0 or 1 3004 - symmetric - PETSC_TRUE or PETSC_FALSE indicating the matrix data structure should be 3005 symmetrized 3006 3007 Output Parameters: 3008 + n - number of columns in the (possibly compressed) matrix 3009 . ia - the column pointers 3010 . ja - the row indices 3011 - done - PETSC_TRUE or PETSC_FALSE, indicating whether the values have been returned 3012 3013 .seealso: MatGetRowIJ(), MatRestoreColumnIJ() 3014 @*/ 3015 int MatGetColumnIJ(Mat mat,int shift,PetscTruth symmetric,int *n,int **ia,int** ja,PetscTruth *done) 3016 { 3017 int ierr; 3018 3019 PetscFunctionBegin; 3020 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3021 if (ia) PetscValidIntPointer(ia); 3022 if (ja) PetscValidIntPointer(ja); 3023 PetscValidIntPointer(done); 3024 3025 if (!mat->ops->getcolumnij) *done = PETSC_FALSE; 3026 else { 3027 *done = PETSC_TRUE; 3028 ierr = (*mat->ops->getcolumnij)(mat,shift,symmetric,n,ia,ja,done); CHKERRQ(ierr); 3029 } 3030 PetscFunctionReturn(0); 3031 } 3032 3033 #undef __FUNC__ 3034 #define __FUNC__ "MatRestoreRowIJ" 3035 /*@C 3036 MatRestoreRowIJ - Call after you are completed with the ia,ja indices obtained with 3037 MatGetRowIJ(). EXPERTS ONLY. 3038 3039 Collective on Mat 3040 3041 Input Parameters: 3042 + mat - the matrix 3043 . shift - 1 or zero indicating we want the indices starting at 0 or 1 3044 - symmetric - PETSC_TRUE or PETSC_FALSE indicating the matrix data structure should be 3045 symmetrized 3046 3047 Output Parameters: 3048 + n - size of (possibly compressed) matrix 3049 . ia - the row pointers 3050 . ja - the column indices 3051 - done - PETSC_TRUE or PETSC_FALSE indicated that the values have been returned 3052 3053 .seealso: MatGetRowIJ(), MatRestoreColumnIJ() 3054 @*/ 3055 int MatRestoreRowIJ(Mat mat,int shift,PetscTruth symmetric,int *n,int **ia,int** ja,PetscTruth *done) 3056 { 3057 int ierr; 3058 3059 PetscFunctionBegin; 3060 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3061 if (ia) PetscValidIntPointer(ia); 3062 if (ja) PetscValidIntPointer(ja); 3063 PetscValidIntPointer(done); 3064 3065 if (!mat->ops->restorerowij) *done = PETSC_FALSE; 3066 else { 3067 *done = PETSC_TRUE; 3068 ierr = (*mat->ops->restorerowij)(mat,shift,symmetric,n,ia,ja,done); CHKERRQ(ierr); 3069 } 3070 PetscFunctionReturn(0); 3071 } 3072 3073 #undef __FUNC__ 3074 #define __FUNC__ "MatRestoreColumnIJ" 3075 /*@C 3076 MatRestoreColumnIJ - Call after you are completed with the ia,ja indices obtained with 3077 MatGetColumnIJ(). EXPERTS ONLY. 3078 3079 Collective on Mat 3080 3081 Input Parameters: 3082 + mat - the matrix 3083 . shift - 1 or zero indicating we want the indices starting at 0 or 1 3084 - symmetric - PETSC_TRUE or PETSC_FALSE indicating the matrix data structure should be 3085 symmetrized 3086 3087 Output Parameters: 3088 + n - size of (possibly compressed) matrix 3089 . ia - the column pointers 3090 . ja - the row indices 3091 - done - PETSC_TRUE or PETSC_FALSE indicated that the values have been returned 3092 3093 .seealso: MatGetColumnIJ(), MatRestoreRowIJ() 3094 @*/ 3095 int MatRestoreColumnIJ(Mat mat,int shift,PetscTruth symmetric,int *n,int **ia,int** ja,PetscTruth *done) 3096 { 3097 int ierr; 3098 3099 PetscFunctionBegin; 3100 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3101 if (ia) PetscValidIntPointer(ia); 3102 if (ja) PetscValidIntPointer(ja); 3103 PetscValidIntPointer(done); 3104 3105 if (!mat->ops->restorecolumnij) *done = PETSC_FALSE; 3106 else { 3107 *done = PETSC_TRUE; 3108 ierr = (*mat->ops->restorecolumnij)(mat,shift,symmetric,n,ia,ja,done); CHKERRQ(ierr); 3109 } 3110 PetscFunctionReturn(0); 3111 } 3112 3113 #undef __FUNC__ 3114 #define __FUNC__ "MatColoringPatch" 3115 /*@C 3116 MatColoringPatch - EXPERTS ONLY, used inside matrix coloring routines that 3117 use MatGetRowIJ() and/or MatGetColumnIJ(). 3118 3119 Collective on Mat 3120 3121 Input Parameters: 3122 + mat - the matrix 3123 . n - number of colors 3124 - colorarray - array indicating color for each column 3125 3126 Output Parameters: 3127 . iscoloring - coloring generated using colorarray information 3128 3129 .seealso: MatGetRowIJ(), MatGetColumnIJ() 3130 3131 .keywords: mat, coloring, patch 3132 @*/ 3133 int MatColoringPatch(Mat mat,int n,int *colorarray,ISColoring *iscoloring) 3134 { 3135 int ierr; 3136 3137 PetscFunctionBegin; 3138 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3139 PetscValidIntPointer(colorarray); 3140 3141 if (!mat->ops->coloringpatch) {SETERRQ(PETSC_ERR_SUP,0,"");} 3142 else { 3143 ierr = (*mat->ops->coloringpatch)(mat,n,colorarray,iscoloring); CHKERRQ(ierr); 3144 } 3145 PetscFunctionReturn(0); 3146 } 3147 3148 3149 #undef __FUNC__ 3150 #define __FUNC__ "MatSetUnfactored" 3151 /*@ 3152 MatSetUnfactored - Resets a factored matrix to be treated as unfactored. 3153 3154 Collective on Mat 3155 3156 Input Parameter: 3157 . mat - the factored matrix to be reset 3158 3159 Notes: 3160 This routine should be used only with factored matrices formed by in-place 3161 factorization via ILU(0) (or by in-place LU factorization for the MATSEQDENSE 3162 format). This option can save memory, for example, when solving nonlinear 3163 systems with a matrix-free Newton-Krylov method and a matrix-based, in-place 3164 ILU(0) preconditioner. 3165 3166 Note that one can specify in-place ILU(0) factorization by calling 3167 $ PCType(pc,PCILU); 3168 $ PCILUSeUseInPlace(pc); 3169 or by using the options -pc_type ilu -pc_ilu_in_place 3170 3171 In-place factorization ILU(0) can also be used as a local 3172 solver for the blocks within the block Jacobi or additive Schwarz 3173 methods (runtime option: -sub_pc_ilu_in_place). See the discussion 3174 of these preconditioners in the users manual for details on setting 3175 local solver options. 3176 3177 .seealso: PCILUSetUseInPlace(), PCLUSetUseInPlace() 3178 3179 .keywords: matrix-free, in-place ILU, in-place LU 3180 @*/ 3181 int MatSetUnfactored(Mat mat) 3182 { 3183 int ierr; 3184 3185 PetscFunctionBegin; 3186 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3187 mat->factor = 0; 3188 if (!mat->ops->setunfactored) PetscFunctionReturn(0); 3189 ierr = (*mat->ops->setunfactored)(mat); CHKERRQ(ierr); 3190 PetscFunctionReturn(0); 3191 } 3192 3193 #undef __FUNC__ 3194 #define __FUNC__ "MatGetType" 3195 /*@C 3196 MatGetType - Gets the matrix type and name (as a string) from the matrix. 3197 3198 Not Collective 3199 3200 Input Parameter: 3201 . mat - the matrix 3202 3203 Output Parameter: 3204 + type - the matrix type (or use PETSC_NULL) 3205 - name - name of matrix type (or use PETSC_NULL) 3206 3207 .keywords: matrix, get, type, name 3208 @*/ 3209 int MatGetType(Mat mat,MatType *type,char **name) 3210 { 3211 int itype = (int)mat->type; 3212 char *matname[10]; 3213 3214 PetscFunctionBegin; 3215 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3216 3217 if (type) *type = (MatType) mat->type; 3218 if (name) { 3219 /* Note: Be sure that this list corresponds to the enum in mat.h */ 3220 matname[0] = "MATSEQDENSE"; 3221 matname[1] = "MATSEQAIJ"; 3222 matname[2] = "MATMPIAIJ"; 3223 matname[3] = "MATSHELL"; 3224 matname[4] = "MATMPIROWBS"; 3225 matname[5] = "MATSEQBDIAG"; 3226 matname[6] = "MATMPIBDIAG"; 3227 matname[7] = "MATMPIDENSE"; 3228 matname[8] = "MATSEQBAIJ"; 3229 matname[9] = "MATMPIBAIJ"; 3230 3231 if (itype < 0 || itype > 9) *name = "Unknown matrix type"; 3232 else *name = matname[itype]; 3233 } 3234 PetscFunctionReturn(0); 3235 } 3236 3237 /*MC 3238 MatGetArrayF90 - Accesses a matrix array from Fortran90. 3239 3240 Synopsis: 3241 MatGetArrayF90(Mat x,{Scalar, pointer :: xx_v(:)},integer ierr) 3242 3243 Not collective 3244 3245 Input Parameter: 3246 . x - matrix 3247 3248 Output Parameters: 3249 + xx_v - the Fortran90 pointer to the array 3250 - ierr - error code 3251 3252 Example of Usage: 3253 .vb 3254 Scalar, pointer xx_v(:) 3255 .... 3256 call MatGetArrayF90(x,xx_v,ierr) 3257 a = xx_v(3) 3258 call MatRestoreArrayF90(x,xx_v,ierr) 3259 .ve 3260 3261 Notes: 3262 Not yet supported for all F90 compilers 3263 3264 .seealso: MatRestoreArrayF90(), MatGetArray(), MatRestoreArray() 3265 3266 .keywords: matrix, array, f90 3267 M*/ 3268 3269 /*MC 3270 MatRestoreArrayF90 - Restores a matrix array that has been 3271 accessed with MatGetArrayF90(). 3272 3273 Synopsis: 3274 MatRestoreArrayF90(Mat x,{Scalar, pointer :: xx_v(:)},integer ierr) 3275 3276 Not collective 3277 3278 Input Parameters: 3279 + x - matrix 3280 - xx_v - the Fortran90 pointer to the array 3281 3282 Output Parameter: 3283 . ierr - error code 3284 3285 Example of Usage: 3286 .vb 3287 Scalar, pointer xx_v(:) 3288 .... 3289 call MatGetArrayF90(x,xx_v,ierr) 3290 a = xx_v(3) 3291 call MatRestoreArrayF90(x,xx_v,ierr) 3292 .ve 3293 3294 Notes: 3295 Not yet supported for all F90 compilers 3296 3297 .seealso: MatGetArrayF90(), MatGetArray(), MatRestoreArray() 3298 3299 .keywords: matrix, array, f90 3300 M*/ 3301 3302 3303 #undef __FUNC__ 3304 #define __FUNC__ "MatGetSubMatrix" 3305 /*@ 3306 MatGetSubMatrix - Gets a single submatrix on the same number of processors 3307 as the original matrix. 3308 3309 Collective on Mat 3310 3311 Input Parameters: 3312 + mat - the original matrix 3313 . isrow - rows this processor should obtain 3314 . iscol - columns for all processors you wish to keep 3315 . csize - number of columns "local" to this processor (does nothing for sequential 3316 matrices). This should match the result from VecGetLocalSize(x,...) if you 3317 plan to use the matrix in a A*x 3318 - cll - either MAT_INITIAL_MATRIX or MAT_REUSE_MATRIX 3319 3320 Output Parameter: 3321 . newmat - the new submatrix, of the same type as the old 3322 3323 .seealso: MatGetSubMatrices() 3324 3325 @*/ 3326 int MatGetSubMatrix(Mat mat,IS isrow,IS iscol,int csize,MatReuse cll,Mat *newmat) 3327 { 3328 int ierr, size; 3329 Mat *local; 3330 3331 PetscFunctionBegin; 3332 MPI_Comm_size(mat->comm,&size); 3333 3334 /* if original matrix is on just one processor then use submatrix generated */ 3335 if (!mat->ops->getsubmatrix && size == 1 && cll == MAT_REUSE_MATRIX) { 3336 ierr = MatGetSubMatrices(mat,1,&isrow,&iscol,MAT_REUSE_MATRIX,&newmat);CHKERRQ(ierr); 3337 PetscFunctionReturn(0); 3338 } else if (!mat->ops->getsubmatrix && size == 1) { 3339 ierr = MatGetSubMatrices(mat,1,&isrow,&iscol,MAT_INITIAL_MATRIX,&local);CHKERRQ(ierr); 3340 *newmat = *local; 3341 PetscFree(local); 3342 PetscFunctionReturn(0); 3343 } 3344 3345 if (!mat->ops->getsubmatrix) SETERRQ(PETSC_ERR_SUP,0,"Not currently implemented"); 3346 ierr = (*mat->ops->getsubmatrix)(mat,isrow,iscol,csize,cll,newmat);CHKERRQ(ierr); 3347 PetscFunctionReturn(0); 3348 } 3349 3350 #undef __FUNC__ 3351 #define __FUNC__ "MatGetMaps" 3352 /*@C 3353 MatGetMaps - Returns the maps associated with the matrix. 3354 3355 Not Collective 3356 3357 Input Parameter: 3358 . mat - the matrix 3359 3360 Output Parameters: 3361 + rmap - the row (right) map 3362 - cmap - the column (left) map 3363 3364 .keywords: matrix, get, map 3365 @*/ 3366 int MatGetMaps(Mat mat,Map *rmap,Map *cmap) 3367 { 3368 int ierr; 3369 3370 PetscFunctionBegin; 3371 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3372 ierr = (*mat->ops->getmaps)(mat,rmap,cmap);CHKERRQ(ierr); 3373 PetscFunctionReturn(0); 3374 } 3375 3376 /* 3377 Version that works for all PETSc matrices 3378 */ 3379 #undef __FUNC__ 3380 #define __FUNC__ "MatGetMaps_Petsc" 3381 int MatGetMaps_Petsc(Mat mat,Map *rmap,Map *cmap) 3382 { 3383 PetscFunctionBegin; 3384 if (rmap) *rmap = mat->rmap; 3385 if (cmap) *cmap = mat->cmap; 3386 PetscFunctionReturn(0); 3387 } 3388