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