1 2 /* 3 The basic KSP routines, Create, View etc. are here. 4 */ 5 #include <petsc-private/kspimpl.h> /*I "petscksp.h" I*/ 6 7 /* Logging support */ 8 PetscClassId KSP_CLASSID; 9 PetscLogEvent KSP_GMRESOrthogonalization, KSP_SetUp, KSP_Solve; 10 11 /* 12 Contains the list of registered KSP routines 13 */ 14 PetscFList KSPList = 0; 15 PetscBool KSPRegisterAllCalled = PETSC_FALSE; 16 17 #undef __FUNCT__ 18 #define __FUNCT__ "KSPLoad" 19 /*@C 20 KSPLoad - Loads a KSP that has been stored in binary with KSPView(). 21 22 Collective on PetscViewer 23 24 Input Parameters: 25 + newdm - the newly loaded KSP, this needs to have been created with KSPCreate() or 26 some related function before a call to KSPLoad(). 27 - viewer - binary file viewer, obtained from PetscViewerBinaryOpen() 28 29 Level: intermediate 30 31 Notes: 32 The type is determined by the data in the file, any type set into the KSP before this call is ignored. 33 34 Notes for advanced users: 35 Most users should not need to know the details of the binary storage 36 format, since KSPLoad() and KSPView() completely hide these details. 37 But for anyone who's interested, the standard binary matrix storage 38 format is 39 .vb 40 has not yet been determined 41 .ve 42 43 .seealso: PetscViewerBinaryOpen(), KSPView(), MatLoad(), VecLoad() 44 @*/ 45 PetscErrorCode KSPLoad(KSP newdm, PetscViewer viewer) 46 { 47 PetscErrorCode ierr; 48 PetscBool isbinary; 49 PetscInt classid; 50 char type[256]; 51 PC pc; 52 53 PetscFunctionBegin; 54 PetscValidHeaderSpecific(newdm,KSP_CLASSID,1); 55 PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,2); 56 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);CHKERRQ(ierr); 57 if (!isbinary) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Invalid viewer; open viewer with PetscViewerBinaryOpen()"); 58 59 ierr = PetscViewerBinaryRead(viewer,&classid,1,PETSC_INT);CHKERRQ(ierr); 60 if (classid != KSP_FILE_CLASSID) SETERRQ(((PetscObject)newdm)->comm,PETSC_ERR_ARG_WRONG,"Not KSP next in file"); 61 ierr = PetscViewerBinaryRead(viewer,type,256,PETSC_CHAR);CHKERRQ(ierr); 62 ierr = KSPSetType(newdm, type);CHKERRQ(ierr); 63 ierr = (*newdm->ops->load)(newdm,viewer);CHKERRQ(ierr); 64 ierr = KSPGetPC(newdm,&pc);CHKERRQ(ierr); 65 ierr = PCLoad(pc,viewer);CHKERRQ(ierr); 66 PetscFunctionReturn(0); 67 } 68 69 #undef __FUNCT__ 70 #define __FUNCT__ "KSPView" 71 /*@C 72 KSPView - Prints the KSP data structure. 73 74 Collective on KSP 75 76 Input Parameters: 77 + ksp - the Krylov space context 78 - viewer - visualization context 79 80 Options Database Keys: 81 . -ksp_view - print the ksp data structure at the end of a KSPSolve call 82 83 Note: 84 The available visualization contexts include 85 + PETSC_VIEWER_STDOUT_SELF - standard output (default) 86 - PETSC_VIEWER_STDOUT_WORLD - synchronized standard 87 output where only the first processor opens 88 the file. All other processors send their 89 data to the first processor to print. 90 91 The user can open an alternative visualization context with 92 PetscViewerASCIIOpen() - output to a specified file. 93 94 Level: beginner 95 96 .keywords: KSP, view 97 98 .seealso: PCView(), PetscViewerASCIIOpen() 99 @*/ 100 PetscErrorCode KSPView(KSP ksp,PetscViewer viewer) 101 { 102 PetscErrorCode ierr; 103 PetscBool iascii,isbinary; 104 105 PetscFunctionBegin; 106 PetscValidHeaderSpecific(ksp,KSP_CLASSID,1); 107 if (!viewer) viewer = PETSC_VIEWER_STDOUT_(((PetscObject)ksp)->comm); 108 PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,2); 109 PetscCheckSameComm(ksp,1,viewer,2); 110 111 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr); 112 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);CHKERRQ(ierr); 113 if (iascii) { 114 ierr = PetscObjectPrintClassNamePrefixType((PetscObject)ksp,viewer,"KSP Object");CHKERRQ(ierr); 115 if (ksp->ops->view) { 116 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 117 ierr = (*ksp->ops->view)(ksp,viewer);CHKERRQ(ierr); 118 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 119 } 120 if (ksp->guess_zero) {ierr = PetscViewerASCIIPrintf(viewer," maximum iterations=%D, initial guess is zero\n",ksp->max_it);CHKERRQ(ierr);} 121 else {ierr = PetscViewerASCIIPrintf(viewer," maximum iterations=%D\n", ksp->max_it);CHKERRQ(ierr);} 122 if (ksp->guess_knoll) {ierr = PetscViewerASCIIPrintf(viewer," using preconditioner applied to right hand side for initial guess\n");CHKERRQ(ierr);} 123 ierr = PetscViewerASCIIPrintf(viewer," tolerances: relative=%G, absolute=%G, divergence=%G\n",ksp->rtol,ksp->abstol,ksp->divtol);CHKERRQ(ierr); 124 if (ksp->pc_side == PC_RIGHT) {ierr = PetscViewerASCIIPrintf(viewer," right preconditioning\n");CHKERRQ(ierr);} 125 else if (ksp->pc_side == PC_SYMMETRIC) {ierr = PetscViewerASCIIPrintf(viewer," symmetric preconditioning\n");CHKERRQ(ierr);} 126 else {ierr = PetscViewerASCIIPrintf(viewer," left preconditioning\n");CHKERRQ(ierr);} 127 if (ksp->guess) {ierr = PetscViewerASCIIPrintf(viewer," using Fischers initial guess method %D with size %D\n",ksp->guess->method,ksp->guess->maxl);CHKERRQ(ierr);} 128 if (ksp->dscale) {ierr = PetscViewerASCIIPrintf(viewer," diagonally scaled system\n");CHKERRQ(ierr);} 129 if (ksp->nullsp) {ierr = PetscViewerASCIIPrintf(viewer," has attached null space\n");CHKERRQ(ierr);} 130 if (!ksp->guess_zero) {ierr = PetscViewerASCIIPrintf(viewer," using nonzero initial guess\n");CHKERRQ(ierr);} 131 ierr = PetscViewerASCIIPrintf(viewer," using %s norm type for convergence test\n",KSPNormTypes[ksp->normtype]);CHKERRQ(ierr); 132 } else if (isbinary) { 133 PetscInt classid = KSP_FILE_CLASSID; 134 MPI_Comm comm; 135 PetscMPIInt rank; 136 char type[256]; 137 138 ierr = PetscObjectGetComm((PetscObject)ksp,&comm);CHKERRQ(ierr); 139 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr); 140 if (!rank) { 141 ierr = PetscViewerBinaryWrite(viewer,&classid,1,PETSC_INT,PETSC_FALSE);CHKERRQ(ierr); 142 ierr = PetscStrncpy(type,((PetscObject)ksp)->type_name,256);CHKERRQ(ierr); 143 ierr = PetscViewerBinaryWrite(viewer,type,256,PETSC_CHAR,PETSC_FALSE);CHKERRQ(ierr); 144 } 145 if (ksp->ops->view) { 146 ierr = (*ksp->ops->view)(ksp,viewer);CHKERRQ(ierr); 147 } 148 } else { 149 if (ksp->ops->view) { 150 ierr = (*ksp->ops->view)(ksp,viewer);CHKERRQ(ierr); 151 } 152 } 153 if (!ksp->pc) {ierr = KSPGetPC(ksp,&ksp->pc);CHKERRQ(ierr);} 154 ierr = PCView(ksp->pc,viewer);CHKERRQ(ierr); 155 PetscFunctionReturn(0); 156 } 157 158 159 #undef __FUNCT__ 160 #define __FUNCT__ "KSPSetNormType" 161 /*@ 162 KSPSetNormType - Sets the norm that is used for convergence testing. 163 164 Logically Collective on KSP 165 166 Input Parameter: 167 + ksp - Krylov solver context 168 - normtype - one of 169 $ KSP_NORM_NONE - skips computing the norm, this should only be used if you are using 170 $ the Krylov method as a smoother with a fixed small number of iterations. 171 $ Implicitly sets KSPSkipConverged as KSP convergence test. 172 $ Supported only by CG, Richardson, Bi-CG-stab, CR, and CGS methods. 173 $ KSP_NORM_PRECONDITIONED - the default for left preconditioned solves, uses the l2 norm 174 $ of the preconditioned residual 175 $ KSP_NORM_UNPRECONDITIONED - uses the l2 norm of the true b - Ax residual, supported only by 176 $ CG, CHEBYSHEV, and RICHARDSON, automatically true for right (see KSPSetPCSide()) 177 $ preconditioning.. 178 $ KSP_NORM_NATURAL - supported by KSPCG, KSPCR, KSPCGNE, KSPCGS 179 180 181 Options Database Key: 182 . -ksp_norm_type <none,preconditioned,unpreconditioned,natural> 183 184 Notes: 185 Currently only works with the CG, Richardson, Bi-CG-stab, CR, and CGS methods. 186 187 Level: advanced 188 189 .keywords: KSP, create, context, norms 190 191 .seealso: KSPSetUp(), KSPSolve(), KSPDestroy(), KSPSkipConverged(), KSPSetCheckNormIteration() 192 @*/ 193 PetscErrorCode KSPSetNormType(KSP ksp,KSPNormType normtype) 194 { 195 PetscErrorCode ierr; 196 197 PetscFunctionBegin; 198 PetscValidHeaderSpecific(ksp,KSP_CLASSID,1); 199 PetscValidLogicalCollectiveEnum(ksp,normtype,2); 200 ksp->normtype = normtype; 201 if (normtype == KSP_NORM_NONE) { 202 ierr = KSPSetConvergenceTest(ksp,KSPSkipConverged,0,0);CHKERRQ(ierr); 203 ierr = PetscInfo(ksp,"Warning: setting KSPNormType to skip computing the norm\n\ 204 KSP convergence test is implicitly set to KSPSkipConverged\n");CHKERRQ(ierr); 205 } 206 PetscFunctionReturn(0); 207 } 208 209 #undef __FUNCT__ 210 #define __FUNCT__ "KSPSetCheckNormIteration" 211 /*@ 212 KSPSetCheckNormIteration - Sets the first iteration at which the norm of the residual will be 213 computed and used in the convergence test. 214 215 Logically Collective on KSP 216 217 Input Parameter: 218 + ksp - Krylov solver context 219 - it - use -1 to check at all iterations 220 221 Notes: 222 Currently only works with KSPCG, KSPBCGS and KSPIBCGS 223 224 Use KSPSetNormType(ksp,KSP_NORM_NONE) to never check the norm 225 226 On steps where the norm is not computed, the previous norm is still in the variable, so if you run with, for example, 227 -ksp_monitor the residual norm will appear to be unchanged for several iterations (though it is not really unchanged). 228 Level: advanced 229 230 .keywords: KSP, create, context, norms 231 232 .seealso: KSPSetUp(), KSPSolve(), KSPDestroy(), KSPSkipConverged(), KSPSetNormType() 233 @*/ 234 PetscErrorCode KSPSetCheckNormIteration(KSP ksp,PetscInt it) 235 { 236 PetscFunctionBegin; 237 PetscValidHeaderSpecific(ksp,KSP_CLASSID,1); 238 PetscValidLogicalCollectiveInt(ksp,it,2); 239 ksp->chknorm = it; 240 PetscFunctionReturn(0); 241 } 242 243 #undef __FUNCT__ 244 #define __FUNCT__ "KSPSetLagNorm" 245 /*@ 246 KSPSetLagNorm - Lags the residual norm calculation so that it is computed as part of the MPI_Allreduce() for 247 computing the inner products for the next iteration. This can reduce communication costs at the expense of doing 248 one additional iteration. 249 250 251 Logically Collective on KSP 252 253 Input Parameter: 254 + ksp - Krylov solver context 255 - flg - PETSC_TRUE or PETSC_FALSE 256 257 Options Database Keys: 258 . -ksp_lag_norm - lag the calculated residual norm 259 260 Notes: 261 Currently only works with KSPIBCGS. 262 263 Use KSPSetNormType(ksp,KSP_NORM_NONE) to never check the norm 264 265 If you lag the norm and run with, for example, -ksp_monitor, the residual norm reported will be the lagged one. 266 Level: advanced 267 268 .keywords: KSP, create, context, norms 269 270 .seealso: KSPSetUp(), KSPSolve(), KSPDestroy(), KSPSkipConverged(), KSPSetNormType(), KSPSetCheckNormIteration() 271 @*/ 272 PetscErrorCode KSPSetLagNorm(KSP ksp,PetscBool flg) 273 { 274 PetscFunctionBegin; 275 PetscValidHeaderSpecific(ksp,KSP_CLASSID,1); 276 PetscValidLogicalCollectiveBool(ksp,flg,2); 277 ksp->lagnorm = flg; 278 PetscFunctionReturn(0); 279 } 280 281 #undef __FUNCT__ 282 #define __FUNCT__ "KSPSetSupportedNorm" 283 /*@ 284 KSPSetSupportedNorm - Sets a norm and preconditioner side supported by a KSP 285 286 Logically Collective 287 288 Input Arguments: 289 + ksp - Krylov method 290 . normtype - supported norm type 291 . pcside - preconditioner side that can be used with this norm 292 - preference - integer preference for this combination, larger values have higher priority 293 294 Level: developer 295 296 Notes: 297 This function should be called from the implementation files KSPCreate_XXX() to declare 298 which norms and preconditioner sides are supported. Users should not need to call this 299 function. 300 301 KSP_NORM_NONE is supported by default with all KSP methods and any PC side. If a KSP explicitly does not support 302 KSP_NORM_NONE, it should set this by setting priority=0. 303 304 .seealso: KSPSetNormType(), KSPSetPCSide() 305 @*/ 306 PetscErrorCode KSPSetSupportedNorm(KSP ksp,KSPNormType normtype,PCSide pcside,PetscInt priority) 307 { 308 309 PetscFunctionBegin; 310 PetscValidHeaderSpecific(ksp,KSP_CLASSID,1); 311 ksp->normsupporttable[normtype][pcside] = priority; 312 PetscFunctionReturn(0); 313 } 314 315 #undef __FUNCT__ 316 #define __FUNCT__ "KSPNormSupportTableReset_Private" 317 PetscErrorCode KSPNormSupportTableReset_Private(KSP ksp) 318 { 319 PetscErrorCode ierr; 320 321 PetscFunctionBegin; 322 ierr = PetscMemzero(ksp->normsupporttable,sizeof(ksp->normsupporttable));CHKERRQ(ierr); 323 ierr = KSPSetSupportedNorm(ksp,KSP_NORM_NONE,PC_LEFT,1);CHKERRQ(ierr); 324 ierr = KSPSetSupportedNorm(ksp,KSP_NORM_NONE,PC_RIGHT,1);CHKERRQ(ierr); 325 PetscFunctionReturn(0); 326 } 327 328 #undef __FUNCT__ 329 #define __FUNCT__ "KSPSetUpNorms_Private" 330 PetscErrorCode KSPSetUpNorms_Private(KSP ksp,KSPNormType *normtype,PCSide *pcside) 331 { 332 PetscInt i,j,best,ibest = 0,jbest = 0; 333 334 PetscFunctionBegin; 335 best = 0; 336 for (i=0; i<KSP_NORM_MAX; i++) { 337 for (j=0; j<PC_SIDE_MAX; j++) { 338 if ((ksp->normtype == KSP_NORM_DEFAULT || ksp->normtype == i) 339 && (ksp->pc_side == PC_SIDE_DEFAULT || ksp->pc_side == j) 340 && (ksp->normsupporttable[i][j] > best)) { 341 if (ksp->normtype == KSP_NORM_DEFAULT && i == KSP_NORM_NONE && ksp->normsupporttable[i][j] <= 1) 342 continue; /* Skip because we don't want to default to no norms unless set by the KSP (preonly). */ 343 best = ksp->normsupporttable[i][j]; 344 ibest = i; 345 jbest = j; 346 } 347 } 348 } 349 if (best < 1) { 350 if (ksp->normtype == KSP_NORM_DEFAULT && ksp->pc_side == PC_SIDE_DEFAULT) SETERRQ1(((PetscObject)ksp)->comm,PETSC_ERR_PLIB,"The %s KSP implementation did not call KSPSetSupportedNorm()",((PetscObject)ksp)->type_name); 351 if (ksp->normtype == KSP_NORM_DEFAULT) SETERRQ2(((PetscObject)ksp)->comm,PETSC_ERR_SUP,"KSP %s does not support %s",((PetscObject)ksp)->type_name,PCSides[ksp->pc_side]); 352 if (ksp->pc_side == PC_SIDE_DEFAULT) SETERRQ2(((PetscObject)ksp)->comm,PETSC_ERR_SUP,"KSP %s does not support %s",((PetscObject)ksp)->type_name,KSPNormTypes[ksp->normtype]); 353 SETERRQ3(((PetscObject)ksp)->comm,PETSC_ERR_SUP,"KSP %s does not support %s with %s",((PetscObject)ksp)->type_name,KSPNormTypes[ksp->normtype],PCSides[ksp->pc_side]); 354 } 355 *normtype = (KSPNormType)ibest; 356 *pcside = (PCSide)jbest; 357 PetscFunctionReturn(0); 358 } 359 360 #undef __FUNCT__ 361 #define __FUNCT__ "KSPGetNormType" 362 /*@ 363 KSPGetNormType - Gets the norm that is used for convergence testing. 364 365 Not Collective 366 367 Input Parameter: 368 . ksp - Krylov solver context 369 370 Output Parameter: 371 . normtype - norm that is used for convergence testing 372 373 Level: advanced 374 375 .keywords: KSP, create, context, norms 376 377 .seealso: KSPNormType, KSPSetNormType(), KSPSkipConverged() 378 @*/ 379 PetscErrorCode KSPGetNormType(KSP ksp, KSPNormType *normtype) 380 { 381 PetscErrorCode ierr; 382 383 PetscFunctionBegin; 384 PetscValidHeaderSpecific(ksp,KSP_CLASSID,1); 385 PetscValidPointer(normtype,2); 386 ierr = KSPSetUpNorms_Private(ksp,&ksp->normtype,&ksp->pc_side);CHKERRQ(ierr); 387 *normtype = ksp->normtype; 388 PetscFunctionReturn(0); 389 } 390 391 #if defined(PETSC_HAVE_AMS) 392 #undef __FUNCT__ 393 #define __FUNCT__ "KSPPublish_Petsc" 394 static PetscErrorCode KSPPublish_Petsc(PetscObject obj) 395 { 396 KSP ksp = (KSP) obj; 397 PetscErrorCode ierr; 398 399 PetscFunctionBegin; 400 ierr = AMS_Memory_add_field(obj->amem,"its",&ksp->its,1,AMS_INT,AMS_READ,AMS_COMMON,AMS_REDUCT_UNDEF);CHKERRQ(ierr); 401 PetscFunctionReturn(0); 402 } 403 #endif 404 405 #undef __FUNCT__ 406 #define __FUNCT__ "KSPSetOperators" 407 /*@ 408 KSPSetOperators - Sets the matrix associated with the linear system 409 and a (possibly) different one associated with the preconditioner. 410 411 Collective on KSP and Mat 412 413 Input Parameters: 414 + ksp - the KSP context 415 . Amat - the matrix associated with the linear system 416 . Pmat - the matrix to be used in constructing the preconditioner, usually the 417 same as Amat. 418 - flag - flag indicating information about the preconditioner matrix structure 419 during successive linear solves. This flag is ignored the first time a 420 linear system is solved, and thus is irrelevant when solving just one linear 421 system. 422 423 Notes: 424 The flag can be used to eliminate unnecessary work in the preconditioner 425 during the repeated solution of linear systems of the same size. The 426 available options are 427 $ SAME_PRECONDITIONER - 428 $ Pmat is identical during successive linear solves. 429 $ This option is intended for folks who are using 430 $ different Amat and Pmat matrices and want to reuse the 431 $ same preconditioner matrix. For example, this option 432 $ saves work by not recomputing incomplete factorization 433 $ for ILU/ICC preconditioners. 434 $ SAME_NONZERO_PATTERN - 435 $ Pmat has the same nonzero structure during 436 $ successive linear solves. 437 $ DIFFERENT_NONZERO_PATTERN - 438 $ Pmat does not have the same nonzero structure. 439 440 All future calls to KSPSetOperators() must use the same size matrices! 441 442 Passing a PETSC_NULL for Amat or Pmat removes the matrix that is currently used. 443 444 If you wish to replace either Amat or Pmat but leave the other one untouched then 445 first call KSPGetOperators() to get the one you wish to keep, call PetscObjectReference() 446 on it and then pass it back in in your call to KSPSetOperators(). 447 448 Caution: 449 If you specify SAME_NONZERO_PATTERN, PETSc believes your assertion 450 and does not check the structure of the matrix. If you erroneously 451 claim that the structure is the same when it actually is not, the new 452 preconditioner will not function correctly. Thus, use this optimization 453 feature carefully! 454 455 If in doubt about whether your preconditioner matrix has changed 456 structure or not, use the flag DIFFERENT_NONZERO_PATTERN. 457 458 Level: beginner 459 460 Alternative usage: If the operators have NOT been set with KSP/PCSetOperators() then the operators 461 are created in PC and returned to the user. In this case, if both operators 462 mat and pmat are requested, two DIFFERENT operators will be returned. If 463 only one is requested both operators in the PC will be the same (i.e. as 464 if one had called KSP/PCSetOperators() with the same argument for both Mats). 465 The user must set the sizes of the returned matrices and their type etc just 466 as if the user created them with MatCreate(). For example, 467 468 $ KSP/PCGetOperators(ksp/pc,&mat,PETSC_NULL,PETSC_NULL); is equivalent to 469 $ set size, type, etc of mat 470 471 $ MatCreate(comm,&mat); 472 $ KSP/PCSetOperators(ksp/pc,mat,mat,SAME_NONZERO_PATTERN); 473 $ PetscObjectDereference((PetscObject)mat); 474 $ set size, type, etc of mat 475 476 and 477 478 $ KSP/PCGetOperators(ksp/pc,&mat,&pmat,PETSC_NULL); is equivalent to 479 $ set size, type, etc of mat and pmat 480 481 $ MatCreate(comm,&mat); 482 $ MatCreate(comm,&pmat); 483 $ KSP/PCSetOperators(ksp/pc,mat,pmat,SAME_NONZERO_PATTERN); 484 $ PetscObjectDereference((PetscObject)mat); 485 $ PetscObjectDereference((PetscObject)pmat); 486 $ set size, type, etc of mat and pmat 487 488 The rational for this support is so that when creating a TS, SNES, or KSP the hierarchy 489 of underlying objects (i.e. SNES, KSP, PC, Mat) and their livespans can be completely 490 managed by the top most level object (i.e. the TS, SNES, or KSP). Another way to look 491 at this is when you create a SNES you do not NEED to create a KSP and attach it to 492 the SNES object (the SNES object manages it for you). Similarly when you create a KSP 493 you do not need to attach a PC to it (the KSP object manages the PC object for you). 494 Thus, why should YOU have to create the Mat and attach it to the SNES/KSP/PC, when 495 it can be created for you? 496 497 .keywords: KSP, set, operators, matrix, preconditioner, linear system 498 499 .seealso: KSPSolve(), KSPGetPC(), PCGetOperators(), PCSetOperators(), KSPGetOperators() 500 @*/ 501 PetscErrorCode KSPSetOperators(KSP ksp,Mat Amat,Mat Pmat,MatStructure flag) 502 { 503 MatNullSpace nullsp; 504 PetscErrorCode ierr; 505 506 PetscFunctionBegin; 507 PetscValidHeaderSpecific(ksp,KSP_CLASSID,1); 508 if (Amat) PetscValidHeaderSpecific(Amat,MAT_CLASSID,2); 509 if (Pmat) PetscValidHeaderSpecific(Pmat,MAT_CLASSID,3); 510 if (Amat) PetscCheckSameComm(ksp,1,Amat,2); 511 if (Pmat) PetscCheckSameComm(ksp,1,Pmat,3); 512 if (!ksp->pc) {ierr = KSPGetPC(ksp,&ksp->pc);CHKERRQ(ierr);} 513 ierr = PCSetOperators(ksp->pc,Amat,Pmat,flag);CHKERRQ(ierr); 514 if (ksp->setupstage == KSP_SETUP_NEWRHS) ksp->setupstage = KSP_SETUP_NEWMATRIX; /* so that next solve call will call PCSetUp() on new matrix */ 515 if (ksp->guess) { 516 ierr = KSPFischerGuessReset(ksp->guess);CHKERRQ(ierr); 517 } 518 if (Pmat) { 519 ierr = MatGetNullSpace(Pmat, &nullsp);CHKERRQ(ierr); 520 if (nullsp) { 521 ierr = KSPSetNullSpace(ksp, nullsp);CHKERRQ(ierr); 522 } 523 } 524 PetscFunctionReturn(0); 525 } 526 527 #undef __FUNCT__ 528 #define __FUNCT__ "KSPGetOperators" 529 /*@ 530 KSPGetOperators - Gets the matrix associated with the linear system 531 and a (possibly) different one associated with the preconditioner. 532 533 Collective on KSP and Mat 534 535 Input Parameter: 536 . ksp - the KSP context 537 538 Output Parameters: 539 + Amat - the matrix associated with the linear system 540 . Pmat - the matrix to be used in constructing the preconditioner, usually the same as Amat. 541 - flag - flag indicating information about the preconditioner matrix structure 542 during successive linear solves. This flag is ignored the first time a 543 linear system is solved, and thus is irrelevant when solving just one linear 544 system. 545 546 Level: intermediate 547 548 Notes: DOES NOT increase the reference counts of the matrix, so you should NOT destroy them. 549 550 .keywords: KSP, set, get, operators, matrix, preconditioner, linear system 551 552 .seealso: KSPSolve(), KSPGetPC(), PCGetOperators(), PCSetOperators(), KSPSetOperators(), KSPGetOperatorsSet() 553 @*/ 554 PetscErrorCode KSPGetOperators(KSP ksp,Mat *Amat,Mat *Pmat,MatStructure *flag) 555 { 556 PetscErrorCode ierr; 557 558 PetscFunctionBegin; 559 PetscValidHeaderSpecific(ksp,KSP_CLASSID,1); 560 if (!ksp->pc) {ierr = KSPGetPC(ksp,&ksp->pc);CHKERRQ(ierr);} 561 ierr = PCGetOperators(ksp->pc,Amat,Pmat,flag);CHKERRQ(ierr); 562 PetscFunctionReturn(0); 563 } 564 565 #undef __FUNCT__ 566 #define __FUNCT__ "KSPGetOperatorsSet" 567 /*@C 568 KSPGetOperatorsSet - Determines if the matrix associated with the linear system and 569 possibly a different one associated with the preconditioner have been set in the KSP. 570 571 Not collective, though the results on all processes should be the same 572 573 Input Parameter: 574 . pc - the KSP context 575 576 Output Parameters: 577 + mat - the matrix associated with the linear system was set 578 - pmat - matrix associated with the preconditioner was set, usually the same 579 580 Level: intermediate 581 582 .keywords: KSP, get, operators, matrix, linear system 583 584 .seealso: PCSetOperators(), KSPGetOperators(), KSPSetOperators(), PCGetOperators(), PCGetOperatorsSet() 585 @*/ 586 PetscErrorCode KSPGetOperatorsSet(KSP ksp,PetscBool *mat,PetscBool *pmat) 587 { 588 PetscErrorCode ierr; 589 590 PetscFunctionBegin; 591 PetscValidHeaderSpecific(ksp,KSP_CLASSID,1); 592 if (!ksp->pc) {ierr = KSPGetPC(ksp,&ksp->pc);CHKERRQ(ierr);} 593 ierr = PCGetOperatorsSet(ksp->pc,mat,pmat);CHKERRQ(ierr); 594 PetscFunctionReturn(0); 595 } 596 597 #undef __FUNCT__ 598 #define __FUNCT__ "KSPCreate" 599 /*@ 600 KSPCreate - Creates the default KSP context. 601 602 Collective on MPI_Comm 603 604 Input Parameter: 605 . comm - MPI communicator 606 607 Output Parameter: 608 . ksp - location to put the KSP context 609 610 Notes: 611 The default KSP type is GMRES with a restart of 30, using modified Gram-Schmidt 612 orthogonalization. 613 614 Level: beginner 615 616 .keywords: KSP, create, context 617 618 .seealso: KSPSetUp(), KSPSolve(), KSPDestroy(), KSP 619 @*/ 620 PetscErrorCode KSPCreate(MPI_Comm comm,KSP *inksp) 621 { 622 KSP ksp; 623 PetscErrorCode ierr; 624 void *ctx; 625 626 PetscFunctionBegin; 627 PetscValidPointer(inksp,2); 628 *inksp = 0; 629 #ifndef PETSC_USE_DYNAMIC_LIBRARIES 630 ierr = KSPInitializePackage(PETSC_NULL);CHKERRQ(ierr); 631 #endif 632 633 ierr = PetscHeaderCreate(ksp,_p_KSP,struct _KSPOps,KSP_CLASSID,-1,"KSP","Krylov Method","KSP",comm,KSPDestroy,KSPView);CHKERRQ(ierr); 634 635 ksp->max_it = 10000; 636 ksp->pc_side = PC_SIDE_DEFAULT; 637 ksp->rtol = 1.e-5; 638 ksp->abstol = 1.e-50; 639 ksp->divtol = 1.e4; 640 641 ksp->chknorm = -1; 642 ksp->normtype = KSP_NORM_DEFAULT; 643 ksp->rnorm = 0.0; 644 ksp->its = 0; 645 ksp->guess_zero = PETSC_TRUE; 646 ksp->calc_sings = PETSC_FALSE; 647 ksp->res_hist = PETSC_NULL; 648 ksp->res_hist_alloc = PETSC_NULL; 649 ksp->res_hist_len = 0; 650 ksp->res_hist_max = 0; 651 ksp->res_hist_reset = PETSC_TRUE; 652 ksp->numbermonitors = 0; 653 654 ierr = KSPDefaultConvergedCreate(&ctx);CHKERRQ(ierr); 655 ierr = KSPSetConvergenceTest(ksp,KSPDefaultConverged,ctx,KSPDefaultConvergedDestroy);CHKERRQ(ierr); 656 ksp->ops->buildsolution = KSPDefaultBuildSolution; 657 ksp->ops->buildresidual = KSPDefaultBuildResidual; 658 #if defined(PETSC_HAVE_AMS) 659 ((PetscObject)ksp)->bops->publish = KSPPublish_Petsc; 660 #endif 661 662 ksp->vec_sol = 0; 663 ksp->vec_rhs = 0; 664 ksp->pc = 0; 665 ksp->data = 0; 666 ksp->nwork = 0; 667 ksp->work = 0; 668 ksp->reason = KSP_CONVERGED_ITERATING; 669 ksp->setupstage = KSP_SETUP_NEW; 670 671 ierr = KSPNormSupportTableReset_Private(ksp);CHKERRQ(ierr); 672 673 *inksp = ksp; 674 PetscFunctionReturn(0); 675 } 676 677 #undef __FUNCT__ 678 #define __FUNCT__ "KSPSetType" 679 /*@C 680 KSPSetType - Builds KSP for a particular solver. 681 682 Logically Collective on KSP 683 684 Input Parameters: 685 + ksp - the Krylov space context 686 - type - a known method 687 688 Options Database Key: 689 . -ksp_type <method> - Sets the method; use -help for a list 690 of available methods (for instance, cg or gmres) 691 692 Notes: 693 See "petsc/include/petscksp.h" for available methods (for instance, 694 KSPCG or KSPGMRES). 695 696 Normally, it is best to use the KSPSetFromOptions() command and 697 then set the KSP type from the options database rather than by using 698 this routine. Using the options database provides the user with 699 maximum flexibility in evaluating the many different Krylov methods. 700 The KSPSetType() routine is provided for those situations where it 701 is necessary to set the iterative solver independently of the command 702 line or options database. This might be the case, for example, when 703 the choice of iterative solver changes during the execution of the 704 program, and the user's application is taking responsibility for 705 choosing the appropriate method. In other words, this routine is 706 not for beginners. 707 708 Level: intermediate 709 710 .keywords: KSP, set, method 711 712 .seealso: PCSetType(), KSPType 713 714 @*/ 715 PetscErrorCode KSPSetType(KSP ksp, KSPType type) 716 { 717 PetscErrorCode ierr,(*r)(KSP); 718 PetscBool match; 719 720 PetscFunctionBegin; 721 PetscValidHeaderSpecific(ksp,KSP_CLASSID,1); 722 PetscValidCharPointer(type,2); 723 724 ierr = PetscObjectTypeCompare((PetscObject)ksp,type,&match);CHKERRQ(ierr); 725 if (match) PetscFunctionReturn(0); 726 727 ierr = PetscFListFind(KSPList,((PetscObject)ksp)->comm,type,PETSC_TRUE,(void (**)(void)) &r);CHKERRQ(ierr); 728 if (!r) SETERRQ1(((PetscObject)ksp)->comm,PETSC_ERR_ARG_UNKNOWN_TYPE,"Unable to find requested KSP type %s",type); 729 /* Destroy the previous private KSP context */ 730 if (ksp->ops->destroy) { 731 ierr = (*ksp->ops->destroy)(ksp);CHKERRQ(ierr); 732 ksp->ops->destroy = PETSC_NULL; 733 } 734 /* Reinitialize function pointers in KSPOps structure */ 735 ierr = PetscMemzero(ksp->ops,sizeof(struct _KSPOps));CHKERRQ(ierr); 736 ksp->ops->buildsolution = KSPDefaultBuildSolution; 737 ksp->ops->buildresidual = KSPDefaultBuildResidual; 738 ierr = KSPNormSupportTableReset_Private(ksp);CHKERRQ(ierr); 739 /* Call the KSPCreate_XXX routine for this particular Krylov solver */ 740 ksp->setupstage = KSP_SETUP_NEW; 741 ierr = PetscObjectChangeTypeName((PetscObject)ksp,type);CHKERRQ(ierr); 742 ierr = (*r)(ksp);CHKERRQ(ierr); 743 #if defined(PETSC_HAVE_AMS) 744 if (PetscAMSPublishAll) { 745 ierr = PetscObjectAMSPublish((PetscObject)ksp);CHKERRQ(ierr); 746 } 747 #endif 748 PetscFunctionReturn(0); 749 } 750 751 #undef __FUNCT__ 752 #define __FUNCT__ "KSPRegisterDestroy" 753 /*@ 754 KSPRegisterDestroy - Frees the list of KSP methods that were 755 registered by KSPRegisterDynamic(). 756 757 Not Collective 758 759 Level: advanced 760 761 .keywords: KSP, register, destroy 762 763 .seealso: KSPRegisterDynamic(), KSPRegisterAll() 764 @*/ 765 PetscErrorCode KSPRegisterDestroy(void) 766 { 767 PetscErrorCode ierr; 768 769 PetscFunctionBegin; 770 ierr = PetscFListDestroy(&KSPList);CHKERRQ(ierr); 771 KSPRegisterAllCalled = PETSC_FALSE; 772 PetscFunctionReturn(0); 773 } 774 775 #undef __FUNCT__ 776 #define __FUNCT__ "KSPGetType" 777 /*@C 778 KSPGetType - Gets the KSP type as a string from the KSP object. 779 780 Not Collective 781 782 Input Parameter: 783 . ksp - Krylov context 784 785 Output Parameter: 786 . name - name of KSP method 787 788 Level: intermediate 789 790 .keywords: KSP, get, method, name 791 792 .seealso: KSPSetType() 793 @*/ 794 PetscErrorCode KSPGetType(KSP ksp,KSPType *type) 795 { 796 PetscFunctionBegin; 797 PetscValidHeaderSpecific(ksp,KSP_CLASSID,1); 798 PetscValidPointer(type,2); 799 *type = ((PetscObject)ksp)->type_name; 800 PetscFunctionReturn(0); 801 } 802 803 #undef __FUNCT__ 804 #define __FUNCT__ "KSPRegister" 805 /*@C 806 KSPRegister - See KSPRegisterDynamic() 807 808 Level: advanced 809 @*/ 810 PetscErrorCode KSPRegister(const char sname[],const char path[],const char name[],PetscErrorCode (*function)(KSP)) 811 { 812 PetscErrorCode ierr; 813 char fullname[PETSC_MAX_PATH_LEN]; 814 815 PetscFunctionBegin; 816 ierr = PetscFListConcat(path,name,fullname);CHKERRQ(ierr); 817 ierr = PetscFListAdd(&KSPList,sname,fullname,(void (*)(void))function);CHKERRQ(ierr); 818 PetscFunctionReturn(0); 819 } 820 821 #undef __FUNCT__ 822 #define __FUNCT__ "KSPSetNullSpace" 823 /*@ 824 KSPSetNullSpace - Sets the null space of the operator 825 826 Logically Collective on KSP 827 828 Input Parameters: 829 + ksp - the Krylov space object 830 - nullsp - the null space of the operator 831 832 Level: advanced 833 834 .seealso: KSPSetOperators(), MatNullSpaceCreate(), KSPGetNullSpace(), MatSetNullSpace() 835 @*/ 836 PetscErrorCode KSPSetNullSpace(KSP ksp,MatNullSpace nullsp) 837 { 838 PetscErrorCode ierr; 839 840 PetscFunctionBegin; 841 PetscValidHeaderSpecific(ksp,KSP_CLASSID,1); 842 PetscValidHeaderSpecific(nullsp,MAT_NULLSPACE_CLASSID,2); 843 ierr = PetscObjectReference((PetscObject)nullsp);CHKERRQ(ierr); 844 if (ksp->nullsp) { ierr = MatNullSpaceDestroy(&ksp->nullsp);CHKERRQ(ierr); } 845 ksp->nullsp = nullsp; 846 PetscFunctionReturn(0); 847 } 848 849 #undef __FUNCT__ 850 #define __FUNCT__ "KSPGetNullSpace" 851 /*@ 852 KSPGetNullSpace - Gets the null space of the operator 853 854 Not Collective 855 856 Input Parameters: 857 + ksp - the Krylov space object 858 - nullsp - the null space of the operator 859 860 Level: advanced 861 862 .seealso: KSPSetOperators(), MatNullSpaceCreate(), KSPSetNullSpace() 863 @*/ 864 PetscErrorCode KSPGetNullSpace(KSP ksp,MatNullSpace *nullsp) 865 { 866 PetscFunctionBegin; 867 PetscValidHeaderSpecific(ksp,KSP_CLASSID,1); 868 PetscValidPointer(nullsp,2); 869 *nullsp = ksp->nullsp; 870 PetscFunctionReturn(0); 871 } 872 873