1 2 #ifndef lint 3 static char vcid[] = "$Id: matrix.c,v 1.227 1997/02/22 02:24:34 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 $ and possibly others. 1885 1886 Notes: 1887 Some options are relevant only for particular matrix types and 1888 are thus ignored by others. Other options are not supported by 1889 certain matrix types and will generate an error message if set. 1890 1891 If using a Fortran 77 module to compute a matrix, one may need to 1892 use the column-oriented option (or convert to the row-oriented 1893 format). 1894 1895 MAT_NO_NEW_NONZERO_LOCATIONS indicates that any add or insertion 1896 that will generate a new entry in the nonzero structure is ignored. 1897 Thus, if memory has not alredy been allocated for this particular 1898 data, then the insertion is ignored. For dense matrices, in which 1899 the entire array is allocated, no entries are ever ignored. 1900 1901 MAT_IGNORE_OFF_PROCESSOR_ENTRIES indicates entries destined for 1902 other processors are dropped, rather than stashed. 1903 1904 .keywords: matrix, option, row-oriented, column-oriented, sorted, nonzero 1905 @*/ 1906 int MatSetOption(Mat mat,MatOption op) 1907 { 1908 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1909 if (mat->ops.setoption) return (*mat->ops.setoption)(mat,op); 1910 return 0; 1911 } 1912 1913 #undef __FUNC__ 1914 #define __FUNC__ "MatZeroEntries" 1915 /*@ 1916 MatZeroEntries - Zeros all entries of a matrix. For sparse matrices 1917 this routine retains the old nonzero structure. 1918 1919 Input Parameters: 1920 . mat - the matrix 1921 1922 .keywords: matrix, zero, entries 1923 1924 .seealso: MatZeroRows() 1925 @*/ 1926 int MatZeroEntries(Mat mat) 1927 { 1928 int ierr; 1929 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1930 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 1931 if (!mat->ops.zeroentries) SETERRQ(PETSC_ERR_SUP,0,""); 1932 1933 PLogEventBegin(MAT_ZeroEntries,mat,0,0,0); 1934 ierr = (*mat->ops.zeroentries)(mat); CHKERRQ(ierr); 1935 PLogEventEnd(MAT_ZeroEntries,mat,0,0,0); 1936 return 0; 1937 } 1938 1939 #undef __FUNC__ 1940 #define __FUNC__ "MatZeroRows" 1941 /*@ 1942 MatZeroRows - Zeros all entries (except possibly the main diagonal) 1943 of a set of rows of a matrix. 1944 1945 Input Parameters: 1946 . mat - the matrix 1947 . is - index set of rows to remove 1948 . diag - pointer to value put in all diagonals of eliminated rows. 1949 Note that diag is not a pointer to an array, but merely a 1950 pointer to a single value. 1951 1952 Notes: 1953 For the AIJ matrix formats this removes the old nonzero structure, 1954 but does not release memory. For the dense and block diagonal 1955 formats this does not alter the nonzero structure. 1956 1957 The user can set a value in the diagonal entry (or for the AIJ and 1958 row formats can optionally remove the main diagonal entry from the 1959 nonzero structure as well, by passing a null pointer as the final 1960 argument). 1961 1962 .keywords: matrix, zero, rows, boundary conditions 1963 1964 .seealso: MatZeroEntries(), 1965 @*/ 1966 int MatZeroRows(Mat mat,IS is, Scalar *diag) 1967 { 1968 int ierr; 1969 1970 PetscValidHeaderSpecific(mat,MAT_COOKIE); 1971 PetscValidHeaderSpecific(is,IS_COOKIE); 1972 if (diag) PetscValidScalarPointer(diag); 1973 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 1974 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 1975 if (!mat->ops.zerorows) SETERRQ(PETSC_ERR_SUP,0,""); 1976 1977 ierr = (*mat->ops.zerorows)(mat,is,diag); CHKERRQ(ierr); 1978 ierr = MatView_Private(mat); CHKERRQ(ierr); 1979 return 0; 1980 } 1981 1982 #undef __FUNC__ 1983 #define __FUNC__ "MatZeroRowsLocal" 1984 /*@ 1985 MatZeroRowsLocal - Zeros all entries (except possibly the main diagonal) 1986 of a set of rows of a matrix; using local numbering of rows. 1987 1988 Input Parameters: 1989 . mat - the matrix 1990 . is - index set of rows to remove 1991 . diag - pointer to value put in all diagonals of eliminated rows. 1992 Note that diag is not a pointer to an array, but merely a 1993 pointer to a single value. 1994 1995 Notes: 1996 For the AIJ matrix formats this removes the old nonzero structure, 1997 but does not release memory. For the dense and block diagonal 1998 formats this does not alter the nonzero structure. 1999 2000 The user can set a value in the diagonal entry (or for the AIJ and 2001 row formats can optionally remove the main diagonal entry from the 2002 nonzero structure as well, by passing a null pointer as the final 2003 argument). 2004 2005 .keywords: matrix, zero, rows, boundary conditions 2006 2007 .seealso: MatZeroEntries(), 2008 @*/ 2009 int MatZeroRowsLocal(Mat mat,IS is, Scalar *diag) 2010 { 2011 int ierr; 2012 IS newis; 2013 2014 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2015 PetscValidHeaderSpecific(is,IS_COOKIE); 2016 if (diag) PetscValidScalarPointer(diag); 2017 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2018 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2019 if (!mat->ops.zerorows) SETERRQ(PETSC_ERR_SUP,0,""); 2020 2021 ierr = ISLocalToGlobalMappingApplyIS(mat->mapping,is,&newis); CHKERRQ(ierr); 2022 ierr = (*mat->ops.zerorows)(mat,newis,diag); CHKERRQ(ierr); 2023 ierr = ISDestroy(newis); 2024 return 0; 2025 } 2026 2027 #undef __FUNC__ 2028 #define __FUNC__ "MatGetSize" /* ADIC Ignore */ 2029 /*@ 2030 MatGetSize - Returns the numbers of rows and columns in a matrix. 2031 2032 Input Parameter: 2033 . mat - the matrix 2034 2035 Output Parameters: 2036 . m - the number of global rows 2037 . n - the number of global columns 2038 2039 .keywords: matrix, dimension, size, rows, columns, global, get 2040 2041 .seealso: MatGetLocalSize() 2042 @*/ 2043 int MatGetSize(Mat mat,int *m,int* n) 2044 { 2045 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2046 PetscValidIntPointer(m); 2047 PetscValidIntPointer(n); 2048 return (*mat->ops.getsize)(mat,m,n); 2049 } 2050 2051 #undef __FUNC__ 2052 #define __FUNC__ "MatGetLocalSize" /* ADIC Ignore */ 2053 /*@ 2054 MatGetLocalSize - Returns the number of rows and columns in a matrix 2055 stored locally. This information may be implementation dependent, so 2056 use with care. 2057 2058 Input Parameters: 2059 . mat - the matrix 2060 2061 Output Parameters: 2062 . m - the number of local rows 2063 . n - the number of local columns 2064 2065 .keywords: matrix, dimension, size, local, rows, columns, get 2066 2067 .seealso: MatGetSize() 2068 @*/ 2069 int MatGetLocalSize(Mat mat,int *m,int* n) 2070 { 2071 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2072 PetscValidIntPointer(m); 2073 PetscValidIntPointer(n); 2074 return (*mat->ops.getlocalsize)(mat,m,n); 2075 } 2076 2077 #undef __FUNC__ 2078 #define __FUNC__ "MatGetOwnershipRange" /* ADIC Ignore */ 2079 /*@ 2080 MatGetOwnershipRange - Returns the range of matrix rows owned by 2081 this processor, assuming that the matrix is laid out with the first 2082 n1 rows on the first processor, the next n2 rows on the second, etc. 2083 For certain parallel layouts this range may not be well defined. 2084 2085 Input Parameters: 2086 . mat - the matrix 2087 2088 Output Parameters: 2089 . m - the first local row 2090 . n - one more then the last local row 2091 2092 .keywords: matrix, get, range, ownership 2093 @*/ 2094 int MatGetOwnershipRange(Mat mat,int *m,int* n) 2095 { 2096 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2097 PetscValidIntPointer(m); 2098 PetscValidIntPointer(n); 2099 if (!mat->ops.getownershiprange) SETERRQ(PETSC_ERR_SUP,0,""); 2100 return (*mat->ops.getownershiprange)(mat,m,n); 2101 } 2102 2103 #undef __FUNC__ 2104 #define __FUNC__ "MatILUFactorSymbolic" 2105 /*@ 2106 MatILUFactorSymbolic - Performs symbolic ILU factorization of a matrix. 2107 Uses levels of fill only, not drop tolerance. Use MatLUFactorNumeric() 2108 to complete the factorization. 2109 2110 Input Parameters: 2111 . mat - the matrix 2112 . row - row permutation 2113 . column - column permutation 2114 . fill - number of levels of fill 2115 . f - expected fill as ratio of the original number of nonzeros, 2116 for example 3.0; choosing this parameter well can result in 2117 more efficient use of time and space. 2118 2119 Output Parameters: 2120 . fact - new matrix that has been symbolically factored 2121 2122 Options Database Key: 2123 $ -mat_ilu_fill <f>, where f is the fill ratio 2124 2125 Notes: 2126 See the file $(PETSC_DIR)/Performace for additional information about 2127 choosing the fill factor for better efficiency. 2128 2129 .keywords: matrix, factor, incomplete, ILU, symbolic, fill 2130 2131 .seealso: MatLUFactorSymbolic(), MatLUFactorNumeric() 2132 @*/ 2133 int MatILUFactorSymbolic(Mat mat,IS row,IS col,double f,int fill,Mat *fact) 2134 { 2135 int ierr,flg; 2136 2137 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2138 PetscValidPointer(fact); 2139 if (fill < 0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,0,"Levels of fill negative"); 2140 if (!mat->ops.ilufactorsymbolic) SETERRQ(PETSC_ERR_SUP,0,""); 2141 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2142 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2143 2144 ierr = OptionsGetDouble(PETSC_NULL,"-mat_ilu_fill",&f,&flg); CHKERRQ(ierr); 2145 PLogEventBegin(MAT_ILUFactorSymbolic,mat,row,col,0); 2146 ierr = (*mat->ops.ilufactorsymbolic)(mat,row,col,f,fill,fact); CHKERRQ(ierr); 2147 PLogEventEnd(MAT_ILUFactorSymbolic,mat,row,col,0); 2148 return 0; 2149 } 2150 2151 #undef __FUNC__ 2152 #define __FUNC__ "MatIncompleteCholeskyFactorSymbolic" 2153 /*@ 2154 MatIncompleteCholeskyFactorSymbolic - Performs symbolic incomplete 2155 Cholesky factorization for a symmetric matrix. Use 2156 MatCholeskyFactorNumeric() to complete the factorization. 2157 2158 Input Parameters: 2159 . mat - the matrix 2160 . perm - row and column permutation 2161 . fill - levels of fill 2162 . f - expected fill as ratio of original fill 2163 2164 Output Parameter: 2165 . fact - the factored matrix 2166 2167 Note: Currently only no-fill factorization is supported. 2168 2169 .keywords: matrix, factor, incomplete, ICC, Cholesky, symbolic, fill 2170 2171 .seealso: MatCholeskyFactorNumeric(), MatCholeskyFactor() 2172 @*/ 2173 int MatIncompleteCholeskyFactorSymbolic(Mat mat,IS perm,double f,int fill, 2174 Mat *fact) 2175 { 2176 int ierr; 2177 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2178 PetscValidPointer(fact); 2179 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2180 if (fill < 0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,0,"Fill negative"); 2181 if (!mat->ops.incompletecholeskyfactorsymbolic) SETERRQ(PETSC_ERR_SUP,0,""); 2182 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2183 2184 PLogEventBegin(MAT_IncompleteCholeskyFactorSymbolic,mat,perm,0,0); 2185 ierr = (*mat->ops.incompletecholeskyfactorsymbolic)(mat,perm,f,fill,fact);CHKERRQ(ierr); 2186 PLogEventEnd(MAT_IncompleteCholeskyFactorSymbolic,mat,perm,0,0); 2187 return 0; 2188 } 2189 2190 #undef __FUNC__ 2191 #define __FUNC__ "MatGetArray" /* ADIC Ignore */ 2192 /*@C 2193 MatGetArray - Returns a pointer to the element values in the matrix. 2194 This routine is implementation dependent, and may not even work for 2195 certain matrix types. You MUST call MatRestoreArray() when you no 2196 longer need to access the array. 2197 2198 Input Parameter: 2199 . mat - the matrix 2200 2201 Output Parameter: 2202 . v - the location of the values 2203 2204 Fortran Note: 2205 The Fortran interface is slightly different from that given below. 2206 See the Fortran chapter of the users manual and 2207 petsc/src/mat/examples for details. 2208 2209 .keywords: matrix, array, elements, values 2210 2211 .seealso: MatRestoreArray() 2212 @*/ 2213 int MatGetArray(Mat mat,Scalar **v) 2214 { 2215 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2216 PetscValidPointer(v); 2217 if (!mat->ops.getarray) SETERRQ(PETSC_ERR_SUP,0,""); 2218 return (*mat->ops.getarray)(mat,v); 2219 } 2220 2221 #undef __FUNC__ 2222 #define __FUNC__ "MatRestoreArray" /* ADIC Ignore */ 2223 /*@C 2224 MatRestoreArray - Restores the matrix after MatGetArray has been called. 2225 2226 Input Parameter: 2227 . mat - the matrix 2228 . v - the location of the values 2229 2230 Fortran Note: 2231 The Fortran interface is slightly different from that given below. 2232 See the users manual and petsc/src/mat/examples for details. 2233 2234 .keywords: matrix, array, elements, values, restore 2235 2236 .seealso: MatGetArray() 2237 @*/ 2238 int MatRestoreArray(Mat mat,Scalar **v) 2239 { 2240 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2241 PetscValidPointer(v); 2242 if (!mat->ops.restorearray) SETERRQ(PETSC_ERR_SUP,0,""); 2243 return (*mat->ops.restorearray)(mat,v); 2244 } 2245 2246 #undef __FUNC__ 2247 #define __FUNC__ "MatGetSubMatrices" 2248 /*@C 2249 MatGetSubMatrices - Extracts several submatrices from a matrix. If submat 2250 points to an array of valid matrices, it may be reused. 2251 2252 Input Parameters: 2253 . mat - the matrix 2254 . n - the number of submatrixes to be extracted 2255 . irow, icol - index sets of rows and columns to extract 2256 . scall - either MAT_INITIAL_MATRIX or MAT_REUSE_MATRIX 2257 2258 Output Parameter: 2259 . submat - the array of submatrices 2260 2261 Limitations: 2262 Currently, MatGetSubMatrices() can extract only sequential submatrices 2263 (from both sequential and parallel matrices). 2264 2265 Notes: 2266 When extracting submatrices from a parallel matrix, each processor can 2267 form a different submatrix by setting the rows and columns of its 2268 individual index sets according to the local submatrix desired. 2269 2270 When finished using the submatrices, the user should destroy 2271 them with MatDestroySubMatrices(). 2272 2273 .keywords: matrix, get, submatrix, submatrices 2274 2275 .seealso: MatDestroyMatrices() 2276 @*/ 2277 int MatGetSubMatrices(Mat mat,int n,IS *irow,IS *icol,MatGetSubMatrixCall scall, 2278 Mat **submat) 2279 { 2280 int ierr; 2281 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2282 if (!mat->ops.getsubmatrices) SETERRQ(PETSC_ERR_SUP,0,""); 2283 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2284 2285 PLogEventBegin(MAT_GetSubMatrices,mat,0,0,0); 2286 ierr = (*mat->ops.getsubmatrices)(mat,n,irow,icol,scall,submat); CHKERRQ(ierr); 2287 PLogEventEnd(MAT_GetSubMatrices,mat,0,0,0); 2288 2289 return 0; 2290 } 2291 2292 #undef __FUNC__ 2293 #define __FUNC__ "MatDestroyMatrices" /* ADIC Ignore */ 2294 /*@C 2295 MatDestroyMatrices - Destroys a set of matrices obtained with MatGetSubMatrices(). 2296 2297 Input Parameters: 2298 . n - the number of local matrices 2299 . mat - the matrices 2300 2301 .keywords: matrix, destroy, submatrix, submatrices 2302 2303 .seealso: MatGetSubMatrices() 2304 @*/ 2305 int MatDestroyMatrices(int n,Mat **mat) 2306 { 2307 int ierr,i; 2308 2309 if (n < 0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,1,"Trying to destroy negative number of matrices"); 2310 PetscValidPointer(mat); 2311 for ( i=0; i<n; i++ ) { 2312 ierr = MatDestroy((*mat)[i]); CHKERRQ(ierr); 2313 } 2314 if (n) PetscFree(*mat); 2315 return 0; 2316 } 2317 2318 #undef __FUNC__ 2319 #define __FUNC__ "MatIncreaseOverlap" 2320 /*@ 2321 MatIncreaseOverlap - Given a set of submatrices indicated by index sets, 2322 replaces the index sets by larger ones that represent submatrices with 2323 additional overlap. 2324 2325 Input Parameters: 2326 . mat - the matrix 2327 . n - the number of index sets 2328 . is - the array of pointers to index sets 2329 . ov - the additional overlap requested 2330 2331 .keywords: matrix, overlap, Schwarz 2332 2333 .seealso: MatGetSubMatrices() 2334 @*/ 2335 int MatIncreaseOverlap(Mat mat,int n, IS *is,int ov) 2336 { 2337 int ierr; 2338 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2339 if (!mat->assembled) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for unassembled matrix"); 2340 if (mat->factor) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,0,"Not for factored matrix"); 2341 2342 if (ov == 0) return 0; 2343 if (!mat->ops.increaseoverlap) SETERRQ(PETSC_ERR_SUP,0,""); 2344 PLogEventBegin(MAT_IncreaseOverlap,mat,0,0,0); 2345 ierr = (*mat->ops.increaseoverlap)(mat,n,is,ov); CHKERRQ(ierr); 2346 PLogEventEnd(MAT_IncreaseOverlap,mat,0,0,0); 2347 return 0; 2348 } 2349 2350 #undef __FUNC__ 2351 #define __FUNC__ "MatPrintHelp" /* ADIC Ignore */ 2352 /*@ 2353 MatPrintHelp - Prints all the options for the matrix. 2354 2355 Input Parameter: 2356 . mat - the matrix 2357 2358 Options Database Keys: 2359 $ -help, -h 2360 2361 .keywords: mat, help 2362 2363 .seealso: MatCreate(), MatCreateXXX() 2364 @*/ 2365 int MatPrintHelp(Mat mat) 2366 { 2367 static int called = 0; 2368 MPI_Comm comm; 2369 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2370 2371 comm = mat->comm; 2372 if (!called) { 2373 PetscPrintf(comm,"General matrix options:\n"); 2374 PetscPrintf(comm," -mat_view_info: view basic matrix info during MatAssemblyEnd()\n"); 2375 PetscPrintf(comm," -mat_view_info_detailed: view detailed matrix info during MatAssemblyEnd()\n"); 2376 PetscPrintf(comm," -mat_view_draw: draw nonzero matrix structure during MatAssemblyEnd()\n"); 2377 PetscPrintf(comm," -draw_pause <sec>: set seconds of display pause\n"); 2378 PetscPrintf(comm," -display <name>: set alternate display\n"); 2379 called = 1; 2380 } 2381 if (mat->ops.printhelp) (*mat->ops.printhelp)(mat); 2382 return 0; 2383 } 2384 2385 #undef __FUNC__ 2386 #define __FUNC__ "MatGetBlockSize" /* ADIC Ignore */ 2387 /*@ 2388 MatGetBlockSize - Returns the matrix block size; useful especially for the 2389 block row and block diagonal formats. 2390 2391 Input Parameter: 2392 . mat - the matrix 2393 2394 Output Parameter: 2395 . bs - block size 2396 2397 Notes: 2398 $ block diagonal formats: MATSEQBDIAG, MATMPIBDIAG 2399 $ block row formats: MATSEQBAIJ, MATMPIBAIJ 2400 2401 .keywords: matrix, get, block, size 2402 2403 .seealso: MatCreateSeqBAIJ(), MatCreateMPIBAIJ(), MatCreateSeqBDiag(), MatCreateMPIBDiag() 2404 @*/ 2405 int MatGetBlockSize(Mat mat,int *bs) 2406 { 2407 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2408 PetscValidIntPointer(bs); 2409 if (!mat->ops.getblocksize) SETERRQ(PETSC_ERR_SUP,0,""); 2410 return (*mat->ops.getblocksize)(mat,bs); 2411 } 2412 2413 #undef __FUNC__ 2414 #define __FUNC__ "MatGetRowIJ" 2415 /*@C 2416 MatGetRowIJ - Returns the compress row storage i and j indices for sequential matrices. 2417 EXPERTS ONLY. 2418 2419 Input Parameters: 2420 . mat - the matrix 2421 . shift - 1 or zero indicating we want the indices starting at 0 or 1 2422 . symmetric - PETSC_TRUE or PETSC_FALSE indicating the matrix data structure should be 2423 symmetrized 2424 2425 Output Parameters: 2426 . n - number of rows and columns in the (possibly compressed) matrix 2427 . ia - the row indices 2428 . ja - the column indices 2429 . done - PETSC_TRUE or PETSC_FALSE indicated that the values have been returned 2430 @*/ 2431 int MatGetRowIJ(Mat mat,int shift,PetscTruth symmetric,int *n,int **ia,int** ja,PetscTruth *done) 2432 { 2433 int ierr; 2434 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2435 if (ia) PetscValidIntPointer(ia); 2436 if (ja) PetscValidIntPointer(ja); 2437 PetscValidIntPointer(done); 2438 if (!mat->ops.getrowij) *done = PETSC_FALSE; 2439 else { 2440 *done = PETSC_TRUE; 2441 ierr = (*mat->ops.getrowij)(mat,shift,symmetric,n,ia,ja,done); CHKERRQ(ierr); 2442 } 2443 return 0; 2444 } 2445 2446 #undef __FUNC__ 2447 #define __FUNC__ "MatGetColumnIJ" /* ADIC Ignore */ 2448 /*@C 2449 MatGetColumnIJ - Returns the compress Column storage i and j indices for sequential matrices. 2450 EXPERTS ONLY. 2451 2452 Input Parameters: 2453 . mat - the matrix 2454 . shift - 1 or zero indicating we want the indices starting at 0 or 1 2455 . symmetric - PETSC_TRUE or PETSC_FALSE indicating the matrix data structure should be 2456 symmetrized 2457 2458 Output Parameters: 2459 . n - number of Columns and columns in the (possibly compressed) matrix 2460 . ia - the Column indices 2461 . ja - the column indices 2462 . done - PETSC_TRUE or PETSC_FALSE indicated that the values have been returned 2463 @*/ 2464 int MatGetColumnIJ(Mat mat,int shift,PetscTruth symmetric,int *n,int **ia,int** ja,PetscTruth *done) 2465 { 2466 int ierr; 2467 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2468 if (ia) PetscValidIntPointer(ia); 2469 if (ja) PetscValidIntPointer(ja); 2470 PetscValidIntPointer(done); 2471 2472 if (!mat->ops.getcolumnij) *done = PETSC_FALSE; 2473 else { 2474 *done = PETSC_TRUE; 2475 ierr = (*mat->ops.getcolumnij)(mat,shift,symmetric,n,ia,ja,done); CHKERRQ(ierr); 2476 } 2477 return 0; 2478 } 2479 2480 #undef __FUNC__ 2481 #define __FUNC__ "MatRestoreRowIJ" /* ADIC Ignore */ 2482 /*@C 2483 MatRestoreRowIJ - Call after you are completed with the ia,ja indices obtained with 2484 MatGetRowIJ(). EXPERTS ONLY. 2485 2486 Input Parameters: 2487 . mat - the matrix 2488 . shift - 1 or zero indicating we want the indices starting at 0 or 1 2489 . symmetric - PETSC_TRUE or PETSC_FALSE indicating the matrix data structure should be 2490 symmetrized 2491 2492 Output Parameters: 2493 . n - size of (possibly compressed) matrix 2494 . ia - the row indices 2495 . ja - the column indices 2496 . done - PETSC_TRUE or PETSC_FALSE indicated that the values have been returned 2497 @*/ 2498 int MatRestoreRowIJ(Mat mat,int shift,PetscTruth symmetric,int *n,int **ia,int** ja,PetscTruth *done) 2499 { 2500 int ierr; 2501 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2502 if (ia) PetscValidIntPointer(ia); 2503 if (ja) PetscValidIntPointer(ja); 2504 PetscValidIntPointer(done); 2505 2506 if (!mat->ops.restorerowij) *done = PETSC_FALSE; 2507 else { 2508 *done = PETSC_TRUE; 2509 ierr = (*mat->ops.restorerowij)(mat,shift,symmetric,n,ia,ja,done); CHKERRQ(ierr); 2510 } 2511 return 0; 2512 } 2513 2514 #undef __FUNC__ 2515 #define __FUNC__ "MatRestoreColumnIJ" /* ADIC Ignore */ 2516 /*@C 2517 MatRestoreColumnIJ - Call after you are completed with the ia,ja indices obtained with 2518 MatGetColumnIJ(). EXPERTS ONLY. 2519 2520 Input Parameters: 2521 . mat - the matrix 2522 . shift - 1 or zero indicating we want the indices starting at 0 or 1 2523 . symmetric - PETSC_TRUE or PETSC_FALSE indicating the matrix data structure should be 2524 symmetrized 2525 2526 Output Parameters: 2527 . n - size of (possibly compressed) matrix 2528 . ia - the Column indices 2529 . ja - the column indices 2530 . done - PETSC_TRUE or PETSC_FALSE indicated that the values have been returned 2531 @*/ 2532 int MatRestoreColumnIJ(Mat mat,int shift,PetscTruth symmetric,int *n,int **ia,int** ja,PetscTruth *done) 2533 { 2534 int ierr; 2535 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2536 if (ia) PetscValidIntPointer(ia); 2537 if (ja) PetscValidIntPointer(ja); 2538 PetscValidIntPointer(done); 2539 2540 if (!mat->ops.restorecolumnij) *done = PETSC_FALSE; 2541 else { 2542 *done = PETSC_TRUE; 2543 ierr = (*mat->ops.restorecolumnij)(mat,shift,symmetric,n,ia,ja,done); CHKERRQ(ierr); 2544 } 2545 return 0; 2546 } 2547 2548 #undef __FUNC__ 2549 #define __FUNC__ "MatColoringPatch" 2550 /*@C 2551 MatColoringPatch - EXPERTS ONLY, used inside matrix coloring routines that 2552 use matGetRowIJ() and/or MatGetColumnIJ(). 2553 2554 Input Parameters: 2555 . mat - the matrix 2556 . n - number of colors 2557 . colorarray - array indicating color for each column 2558 2559 Output Parameters: 2560 . iscoloring - coloring generated using colorarray information 2561 2562 @*/ 2563 int MatColoringPatch(Mat mat,int n,int *colorarray,ISColoring *iscoloring) 2564 { 2565 int ierr; 2566 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2567 PetscValidIntPointer(colorarray); 2568 2569 if (!mat->ops.coloringpatch) {SETERRQ(PETSC_ERR_SUP,0,"");} 2570 else { 2571 ierr = (*mat->ops.coloringpatch)(mat,n,colorarray,iscoloring); CHKERRQ(ierr); 2572 } 2573 return 0; 2574 } 2575 2576 2577 /*@ 2578 MatSetUnfactored - Resets a factored matrix to be treated as unfactored. 2579 2580 Input Paramter: 2581 . mat - the factored matrix to be reset 2582 2583 Notes: 2584 This routine should be used only with factored matrices formed by in-place 2585 factorization via ILU(0) (or by in-place LU factorization for the MATSEQDENSE 2586 format). This option can save memory, for example, when solving nonlinear 2587 systems with a matrix-free Newton-Krylov method and a matrix-based, in-place 2588 ILU(0) preconditioner. 2589 2590 Note that one can specify in-place ILU(0) factorization by calling 2591 $ PCType(pc,PCILU); 2592 $ PCILUSeUseInPlace(pc); 2593 or by using the options -pc_type ilu -pc_ilu_in_place 2594 2595 In-place factorization ILU(0) can also be used as a local 2596 solver for the blocks within the block Jacobi or additive Schwarz 2597 methods (runtime option: -sub_pc_ilu_in_place). See the discussion 2598 of these preconditioners in the users manual for details on setting 2599 local solver options. 2600 2601 .seealso: PCILUSetUseInPlace(), PCLUSetUseInPlace() 2602 2603 .keywords: matrix-free, in-place ILU, in-place LU 2604 @*/ 2605 int MatSetUnfactored(Mat mat) 2606 { 2607 int ierr; 2608 2609 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2610 mat->factor = 0; 2611 if (!mat->ops.setunfactored) return 0; 2612 ierr = (*mat->ops.setunfactored)(mat); CHKERRQ(ierr); 2613 return 0; 2614 } 2615 2616 #undef __FUNC__ 2617 #define __FUNC__ "MatGetType" /* ADIC Ignore */ 2618 /*@C 2619 MatGetType - Gets the matrix type and name (as a string) from the matrix. 2620 2621 Input Parameter: 2622 . mat - the matrix 2623 2624 Output Parameter: 2625 . type - the matrix type (or use PETSC_NULL) 2626 . name - name of matrix type (or use PETSC_NULL) 2627 2628 .keywords: matrix, get, type, name 2629 @*/ 2630 int MatGetType(Mat mat,MatType *type,char **name) 2631 { 2632 int itype = (int)mat->type; 2633 char *matname[10]; 2634 2635 PetscValidHeaderSpecific(mat,MAT_COOKIE); 2636 2637 if (type) *type = (MatType) mat->type; 2638 if (name) { 2639 /* Note: Be sure that this list corresponds to the enum in mat.h */ 2640 matname[0] = "MATSEQDENSE"; 2641 matname[1] = "MATSEQAIJ"; 2642 matname[2] = "MATMPIAIJ"; 2643 matname[3] = "MATSHELL"; 2644 matname[4] = "MATMPIROWBS"; 2645 matname[5] = "MATSEQBDIAG"; 2646 matname[6] = "MATMPIBDIAG"; 2647 matname[7] = "MATMPIDENSE"; 2648 matname[8] = "MATSEQBAIJ"; 2649 matname[9] = "MATMPIBAIJ"; 2650 2651 if (itype < 0 || itype > 9) *name = "Unknown matrix type"; 2652 else *name = matname[itype]; 2653 } 2654 return 0; 2655 } 2656