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