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