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