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