1 #ifdef PETSC_RCS_HEADER 2 static char vcid[] = "$Id: matrix.c,v 1.313 1999/01/13 21:40:15 bsmith Exp bsmith $"; 3 #endif 4 5 /* 6 This is where the abstract matrix operations are defined 7 */ 8 9 #include "src/mat/matimpl.h" /*I "mat.h" I*/ 10 #include "src/vec/vecimpl.h" 11 12 #undef __FUNC__ 13 #define __FUNC__ "MatGetRow" 14 /*@C 15 MatGetRow - Gets a row of a matrix. You MUST call MatRestoreRow() 16 for each row that you get to ensure that your application does 17 not bleed memory. 18 19 Not Collective 20 21 Input Parameters: 22 + mat - the matrix 23 - row - the row to get 24 25 Output Parameters: 26 + ncols - the number of nonzeros in the row 27 . cols - if nonzero, the column numbers 28 - vals - if nonzero, the values 29 30 Notes: 31 This routine is provided for people who need to have direct access 32 to the structure of a matrix. We hope that we provide enough 33 high-level matrix routines that few users will need it. 34 35 MatGetRow() always returns 0-based column indices, regardless of 36 whether the internal representation is 0-based (default) or 1-based. 37 38 For better efficiency, set cols and/or vals to PETSC_NULL if you do 39 not wish to extract these quantities. 40 41 The user can only examine the values extracted with MatGetRow(); 42 the values cannot be altered. To change the matrix entries, one 43 must use MatSetValues(). 44 45 You can only have one call to MatGetRow() outstanding for a particular 46 matrix at a time. 47 48 Fortran Notes: 49 The calling sequence from Fortran is 50 .vb 51 MatGetRow(matrix,row,ncols,cols,values,ierr) 52 Mat matrix (input) 53 integer row (input) 54 integer ncols (output) 55 integer cols(maxcols) (output) 56 double precision (or double complex) values(maxcols) output 57 .ve 58 where maxcols >= maximum nonzeros in any row of the matrix. 59 60 Caution: 61 Do not try to change the contents of the output arrays (cols and vals). 62 In some cases, this may corrupt the matrix. 63 64 .keywords: matrix, row, get, extract 65 66 .seealso: MatRestoreRow(), MatSetValues() 67 @*/ 68 int MatGetRow(Mat mat,int row,int *ncols,int **cols,Scalar **vals) 69 { 70 int ierr; 71 72 PetscFunctionBegin; 73 PetscValidHeaderSpecific(mat,MAT_COOKIE); 74 PetscValidIntPointer(ncols); 75 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 76 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 77 if (!mat->ops->getrow) SETERRQ(PETSC_ERR_SUP,0,""); 78 PLogEventBegin(MAT_GetRow,mat,0,0,0); 79 ierr = (*mat->ops->getrow)(mat,row,ncols,cols,vals); CHKERRQ(ierr); 80 PLogEventEnd(MAT_GetRow,mat,0,0,0); 81 PetscFunctionReturn(0); 82 } 83 84 #undef __FUNC__ 85 #define __FUNC__ "MatRestoreRow" 86 /*@C 87 MatRestoreRow - Frees any temporary space allocated by MatGetRow(). 88 89 Not Collective 90 91 Input Parameters: 92 + mat - the matrix 93 . row - the row to get 94 . ncols, cols - the number of nonzeros and their columns 95 - vals - if nonzero the column values 96 97 Notes: 98 This routine should be called after you have finished examining the entries. 99 100 Fortran Notes: 101 The calling sequence from Fortran is 102 .vb 103 MatRestoreRow(matrix,row,ncols,cols,values,ierr) 104 Mat matrix (input) 105 integer row (input) 106 integer ncols (output) 107 integer cols(maxcols) (output) 108 double precision (or double complex) values(maxcols) output 109 .ve 110 Where maxcols >= maximum nonzeros in any row of the matrix. 111 112 In Fortran MatRestoreRow() MUST be called after MatGetRow() 113 before another call to MatGetRow() can be made. 114 115 .keywords: matrix, row, restore 116 117 .seealso: MatGetRow() 118 @*/ 119 int MatRestoreRow(Mat mat,int row,int *ncols,int **cols,Scalar **vals) 120 { 121 int ierr; 122 123 PetscFunctionBegin; 124 PetscValidHeaderSpecific(mat,MAT_COOKIE); 125 PetscValidIntPointer(ncols); 126 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 127 if (!mat->ops->restorerow) PetscFunctionReturn(0); 128 ierr = (*mat->ops->restorerow)(mat,row,ncols,cols,vals);CHKERRQ(ierr); 129 PetscFunctionReturn(0); 130 } 131 132 #undef __FUNC__ 133 #define __FUNC__ "MatView" 134 /*@C 135 MatView - Visualizes a matrix object. 136 137 Collective on Mat unless Viewer is VIEWER_STDOUT_SELF 138 139 Input Parameters: 140 + mat - the matrix 141 - ptr - visualization context 142 143 Notes: 144 The available visualization contexts include 145 + VIEWER_STDOUT_SELF - standard output (default) 146 . VIEWER_STDOUT_WORLD - synchronized standard 147 output where only the first processor opens 148 the file. All other processors send their 149 data to the first processor to print. 150 - VIEWER_DRAW_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_DRAW_(mat->comm),VIEWER_FORMAT_DRAW_CONTOUR,0);CHKERRQ(ierr); 1137 } 1138 ierr = MatView(*fact,VIEWER_DRAW_(mat->comm)); CHKERRQ(ierr); 1139 ierr = ViewerFlush(VIEWER_DRAW_(mat->comm)); CHKERRQ(ierr); 1140 if (flg) { 1141 ViewerPopFormat(VIEWER_DRAW_(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_DRAW_(mat->comm),VIEWER_FORMAT_DRAW_CONTOUR,0);CHKERRQ(ierr); 2152 } 2153 ierr = MatView(mat,VIEWER_DRAW_(mat->comm)); CHKERRQ(ierr); 2154 ierr = ViewerFlush(VIEWER_DRAW_(mat->comm)); CHKERRQ(ierr); 2155 if (flg) { 2156 ViewerPopFormat(VIEWER_DRAW_(mat->comm));CHKERRQ(ierr); 2157 } 2158 } 2159 ierr = OptionsHasName(PETSC_NULL,"-mat_view_socket",&flg); CHKERRQ(ierr); 2160 if (flg) { 2161 ierr = MatView(mat,VIEWER_SOCKET_(mat->comm)); CHKERRQ(ierr); 2162 ierr = ViewerFlush(VIEWER_SOCKET_(mat->comm)); CHKERRQ(ierr); 2163 } 2164 PetscFunctionReturn(0); 2165 } 2166 2167 #undef __FUNC__ 2168 #define __FUNC__ "MatAssemblyEnd" 2169 /*@ 2170 MatAssemblyEnd - Completes assembling the matrix. This routine should 2171 be called after MatAssemblyBegin(). 2172 2173 Collective on Mat 2174 2175 Input Parameters: 2176 + mat - the matrix 2177 - type - type of assembly, either MAT_FLUSH_ASSEMBLY or MAT_FINAL_ASSEMBLY 2178 2179 Options Database Keys: 2180 + -mat_view_info - Prints info on matrix at conclusion of MatEndAssembly() 2181 . -mat_view_info_detailed - Prints more detailed info 2182 . -mat_view - Prints matrix in ASCII format 2183 . -mat_view_matlab - Prints matrix in Matlab format 2184 . -mat_view_draw - Draws nonzero structure of matrix, using MatView() and DrawOpenX(). 2185 . -display <name> - Sets display name (default is host) 2186 - -draw_pause <sec> - Sets number of seconds to pause after display 2187 2188 Notes: 2189 MatSetValues() generally caches the values. The matrix is ready to 2190 use only after MatAssemblyBegin() and MatAssemblyEnd() have been called. 2191 Use MAT_FLUSH_ASSEMBLY when switching between ADD_VALUES and INSERT_VALUES 2192 in MatSetValues(); use MAT_FINAL_ASSEMBLY for the final assembly before 2193 using the matrix. 2194 2195 .keywords: matrix, assembly, assemble, end 2196 2197 .seealso: MatAssemblyBegin(), MatSetValues(), DrawOpenX(), MatView() 2198 @*/ 2199 int MatAssemblyEnd(Mat mat,MatAssemblyType type) 2200 { 2201 int ierr; 2202 static int inassm = 0; 2203 2204 PetscFunctionBegin; 2205 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2206 2207 inassm++; 2208 MatAssemblyEnd_InUse++; 2209 PLogEventBegin(MAT_AssemblyEnd,mat,0,0,0); 2210 if (mat->ops->assemblyend) { 2211 ierr = (*mat->ops->assemblyend)(mat,type); CHKERRQ(ierr); 2212 } 2213 2214 /* Flush assembly is not a true assembly */ 2215 if (type != MAT_FLUSH_ASSEMBLY) { 2216 mat->assembled = PETSC_TRUE; mat->num_ass++; 2217 } 2218 mat->insertmode = NOT_SET_VALUES; 2219 PLogEventEnd(MAT_AssemblyEnd,mat,0,0,0); 2220 MatAssemblyEnd_InUse--; 2221 2222 if (inassm == 1 && type != MAT_FLUSH_ASSEMBLY) { 2223 ierr = MatView_Private(mat); CHKERRQ(ierr); 2224 } 2225 inassm--; 2226 PetscFunctionReturn(0); 2227 } 2228 2229 2230 #undef __FUNC__ 2231 #define __FUNC__ "MatCompress" 2232 /*@ 2233 MatCompress - Tries to store the matrix in as little space as 2234 possible. May fail if memory is already fully used, since it 2235 tries to allocate new space. 2236 2237 Collective on Mat 2238 2239 Input Parameters: 2240 . mat - the matrix 2241 2242 .keywords: matrix, compress 2243 @*/ 2244 int MatCompress(Mat mat) 2245 { 2246 int ierr; 2247 2248 PetscFunctionBegin; 2249 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2250 if (mat->ops->compress) {ierr = (*mat->ops->compress)(mat);CHKERRQ(ierr);} 2251 PetscFunctionReturn(0); 2252 } 2253 2254 #undef __FUNC__ 2255 #define __FUNC__ "MatSetOption" 2256 /*@ 2257 MatSetOption - Sets a parameter option for a matrix. Some options 2258 may be specific to certain storage formats. Some options 2259 determine how values will be inserted (or added). Sorted, 2260 row-oriented input will generally assemble the fastest. The default 2261 is row-oriented, nonsorted input. 2262 2263 Collective on Mat 2264 2265 Input Parameters: 2266 + mat - the matrix 2267 - option - the option, one of those listed below (and possibly others), 2268 e.g., MAT_ROWS_SORTED, MAT_NEW_NONZERO_LOCATION_ERROR 2269 2270 Options Describing Matrix Structure: 2271 + MAT_SYMMETRIC - symmetric in terms of both structure and value 2272 - MAT_STRUCTURALLY_SYMMETRIC - symmetric nonzero structure 2273 2274 Options For Use with MatSetValues(): 2275 Insert a logically dense subblock, which can be 2276 + MAT_ROW_ORIENTED - row-oriented 2277 . MAT_COLUMN_ORIENTED - column-oriented 2278 . MAT_ROWS_SORTED - sorted by row 2279 . MAT_ROWS_UNSORTED - not sorted by row 2280 . MAT_COLUMNS_SORTED - sorted by column 2281 - MAT_COLUMNS_UNSORTED - not sorted by column 2282 2283 Not these options reflect the data you pass in with MatSetValues(); it has 2284 nothing to do with how the data is stored internally in the matrix 2285 data structure. 2286 2287 When (re)assembling a matrix, we can restrict the input for 2288 efficiency/debugging purposes. These options include 2289 + MAT_NO_NEW_NONZERO_LOCATIONS - additional insertions will not be 2290 allowed if they generate a new nonzero 2291 . MAT_YES_NEW_NONZERO_LOCATIONS - additional insertions will be allowed 2292 . MAT_NO_NEW_DIAGONALS - additional insertions will not be allowed if 2293 they generate a nonzero in a new diagonal (for block diagonal format only) 2294 . MAT_YES_NEW_DIAGONALS - new diagonals will be allowed (for block diagonal format only) 2295 . MAT_IGNORE_OFF_PROC_ENTRIES - drops off-processor entries 2296 . MAT_NEW_NONZERO_LOCATION_ERROR - generates an error for new matrix entry 2297 - MAT_USE_HASH_TABLE - uses a hash table to speed up matrix assembly 2298 2299 Notes: 2300 Some options are relevant only for particular matrix types and 2301 are thus ignored by others. Other options are not supported by 2302 certain matrix types and will generate an error message if set. 2303 2304 If using a Fortran 77 module to compute a matrix, one may need to 2305 use the column-oriented option (or convert to the row-oriented 2306 format). 2307 2308 MAT_NO_NEW_NONZERO_LOCATIONS indicates that any add or insertion 2309 that would generate a new entry in the nonzero structure is instead 2310 ignored. Thus, if memory has not alredy been allocated for this particular 2311 data, then the insertion is ignored. For dense matrices, in which 2312 the entire array is allocated, no entries are ever ignored. 2313 2314 MAT_NEW_NONZERO_LOCATION_ERROR indicates that any add or insertion 2315 that would generate a new entry in the nonzero structure instead produces 2316 an error. (Currently supported for AIJ and BAIJ formats only.) 2317 This is a useful flag when using SAME_NONZERO_PATTERN in calling 2318 SLESSetOperators() to ensure that the nonzero pattern truely does 2319 remain unchanged. 2320 2321 MAT_NEW_NONZERO_ALLOCATION_ERROR indicates that any add or insertion 2322 that would generate a new entry that has not been preallocated will 2323 instead produce an error. (Currently supported for AIJ and BAIJ formats 2324 only.) This is a useful flag when debugging matrix memory preallocation. 2325 2326 MAT_IGNORE_OFF_PROC_ENTRIES indicates entries destined for 2327 other processors should be dropped, rather than stashed. 2328 This is useful if you know that the "owning" processor is also 2329 always generating the correct matrix entries, so that PETSc need 2330 not transfer duplicate entries generated on another processor. 2331 2332 MAT_USE_HASH_TABLE indicates that a hash table be used to improve the 2333 searches during matrix assembly. When this flag is set, the hash table 2334 is created during the first Matrix Assembly. This hash table is 2335 used the next time through, during MatSetVaules()/MatSetVaulesBlocked() 2336 to improve the searching of indices. MAT_NO_NEW_NONZERO_LOCATIONS flag 2337 should be used with MAT_USE_HASH_TABLE flag. This option is currently 2338 supported by MATMPIBAIJ format only. 2339 2340 .keywords: matrix, option, row-oriented, column-oriented, sorted, nonzero 2341 @*/ 2342 int MatSetOption(Mat mat,MatOption op) 2343 { 2344 int ierr; 2345 2346 PetscFunctionBegin; 2347 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2348 if (mat->ops->setoption) {ierr = (*mat->ops->setoption)(mat,op);CHKERRQ(ierr);} 2349 PetscFunctionReturn(0); 2350 } 2351 2352 #undef __FUNC__ 2353 #define __FUNC__ "MatZeroEntries" 2354 /*@ 2355 MatZeroEntries - Zeros all entries of a matrix. For sparse matrices 2356 this routine retains the old nonzero structure. 2357 2358 Collective on Mat 2359 2360 Input Parameters: 2361 . mat - the matrix 2362 2363 .keywords: matrix, zero, entries 2364 2365 .seealso: MatZeroRows() 2366 @*/ 2367 int MatZeroEntries(Mat mat) 2368 { 2369 int ierr; 2370 2371 PetscFunctionBegin; 2372 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2373 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2374 if (!mat->ops->zeroentries) SETERRQ(PETSC_ERR_SUP,0,""); 2375 2376 PLogEventBegin(MAT_ZeroEntries,mat,0,0,0); 2377 ierr = (*mat->ops->zeroentries)(mat); CHKERRQ(ierr); 2378 PLogEventEnd(MAT_ZeroEntries,mat,0,0,0); 2379 PetscFunctionReturn(0); 2380 } 2381 2382 #undef __FUNC__ 2383 #define __FUNC__ "MatZeroRows" 2384 /*@ 2385 MatZeroRows - Zeros all entries (except possibly the main diagonal) 2386 of a set of rows of a matrix. 2387 2388 Collective on Mat 2389 2390 Input Parameters: 2391 + mat - the matrix 2392 . is - index set of rows to remove 2393 - diag - pointer to value put in all diagonals of eliminated rows. 2394 Note that diag is not a pointer to an array, but merely a 2395 pointer to a single value. 2396 2397 Notes: 2398 For the AIJ matrix formats this removes the old nonzero structure, 2399 but does not release memory. For the dense and block diagonal 2400 formats this does not alter the nonzero structure. 2401 2402 The user can set a value in the diagonal entry (or for the AIJ and 2403 row formats can optionally remove the main diagonal entry from the 2404 nonzero structure as well, by passing a null pointer as the final 2405 argument). 2406 2407 For the parallel case, all processes that share the matrix (i.e., 2408 those in the communicator used for matrix creation) MUST call this 2409 routine, regardless of whether any rows being zeroed are owned by 2410 them. 2411 2412 .keywords: matrix, zero, rows, boundary conditions 2413 2414 .seealso: MatZeroEntries(), 2415 @*/ 2416 int MatZeroRows(Mat mat,IS is, Scalar *diag) 2417 { 2418 int ierr; 2419 2420 PetscFunctionBegin; 2421 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2422 PetscValidHeaderSpecific(is,IS_COOKIE); 2423 if (diag) PetscValidScalarPointer(diag); 2424 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2425 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2426 if (!mat->ops->zerorows) SETERRQ(PETSC_ERR_SUP,0,""); 2427 2428 ierr = (*mat->ops->zerorows)(mat,is,diag); CHKERRQ(ierr); 2429 ierr = MatView_Private(mat); CHKERRQ(ierr); 2430 PetscFunctionReturn(0); 2431 } 2432 2433 #undef __FUNC__ 2434 #define __FUNC__ "MatZeroRowsLocal" 2435 /*@ 2436 MatZeroRowsLocal - Zeros all entries (except possibly the main diagonal) 2437 of a set of rows of a matrix; using local numbering of rows. 2438 2439 Collective on Mat 2440 2441 Input Parameters: 2442 + mat - the matrix 2443 . is - index set of rows to remove 2444 - diag - pointer to value put in all diagonals of eliminated rows. 2445 Note that diag is not a pointer to an array, but merely a 2446 pointer to a single value. 2447 2448 Notes: 2449 For the AIJ matrix formats this removes the old nonzero structure, 2450 but does not release memory. For the dense and block diagonal 2451 formats this does not alter the nonzero structure. 2452 2453 The user can set a value in the diagonal entry (or for the AIJ and 2454 row formats can optionally remove the main diagonal entry from the 2455 nonzero structure as well, by passing a null pointer as the final 2456 argument). 2457 2458 .keywords: matrix, zero, rows, boundary conditions 2459 2460 .seealso: MatZeroEntries(), 2461 @*/ 2462 int MatZeroRowsLocal(Mat mat,IS is, Scalar *diag) 2463 { 2464 int ierr; 2465 IS newis; 2466 2467 PetscFunctionBegin; 2468 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2469 PetscValidHeaderSpecific(is,IS_COOKIE); 2470 if (diag) PetscValidScalarPointer(diag); 2471 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2472 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2473 if (!mat->ops->zerorows) SETERRQ(PETSC_ERR_SUP,0,""); 2474 2475 ierr = ISLocalToGlobalMappingApplyIS(mat->mapping,is,&newis); CHKERRQ(ierr); 2476 ierr = (*mat->ops->zerorows)(mat,newis,diag); CHKERRQ(ierr); 2477 ierr = ISDestroy(newis); 2478 PetscFunctionReturn(0); 2479 } 2480 2481 #undef __FUNC__ 2482 #define __FUNC__ "MatGetSize" 2483 /*@ 2484 MatGetSize - Returns the numbers of rows and columns in a matrix. 2485 2486 Not Collective 2487 2488 Input Parameter: 2489 . mat - the matrix 2490 2491 Output Parameters: 2492 + m - the number of global rows 2493 - n - the number of global columns 2494 2495 .keywords: matrix, dimension, size, rows, columns, global, get 2496 2497 .seealso: MatGetLocalSize() 2498 @*/ 2499 int MatGetSize(Mat mat,int *m,int* n) 2500 { 2501 int ierr; 2502 2503 PetscFunctionBegin; 2504 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2505 ierr = (*mat->ops->getsize)(mat,m,n);CHKERRQ(ierr); 2506 PetscFunctionReturn(0); 2507 } 2508 2509 #undef __FUNC__ 2510 #define __FUNC__ "MatGetLocalSize" 2511 /*@ 2512 MatGetLocalSize - Returns the number of rows and columns in a matrix 2513 stored locally. This information may be implementation dependent, so 2514 use with care. 2515 2516 Not Collective 2517 2518 Input Parameters: 2519 . mat - the matrix 2520 2521 Output Parameters: 2522 + m - the number of local rows 2523 - n - the number of local columns 2524 2525 .keywords: matrix, dimension, size, local, rows, columns, get 2526 2527 .seealso: MatGetSize() 2528 @*/ 2529 int MatGetLocalSize(Mat mat,int *m,int* n) 2530 { 2531 int ierr; 2532 2533 PetscFunctionBegin; 2534 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2535 ierr = (*mat->ops->getlocalsize)(mat,m,n);CHKERRQ(ierr); 2536 PetscFunctionReturn(0); 2537 } 2538 2539 #undef __FUNC__ 2540 #define __FUNC__ "MatGetOwnershipRange" 2541 /*@ 2542 MatGetOwnershipRange - Returns the range of matrix rows owned by 2543 this processor, assuming that the matrix is laid out with the first 2544 n1 rows on the first processor, the next n2 rows on the second, etc. 2545 For certain parallel layouts this range may not be well defined. 2546 2547 Not Collective 2548 2549 Input Parameters: 2550 . mat - the matrix 2551 2552 Output Parameters: 2553 + m - the global index of the first local row 2554 - n - one more than the global index of the last local row 2555 2556 .keywords: matrix, get, range, ownership 2557 @*/ 2558 int MatGetOwnershipRange(Mat mat,int *m,int* n) 2559 { 2560 int ierr; 2561 2562 PetscFunctionBegin; 2563 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2564 PetscValidIntPointer(m); 2565 PetscValidIntPointer(n); 2566 if (!mat->ops->getownershiprange) SETERRQ(PETSC_ERR_SUP,0,""); 2567 ierr = (*mat->ops->getownershiprange)(mat,m,n);CHKERRQ(ierr); 2568 PetscFunctionReturn(0); 2569 } 2570 2571 #undef __FUNC__ 2572 #define __FUNC__ "MatILUFactorSymbolic" 2573 /*@ 2574 MatILUFactorSymbolic - Performs symbolic ILU factorization of a matrix. 2575 Uses levels of fill only, not drop tolerance. Use MatLUFactorNumeric() 2576 to complete the factorization. 2577 2578 Collective on Mat 2579 2580 Input Parameters: 2581 + mat - the matrix 2582 . row - row permutation 2583 . column - column permutation 2584 . fill - number of levels of fill 2585 - f - expected fill as ratio of the original number of nonzeros, 2586 for example 3.0; choosing this parameter well can result in 2587 more efficient use of time and space. Run your code with -log_info 2588 to determine an optimal value to use. 2589 2590 Output Parameters: 2591 . fact - new matrix that has been symbolically factored 2592 2593 Notes: 2594 See the file ${PETSC_DIR}/Performace for additional information about 2595 choosing the fill factor for better efficiency. 2596 2597 .keywords: matrix, factor, incomplete, ILU, symbolic, fill 2598 2599 .seealso: MatLUFactorSymbolic(), MatLUFactorNumeric() 2600 @*/ 2601 int MatILUFactorSymbolic(Mat mat,IS row,IS col,double f,int fill,Mat *fact) 2602 { 2603 int ierr; 2604 2605 PetscFunctionBegin; 2606 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2607 PetscValidPointer(fact); 2608 if (fill < 0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,0,"Levels of fill negative %d",fill); 2609 if (f < 1.0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,0,"Expected fill less than 1.0 %g",f); 2610 if (!mat->ops->ilufactorsymbolic) SETERRQ(PETSC_ERR_SUP,0,"Only MatCreateMPIRowbs() matrices support parallel ILU"); 2611 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2612 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2613 2614 PLogEventBegin(MAT_ILUFactorSymbolic,mat,row,col,0); 2615 ierr = (*mat->ops->ilufactorsymbolic)(mat,row,col,f,fill,fact); CHKERRQ(ierr); 2616 PLogEventEnd(MAT_ILUFactorSymbolic,mat,row,col,0); 2617 PetscFunctionReturn(0); 2618 } 2619 2620 #undef __FUNC__ 2621 #define __FUNC__ "MatIncompleteCholeskyFactorSymbolic" 2622 /*@ 2623 MatIncompleteCholeskyFactorSymbolic - Performs symbolic incomplete 2624 Cholesky factorization for a symmetric matrix. Use 2625 MatCholeskyFactorNumeric() to complete the factorization. 2626 2627 Collective on Mat 2628 2629 Input Parameters: 2630 + mat - the matrix 2631 . perm - row and column permutation 2632 . fill - levels of fill 2633 - f - expected fill as ratio of original fill 2634 2635 Output Parameter: 2636 . fact - the factored matrix 2637 2638 Note: Currently only no-fill factorization is supported. 2639 2640 .keywords: matrix, factor, incomplete, ICC, Cholesky, symbolic, fill 2641 2642 .seealso: MatCholeskyFactorNumeric(), MatCholeskyFactor() 2643 @*/ 2644 int MatIncompleteCholeskyFactorSymbolic(Mat mat,IS perm,double f,int fill,Mat *fact) 2645 { 2646 int ierr; 2647 2648 PetscFunctionBegin; 2649 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2650 PetscValidPointer(fact); 2651 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2652 if (fill < 0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,0,"Fill negative %d",fill); 2653 if (f < 1.0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,0,"Expected fill less than 1.0 %g",f); 2654 if (!mat->ops->incompletecholeskyfactorsymbolic) SETERRQ(PETSC_ERR_SUP,0,"Currently only MatCreateMPIRowbs() matrices support ICC in parallel"); 2655 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2656 2657 PLogEventBegin(MAT_IncompleteCholeskyFactorSymbolic,mat,perm,0,0); 2658 ierr = (*mat->ops->incompletecholeskyfactorsymbolic)(mat,perm,f,fill,fact);CHKERRQ(ierr); 2659 PLogEventEnd(MAT_IncompleteCholeskyFactorSymbolic,mat,perm,0,0); 2660 PetscFunctionReturn(0); 2661 } 2662 2663 #undef __FUNC__ 2664 #define __FUNC__ "MatGetArray" 2665 /*@C 2666 MatGetArray - Returns a pointer to the element values in the matrix. 2667 The result of this routine is dependent on the underlying matrix data 2668 structure, and may not even work for certain matrix types. You MUST 2669 call MatRestoreArray() when you no longer need to access the array. 2670 2671 Not Collective 2672 2673 Input Parameter: 2674 . mat - the matrix 2675 2676 Output Parameter: 2677 . v - the location of the values 2678 2679 Currently returns an array only for the dense formats, giving access to 2680 the local portion of the matrix in the usual Fortran column-oriented format. 2681 2682 Fortran Note: 2683 This routine is used differently from Fortran 2684 $ Mat mat 2685 $ Scalar mat_array(1) 2686 $ PetscOffset i_mat 2687 $ int ierr 2688 $ call MatGetArray(mat,mat_array,i_mat,ierr) 2689 $ 2690 $ Access first local entry in matrix with 2691 $ value = mat_array(i_mat + 1) 2692 $ (note here array is treated as one dimensional) 2693 $ ...... other code 2694 $ call MatRestoreArray(mat,mat_array,i_mat,ierr) 2695 2696 See the Fortran chapter of the users manual and 2697 petsc/src/mat/examples/tests for details. 2698 2699 .keywords: matrix, array, elements, values 2700 2701 .seealso: MatRestoreArray() 2702 @*/ 2703 int MatGetArray(Mat mat,Scalar **v) 2704 { 2705 int ierr; 2706 2707 PetscFunctionBegin; 2708 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2709 PetscValidPointer(v); 2710 if (!mat->ops->getarray) SETERRQ(PETSC_ERR_SUP,0,""); 2711 ierr = (*mat->ops->getarray)(mat,v);CHKERRQ(ierr); 2712 PetscFunctionReturn(0); 2713 } 2714 2715 #undef __FUNC__ 2716 #define __FUNC__ "MatRestoreArray" 2717 /*@C 2718 MatRestoreArray - Restores the matrix after MatGetArray() has been called. 2719 2720 Not Collective 2721 2722 Input Parameter: 2723 + mat - the matrix 2724 - v - the location of the values 2725 2726 Fortran Note: 2727 This routine is used differently from Fortran 2728 $ Mat mat 2729 $ Scalar mat_array(1) 2730 $ PetscOffset i_mat 2731 $ int ierr 2732 $ call MatGetArray(mat,mat_array,i_mat,ierr) 2733 $ 2734 $ Access first local entry in matrix with 2735 $ value = mat_array(i_mat + 1) 2736 $ (note here the array is treated as one dimensional) 2737 $ 2738 $ ...... other code 2739 $ call MatRestoreArray(mat,mat_array,i_mat,ierr) 2740 2741 See the Fortran chapter of the users manual and 2742 petsc/src/mat/examples/tests for details 2743 2744 .keywords: matrix, array, elements, values, restore 2745 2746 .seealso: MatGetArray() 2747 @*/ 2748 int MatRestoreArray(Mat mat,Scalar **v) 2749 { 2750 int ierr; 2751 2752 PetscFunctionBegin; 2753 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2754 PetscValidPointer(v); 2755 if (!mat->ops->restorearray) SETERRQ(PETSC_ERR_SUP,0,""); 2756 ierr = (*mat->ops->restorearray)(mat,v);CHKERRQ(ierr); 2757 PetscFunctionReturn(0); 2758 } 2759 2760 #undef __FUNC__ 2761 #define __FUNC__ "MatGetSubMatrices" 2762 /*@C 2763 MatGetSubMatrices - Extracts several submatrices from a matrix. If submat 2764 points to an array of valid matrices, they may be reused to store the new 2765 submatrices. 2766 2767 Collective on Mat 2768 2769 Input Parameters: 2770 + mat - the matrix 2771 . n - the number of submatrixes to be extracted 2772 . irow, icol - index sets of rows and columns to extract 2773 - scall - either MAT_INITIAL_MATRIX or MAT_REUSE_MATRIX 2774 2775 Output Parameter: 2776 . submat - the array of submatrices 2777 2778 Notes: 2779 MatGetSubMatrices() can extract only sequential submatrices 2780 (from both sequential and parallel matrices). Use MatGetSubMatrix() 2781 to extract a parallel submatrix. 2782 2783 When extracting submatrices from a parallel matrix, each processor can 2784 form a different submatrix by setting the rows and columns of its 2785 individual index sets according to the local submatrix desired. 2786 2787 When finished using the submatrices, the user should destroy 2788 them with MatDestroySubMatrices(). 2789 2790 MAT_REUSE_MATRIX can only be used when the nonzero structure of the 2791 original matrix has not changed from that last call to MatGetSubMatrices() 2792 2793 Fortran Note: 2794 The Fortran interface is slightly different from that given below, it 2795 requires one to pass in as submat a Mat (integer) array of size at least m. 2796 2797 .keywords: matrix, get, submatrix, submatrices 2798 2799 .seealso: MatDestroyMatrices(), MatGetSubMatrix() 2800 @*/ 2801 int MatGetSubMatrices(Mat mat,int n,IS *irow,IS *icol,MatReuse scall,Mat **submat) 2802 { 2803 int ierr; 2804 2805 PetscFunctionBegin; 2806 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2807 if (!mat->ops->getsubmatrices) SETERRQ(PETSC_ERR_SUP,0,""); 2808 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2809 2810 PLogEventBegin(MAT_GetSubMatrices,mat,0,0,0); 2811 ierr = (*mat->ops->getsubmatrices)(mat,n,irow,icol,scall,submat); CHKERRQ(ierr); 2812 PLogEventEnd(MAT_GetSubMatrices,mat,0,0,0); 2813 2814 PetscFunctionReturn(0); 2815 } 2816 2817 #undef __FUNC__ 2818 #define __FUNC__ "MatDestroyMatrices" 2819 /*@C 2820 MatDestroyMatrices - Destroys a set of matrices obtained with MatGetSubMatrices(). 2821 2822 Collective on Mat 2823 2824 Input Parameters: 2825 + n - the number of local matrices 2826 - mat - the matrices 2827 2828 .keywords: matrix, destroy, submatrix, submatrices 2829 2830 .seealso: MatGetSubMatrices() 2831 @*/ 2832 int MatDestroyMatrices(int n,Mat **mat) 2833 { 2834 int ierr,i; 2835 2836 PetscFunctionBegin; 2837 if (n < 0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,1,"Trying to destroy negative number of matrices %d",n); 2838 PetscValidPointer(mat); 2839 for ( i=0; i<n; i++ ) { 2840 ierr = MatDestroy((*mat)[i]); CHKERRQ(ierr); 2841 } 2842 if (n) PetscFree(*mat); 2843 PetscFunctionReturn(0); 2844 } 2845 2846 #undef __FUNC__ 2847 #define __FUNC__ "MatIncreaseOverlap" 2848 /*@ 2849 MatIncreaseOverlap - Given a set of submatrices indicated by index sets, 2850 replaces the index sets by larger ones that represent submatrices with 2851 additional overlap. 2852 2853 Collective on Mat 2854 2855 Input Parameters: 2856 + mat - the matrix 2857 . n - the number of index sets 2858 . is - the array of pointers to index sets 2859 - ov - the additional overlap requested 2860 2861 .keywords: matrix, overlap, Schwarz 2862 2863 .seealso: MatGetSubMatrices() 2864 @*/ 2865 int MatIncreaseOverlap(Mat mat,int n, IS *is,int ov) 2866 { 2867 int ierr; 2868 2869 PetscFunctionBegin; 2870 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2871 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2872 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2873 2874 if (ov == 0) PetscFunctionReturn(0); 2875 if (!mat->ops->increaseoverlap) SETERRQ(PETSC_ERR_SUP,0,""); 2876 PLogEventBegin(MAT_IncreaseOverlap,mat,0,0,0); 2877 ierr = (*mat->ops->increaseoverlap)(mat,n,is,ov); CHKERRQ(ierr); 2878 PLogEventEnd(MAT_IncreaseOverlap,mat,0,0,0); 2879 PetscFunctionReturn(0); 2880 } 2881 2882 #undef __FUNC__ 2883 #define __FUNC__ "MatPrintHelp" 2884 /*@ 2885 MatPrintHelp - Prints all the options for the matrix. 2886 2887 Collective on Mat 2888 2889 Input Parameter: 2890 . mat - the matrix 2891 2892 Options Database Keys: 2893 + -help - Prints matrix options 2894 - -h - Prints matrix options 2895 2896 .keywords: mat, help 2897 2898 .seealso: MatCreate(), MatCreateXXX() 2899 @*/ 2900 int MatPrintHelp(Mat mat) 2901 { 2902 static int called = 0; 2903 int ierr; 2904 MPI_Comm comm; 2905 2906 PetscFunctionBegin; 2907 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2908 2909 comm = mat->comm; 2910 if (!called) { 2911 (*PetscHelpPrintf)(comm,"General matrix options:\n"); 2912 (*PetscHelpPrintf)(comm," -mat_view_info: view basic matrix info during MatAssemblyEnd()\n"); 2913 (*PetscHelpPrintf)(comm," -mat_view_info_detailed: view detailed matrix info during MatAssemblyEnd()\n"); 2914 (*PetscHelpPrintf)(comm," -mat_view_draw: draw nonzero matrix structure during MatAssemblyEnd()\n"); 2915 (*PetscHelpPrintf)(comm," -draw_pause <sec>: set seconds of display pause\n"); 2916 (*PetscHelpPrintf)(comm," -display <name>: set alternate display\n"); 2917 called = 1; 2918 } 2919 if (mat->ops->printhelp) { 2920 ierr = (*mat->ops->printhelp)(mat); CHKERRQ(ierr); 2921 } 2922 PetscFunctionReturn(0); 2923 } 2924 2925 #undef __FUNC__ 2926 #define __FUNC__ "MatGetBlockSize" 2927 /*@ 2928 MatGetBlockSize - Returns the matrix block size; useful especially for the 2929 block row and block diagonal formats. 2930 2931 Not Collective 2932 2933 Input Parameter: 2934 . mat - the matrix 2935 2936 Output Parameter: 2937 . bs - block size 2938 2939 Notes: 2940 Block diagonal formats are MATSEQBDIAG, MATMPIBDIAG. 2941 Block row formats are MATSEQBAIJ, MATMPIBAIJ 2942 2943 .keywords: matrix, get, block, size 2944 2945 .seealso: MatCreateSeqBAIJ(), MatCreateMPIBAIJ(), MatCreateSeqBDiag(), MatCreateMPIBDiag() 2946 @*/ 2947 int MatGetBlockSize(Mat mat,int *bs) 2948 { 2949 int ierr; 2950 2951 PetscFunctionBegin; 2952 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2953 PetscValidIntPointer(bs); 2954 if (!mat->ops->getblocksize) SETERRQ(PETSC_ERR_SUP,0,""); 2955 ierr = (*mat->ops->getblocksize)(mat,bs);CHKERRQ(ierr); 2956 PetscFunctionReturn(0); 2957 } 2958 2959 #undef __FUNC__ 2960 #define __FUNC__ "MatGetRowIJ" 2961 /*@C 2962 MatGetRowIJ - Returns the compressed row storage i and j indices for sequential matrices. 2963 EXPERTS ONLY. 2964 2965 Collective on Mat 2966 2967 Input Parameters: 2968 + mat - the matrix 2969 . shift - 0 or 1 indicating we want the indices starting at 0 or 1 2970 - symmetric - PETSC_TRUE or PETSC_FALSE indicating the matrix data structure should be 2971 symmetrized 2972 2973 Output Parameters: 2974 + n - number of rows in the (possibly compressed) matrix 2975 . ia - the row pointers 2976 . ja - the column indices 2977 - done - PETSC_TRUE or PETSC_FALSE, indicating whether the values have been returned 2978 2979 .seealso: MatGetColumnIJ(), MatRestoreRowIJ() 2980 @*/ 2981 int MatGetRowIJ(Mat mat,int shift,PetscTruth symmetric,int *n,int **ia,int** ja,PetscTruth *done) 2982 { 2983 int ierr; 2984 2985 PetscFunctionBegin; 2986 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2987 if (ia) PetscValidIntPointer(ia); 2988 if (ja) PetscValidIntPointer(ja); 2989 PetscValidIntPointer(done); 2990 if (!mat->ops->getrowij) *done = PETSC_FALSE; 2991 else { 2992 *done = PETSC_TRUE; 2993 ierr = (*mat->ops->getrowij)(mat,shift,symmetric,n,ia,ja,done); CHKERRQ(ierr); 2994 } 2995 PetscFunctionReturn(0); 2996 } 2997 2998 #undef __FUNC__ 2999 #define __FUNC__ "MatGetColumnIJ" 3000 /*@C 3001 MatGetColumnIJ - Returns the compressed column storage i and j indices for sequential matrices. 3002 EXPERTS ONLY. 3003 3004 Collective on Mat 3005 3006 Input Parameters: 3007 + mat - the matrix 3008 . shift - 1 or zero indicating we want the indices starting at 0 or 1 3009 - symmetric - PETSC_TRUE or PETSC_FALSE indicating the matrix data structure should be 3010 symmetrized 3011 3012 Output Parameters: 3013 + n - number of columns in the (possibly compressed) matrix 3014 . ia - the column pointers 3015 . ja - the row indices 3016 - done - PETSC_TRUE or PETSC_FALSE, indicating whether the values have been returned 3017 3018 .seealso: MatGetRowIJ(), MatRestoreColumnIJ() 3019 @*/ 3020 int MatGetColumnIJ(Mat mat,int shift,PetscTruth symmetric,int *n,int **ia,int** ja,PetscTruth *done) 3021 { 3022 int ierr; 3023 3024 PetscFunctionBegin; 3025 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3026 if (ia) PetscValidIntPointer(ia); 3027 if (ja) PetscValidIntPointer(ja); 3028 PetscValidIntPointer(done); 3029 3030 if (!mat->ops->getcolumnij) *done = PETSC_FALSE; 3031 else { 3032 *done = PETSC_TRUE; 3033 ierr = (*mat->ops->getcolumnij)(mat,shift,symmetric,n,ia,ja,done); CHKERRQ(ierr); 3034 } 3035 PetscFunctionReturn(0); 3036 } 3037 3038 #undef __FUNC__ 3039 #define __FUNC__ "MatRestoreRowIJ" 3040 /*@C 3041 MatRestoreRowIJ - Call after you are completed with the ia,ja indices obtained with 3042 MatGetRowIJ(). EXPERTS ONLY. 3043 3044 Collective on Mat 3045 3046 Input Parameters: 3047 + mat - the matrix 3048 . shift - 1 or zero indicating we want the indices starting at 0 or 1 3049 - symmetric - PETSC_TRUE or PETSC_FALSE indicating the matrix data structure should be 3050 symmetrized 3051 3052 Output Parameters: 3053 + n - size of (possibly compressed) matrix 3054 . ia - the row pointers 3055 . ja - the column indices 3056 - done - PETSC_TRUE or PETSC_FALSE indicated that the values have been returned 3057 3058 .seealso: MatGetRowIJ(), MatRestoreColumnIJ() 3059 @*/ 3060 int MatRestoreRowIJ(Mat mat,int shift,PetscTruth symmetric,int *n,int **ia,int** ja,PetscTruth *done) 3061 { 3062 int ierr; 3063 3064 PetscFunctionBegin; 3065 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3066 if (ia) PetscValidIntPointer(ia); 3067 if (ja) PetscValidIntPointer(ja); 3068 PetscValidIntPointer(done); 3069 3070 if (!mat->ops->restorerowij) *done = PETSC_FALSE; 3071 else { 3072 *done = PETSC_TRUE; 3073 ierr = (*mat->ops->restorerowij)(mat,shift,symmetric,n,ia,ja,done); CHKERRQ(ierr); 3074 } 3075 PetscFunctionReturn(0); 3076 } 3077 3078 #undef __FUNC__ 3079 #define __FUNC__ "MatRestoreColumnIJ" 3080 /*@C 3081 MatRestoreColumnIJ - Call after you are completed with the ia,ja indices obtained with 3082 MatGetColumnIJ(). EXPERTS ONLY. 3083 3084 Collective on Mat 3085 3086 Input Parameters: 3087 + mat - the matrix 3088 . shift - 1 or zero indicating we want the indices starting at 0 or 1 3089 - symmetric - PETSC_TRUE or PETSC_FALSE indicating the matrix data structure should be 3090 symmetrized 3091 3092 Output Parameters: 3093 + n - size of (possibly compressed) matrix 3094 . ia - the column pointers 3095 . ja - the row indices 3096 - done - PETSC_TRUE or PETSC_FALSE indicated that the values have been returned 3097 3098 .seealso: MatGetColumnIJ(), MatRestoreRowIJ() 3099 @*/ 3100 int MatRestoreColumnIJ(Mat mat,int shift,PetscTruth symmetric,int *n,int **ia,int** ja,PetscTruth *done) 3101 { 3102 int ierr; 3103 3104 PetscFunctionBegin; 3105 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3106 if (ia) PetscValidIntPointer(ia); 3107 if (ja) PetscValidIntPointer(ja); 3108 PetscValidIntPointer(done); 3109 3110 if (!mat->ops->restorecolumnij) *done = PETSC_FALSE; 3111 else { 3112 *done = PETSC_TRUE; 3113 ierr = (*mat->ops->restorecolumnij)(mat,shift,symmetric,n,ia,ja,done); CHKERRQ(ierr); 3114 } 3115 PetscFunctionReturn(0); 3116 } 3117 3118 #undef __FUNC__ 3119 #define __FUNC__ "MatColoringPatch" 3120 /*@C 3121 MatColoringPatch - EXPERTS ONLY, used inside matrix coloring routines that 3122 use MatGetRowIJ() and/or MatGetColumnIJ(). 3123 3124 Collective on Mat 3125 3126 Input Parameters: 3127 + mat - the matrix 3128 . n - number of colors 3129 - colorarray - array indicating color for each column 3130 3131 Output Parameters: 3132 . iscoloring - coloring generated using colorarray information 3133 3134 .seealso: MatGetRowIJ(), MatGetColumnIJ() 3135 3136 .keywords: mat, coloring, patch 3137 @*/ 3138 int MatColoringPatch(Mat mat,int n,int *colorarray,ISColoring *iscoloring) 3139 { 3140 int ierr; 3141 3142 PetscFunctionBegin; 3143 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3144 PetscValidIntPointer(colorarray); 3145 3146 if (!mat->ops->coloringpatch) {SETERRQ(PETSC_ERR_SUP,0,"");} 3147 else { 3148 ierr = (*mat->ops->coloringpatch)(mat,n,colorarray,iscoloring); CHKERRQ(ierr); 3149 } 3150 PetscFunctionReturn(0); 3151 } 3152 3153 3154 #undef __FUNC__ 3155 #define __FUNC__ "MatSetUnfactored" 3156 /*@ 3157 MatSetUnfactored - Resets a factored matrix to be treated as unfactored. 3158 3159 Collective on Mat 3160 3161 Input Parameter: 3162 . mat - the factored matrix to be reset 3163 3164 Notes: 3165 This routine should be used only with factored matrices formed by in-place 3166 factorization via ILU(0) (or by in-place LU factorization for the MATSEQDENSE 3167 format). This option can save memory, for example, when solving nonlinear 3168 systems with a matrix-free Newton-Krylov method and a matrix-based, in-place 3169 ILU(0) preconditioner. 3170 3171 Note that one can specify in-place ILU(0) factorization by calling 3172 $ PCType(pc,PCILU); 3173 $ PCILUSeUseInPlace(pc); 3174 or by using the options -pc_type ilu -pc_ilu_in_place 3175 3176 In-place factorization ILU(0) can also be used as a local 3177 solver for the blocks within the block Jacobi or additive Schwarz 3178 methods (runtime option: -sub_pc_ilu_in_place). See the discussion 3179 of these preconditioners in the users manual for details on setting 3180 local solver options. 3181 3182 .seealso: PCILUSetUseInPlace(), PCLUSetUseInPlace() 3183 3184 .keywords: matrix-free, in-place ILU, in-place LU 3185 @*/ 3186 int MatSetUnfactored(Mat mat) 3187 { 3188 int ierr; 3189 3190 PetscFunctionBegin; 3191 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3192 mat->factor = 0; 3193 if (!mat->ops->setunfactored) PetscFunctionReturn(0); 3194 ierr = (*mat->ops->setunfactored)(mat); CHKERRQ(ierr); 3195 PetscFunctionReturn(0); 3196 } 3197 3198 #undef __FUNC__ 3199 #define __FUNC__ "MatGetType" 3200 /*@C 3201 MatGetType - Gets the matrix type and name (as a string) from the matrix. 3202 3203 Not Collective 3204 3205 Input Parameter: 3206 . mat - the matrix 3207 3208 Output Parameter: 3209 + type - the matrix type (or use PETSC_NULL) 3210 - name - name of matrix type (or use PETSC_NULL) 3211 3212 .keywords: matrix, get, type, name 3213 @*/ 3214 int MatGetType(Mat mat,MatType *type,char **name) 3215 { 3216 int itype = (int)mat->type; 3217 char *matname[10]; 3218 3219 PetscFunctionBegin; 3220 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3221 3222 if (type) *type = (MatType) mat->type; 3223 if (name) { 3224 /* Note: Be sure that this list corresponds to the enum in mat.h */ 3225 matname[0] = "MATSEQDENSE"; 3226 matname[1] = "MATSEQAIJ"; 3227 matname[2] = "MATMPIAIJ"; 3228 matname[3] = "MATSHELL"; 3229 matname[4] = "MATMPIROWBS"; 3230 matname[5] = "MATSEQBDIAG"; 3231 matname[6] = "MATMPIBDIAG"; 3232 matname[7] = "MATMPIDENSE"; 3233 matname[8] = "MATSEQBAIJ"; 3234 matname[9] = "MATMPIBAIJ"; 3235 3236 if (itype < 0 || itype > 9) *name = "Unknown matrix type"; 3237 else *name = matname[itype]; 3238 } 3239 PetscFunctionReturn(0); 3240 } 3241 3242 /*MC 3243 MatGetArrayF90 - Accesses a matrix array from Fortran90. 3244 3245 Synopsis: 3246 MatGetArrayF90(Mat x,{Scalar, pointer :: xx_v(:)},integer ierr) 3247 3248 Not collective 3249 3250 Input Parameter: 3251 . x - matrix 3252 3253 Output Parameters: 3254 + xx_v - the Fortran90 pointer to the array 3255 - ierr - error code 3256 3257 Example of Usage: 3258 .vb 3259 Scalar, pointer xx_v(:) 3260 .... 3261 call MatGetArrayF90(x,xx_v,ierr) 3262 a = xx_v(3) 3263 call MatRestoreArrayF90(x,xx_v,ierr) 3264 .ve 3265 3266 Notes: 3267 Not yet supported for all F90 compilers 3268 3269 .seealso: MatRestoreArrayF90(), MatGetArray(), MatRestoreArray() 3270 3271 .keywords: matrix, array, f90 3272 M*/ 3273 3274 /*MC 3275 MatRestoreArrayF90 - Restores a matrix array that has been 3276 accessed with MatGetArrayF90(). 3277 3278 Synopsis: 3279 MatRestoreArrayF90(Mat x,{Scalar, pointer :: xx_v(:)},integer ierr) 3280 3281 Not collective 3282 3283 Input Parameters: 3284 + x - matrix 3285 - xx_v - the Fortran90 pointer to the array 3286 3287 Output Parameter: 3288 . ierr - error code 3289 3290 Example of Usage: 3291 .vb 3292 Scalar, pointer xx_v(:) 3293 .... 3294 call MatGetArrayF90(x,xx_v,ierr) 3295 a = xx_v(3) 3296 call MatRestoreArrayF90(x,xx_v,ierr) 3297 .ve 3298 3299 Notes: 3300 Not yet supported for all F90 compilers 3301 3302 .seealso: MatGetArrayF90(), MatGetArray(), MatRestoreArray() 3303 3304 .keywords: matrix, array, f90 3305 M*/ 3306 3307 3308 #undef __FUNC__ 3309 #define __FUNC__ "MatGetSubMatrix" 3310 /*@ 3311 MatGetSubMatrix - Gets a single submatrix on the same number of processors 3312 as the original matrix. 3313 3314 Collective on Mat 3315 3316 Input Parameters: 3317 + mat - the original matrix 3318 . isrow - rows this processor should obtain 3319 . iscol - columns for all processors you wish to keep 3320 . csize - number of columns "local" to this processor (does nothing for sequential 3321 matrices). This should match the result from VecGetLocalSize(x,...) if you 3322 plan to use the matrix in a A*x 3323 - cll - either MAT_INITIAL_MATRIX or MAT_REUSE_MATRIX 3324 3325 Output Parameter: 3326 . newmat - the new submatrix, of the same type as the old 3327 3328 .seealso: MatGetSubMatrices() 3329 3330 @*/ 3331 int MatGetSubMatrix(Mat mat,IS isrow,IS iscol,int csize,MatReuse cll,Mat *newmat) 3332 { 3333 int ierr, size; 3334 Mat *local; 3335 3336 PetscFunctionBegin; 3337 MPI_Comm_size(mat->comm,&size); 3338 3339 /* if original matrix is on just one processor then use submatrix generated */ 3340 if (!mat->ops->getsubmatrix && size == 1 && cll == MAT_REUSE_MATRIX) { 3341 ierr = MatGetSubMatrices(mat,1,&isrow,&iscol,MAT_REUSE_MATRIX,&newmat);CHKERRQ(ierr); 3342 PetscFunctionReturn(0); 3343 } else if (!mat->ops->getsubmatrix && size == 1) { 3344 ierr = MatGetSubMatrices(mat,1,&isrow,&iscol,MAT_INITIAL_MATRIX,&local);CHKERRQ(ierr); 3345 *newmat = *local; 3346 PetscFree(local); 3347 PetscFunctionReturn(0); 3348 } 3349 3350 if (!mat->ops->getsubmatrix) SETERRQ(PETSC_ERR_SUP,0,"Not currently implemented"); 3351 ierr = (*mat->ops->getsubmatrix)(mat,isrow,iscol,csize,cll,newmat);CHKERRQ(ierr); 3352 PetscFunctionReturn(0); 3353 } 3354 3355 #undef __FUNC__ 3356 #define __FUNC__ "MatGetMaps" 3357 /*@C 3358 MatGetMaps - Returns the maps associated with the matrix. 3359 3360 Not Collective 3361 3362 Input Parameter: 3363 . mat - the matrix 3364 3365 Output Parameters: 3366 + rmap - the row (right) map 3367 - cmap - the column (left) map 3368 3369 .keywords: matrix, get, map 3370 @*/ 3371 int MatGetMaps(Mat mat,Map *rmap,Map *cmap) 3372 { 3373 int ierr; 3374 3375 PetscFunctionBegin; 3376 PetscValidHeaderSpecific(mat,MAT_COOKIE); 3377 ierr = (*mat->ops->getmaps)(mat,rmap,cmap);CHKERRQ(ierr); 3378 PetscFunctionReturn(0); 3379 } 3380 3381 /* 3382 Version that works for all PETSc matrices 3383 */ 3384 #undef __FUNC__ 3385 #define __FUNC__ "MatGetMaps_Petsc" 3386 int MatGetMaps_Petsc(Mat mat,Map *rmap,Map *cmap) 3387 { 3388 PetscFunctionBegin; 3389 if (rmap) *rmap = mat->rmap; 3390 if (cmap) *cmap = mat->cmap; 3391 PetscFunctionReturn(0); 3392 } 3393