1 /* 2 Defines the interface functions for the Krylov subspace accelerators. 3 */ 4 #ifndef __PETSCKSP_H 5 #define __PETSCKSP_H 6 #include <petscpc.h> 7 8 PETSC_EXTERN PetscErrorCode KSPInitializePackage(void); 9 10 /*S 11 KSP - Abstract PETSc object that manages all Krylov methods. This is the object that manages the 12 linear solves in PETSc (even those such as direct solvers that do no use Krylov accelerators). 13 14 Level: beginner 15 16 Concepts: Krylov methods 17 18 Notes: 19 When a direct solver is used but no Krylov solver is used the KSP object is still used by with a 20 KSPType of KSPPREONLY (meaning application of the preconditioner is only used as the linear solver). 21 22 .seealso: KSPCreate(), KSPSetType(), KSPType, SNES, TS, PC, KSP, KSPDestroy(), KSPCG, KSPGMRES 23 S*/ 24 typedef struct _p_KSP* KSP; 25 26 /*J 27 KSPType - String with the name of a PETSc Krylov method. 28 29 Level: beginner 30 31 .seealso: KSPSetType(), KSP, KSPRegister(), KSPCreate(), KSPSetFromOptions() 32 J*/ 33 typedef const char* KSPType; 34 #define KSPRICHARDSON "richardson" 35 #define KSPCHEBYSHEV "chebyshev" 36 #define KSPCG "cg" 37 #define KSPGROPPCG "groppcg" 38 #define KSPPIPECG "pipecg" 39 #define KSPPIPECGRR "pipecgrr" 40 #define KSPPIPELCG "pipelcg" 41 #define KSPCGNE "cgne" 42 #define KSPCGNASH "nash" 43 #define KSPCGSTCG "stcg" 44 #define KSPCGGLTR "gltr" 45 #define KSPFCG "fcg" 46 #define KSPPIPEFCG "pipefcg" 47 #define KSPGMRES "gmres" 48 #define KSPPIPEFGMRES "pipefgmres" 49 #define KSPFGMRES "fgmres" 50 #define KSPLGMRES "lgmres" 51 #define KSPDGMRES "dgmres" 52 #define KSPPGMRES "pgmres" 53 #define KSPTCQMR "tcqmr" 54 #define KSPBCGS "bcgs" 55 #define KSPIBCGS "ibcgs" 56 #define KSPFBCGS "fbcgs" 57 #define KSPFBCGSR "fbcgsr" 58 #define KSPBCGSL "bcgsl" 59 #define KSPPIPEBCGS "pipebcgs" 60 #define KSPCGS "cgs" 61 #define KSPTFQMR "tfqmr" 62 #define KSPCR "cr" 63 #define KSPPIPECR "pipecr" 64 #define KSPLSQR "lsqr" 65 #define KSPPREONLY "preonly" 66 #define KSPQCG "qcg" 67 #define KSPBICG "bicg" 68 #define KSPMINRES "minres" 69 #define KSPSYMMLQ "symmlq" 70 #define KSPLCD "lcd" 71 #define KSPPYTHON "python" 72 #define KSPGCR "gcr" 73 #define KSPPIPEGCR "pipegcr" 74 #define KSPTSIRM "tsirm" 75 #define KSPCGLS "cgls" 76 #define KSPFETIDP "fetidp" 77 78 /* Logging support */ 79 PETSC_EXTERN PetscClassId KSP_CLASSID; 80 PETSC_EXTERN PetscClassId KSPGUESS_CLASSID; 81 PETSC_EXTERN PetscClassId DMKSP_CLASSID; 82 83 PETSC_EXTERN PetscErrorCode KSPCreate(MPI_Comm,KSP*); 84 PETSC_EXTERN PetscErrorCode KSPSetType(KSP,KSPType); 85 PETSC_EXTERN PetscErrorCode KSPGetType(KSP,KSPType*); 86 PETSC_EXTERN PetscErrorCode KSPSetUp(KSP); 87 PETSC_EXTERN PetscErrorCode KSPSetUpOnBlocks(KSP); 88 PETSC_EXTERN PetscErrorCode KSPSolve(KSP,Vec,Vec); 89 PETSC_EXTERN PetscErrorCode KSPSolveTranspose(KSP,Vec,Vec); 90 PETSC_EXTERN PetscErrorCode KSPReset(KSP); 91 PETSC_EXTERN PetscErrorCode KSPDestroy(KSP*); 92 PETSC_EXTERN PetscErrorCode KSPSetReusePreconditioner(KSP,PetscBool); 93 PETSC_EXTERN PetscErrorCode KSPSetSkipPCSetFromOptions(KSP,PetscBool); 94 95 PETSC_EXTERN PetscFunctionList KSPList; 96 PETSC_EXTERN PetscFunctionList KSPGuessList; 97 PETSC_EXTERN PetscErrorCode KSPRegister(const char[],PetscErrorCode (*)(KSP)); 98 99 PETSC_EXTERN PetscErrorCode KSPSetPCSide(KSP,PCSide); 100 PETSC_EXTERN PetscErrorCode KSPGetPCSide(KSP,PCSide*); 101 PETSC_EXTERN PetscErrorCode KSPSetTolerances(KSP,PetscReal,PetscReal,PetscReal,PetscInt); 102 PETSC_EXTERN PetscErrorCode KSPGetTolerances(KSP,PetscReal*,PetscReal*,PetscReal*,PetscInt*); 103 PETSC_EXTERN PetscErrorCode KSPSetInitialGuessNonzero(KSP,PetscBool); 104 PETSC_EXTERN PetscErrorCode KSPGetInitialGuessNonzero(KSP,PetscBool*); 105 PETSC_EXTERN PetscErrorCode KSPSetErrorIfNotConverged(KSP,PetscBool); 106 PETSC_EXTERN PetscErrorCode KSPGetErrorIfNotConverged(KSP,PetscBool*); 107 PETSC_EXTERN PetscErrorCode KSPSetComputeEigenvalues(KSP,PetscBool); 108 PETSC_EXTERN PetscErrorCode KSPSetComputeRitz(KSP,PetscBool); 109 PETSC_EXTERN PetscErrorCode KSPGetComputeEigenvalues(KSP,PetscBool*); 110 PETSC_EXTERN PetscErrorCode KSPSetComputeSingularValues(KSP,PetscBool); 111 PETSC_EXTERN PetscErrorCode KSPGetComputeSingularValues(KSP,PetscBool*); 112 PETSC_EXTERN PetscErrorCode KSPGetRhs(KSP,Vec*); 113 PETSC_EXTERN PetscErrorCode KSPGetSolution(KSP,Vec*); 114 PETSC_EXTERN PetscErrorCode KSPGetResidualNorm(KSP,PetscReal*); 115 PETSC_EXTERN PetscErrorCode KSPGetIterationNumber(KSP,PetscInt*); 116 PETSC_EXTERN PetscErrorCode KSPGetTotalIterations(KSP,PetscInt*); 117 PETSC_EXTERN PetscErrorCode KSPCreateVecs(KSP,PetscInt,Vec**,PetscInt,Vec**); 118 PETSC_DEPRECATED("Use KSPCreateVecs()") PETSC_STATIC_INLINE PetscErrorCode KSPGetVecs(KSP ksp,PetscInt n,Vec **x,PetscInt m,Vec **y) {return KSPCreateVecs(ksp,n,x,m,y);} 119 120 PETSC_EXTERN PetscErrorCode KSPSetPreSolve(KSP,PetscErrorCode (*)(KSP,Vec,Vec,void*),void*); 121 PETSC_EXTERN PetscErrorCode KSPSetPostSolve(KSP,PetscErrorCode (*)(KSP,Vec,Vec,void*),void*); 122 123 PETSC_EXTERN PetscErrorCode KSPSetPC(KSP,PC); 124 PETSC_EXTERN PetscErrorCode KSPGetPC(KSP,PC*); 125 126 PETSC_EXTERN PetscErrorCode KSPMonitor(KSP,PetscInt,PetscReal); 127 PETSC_EXTERN PetscErrorCode KSPMonitorSet(KSP,PetscErrorCode (*)(KSP,PetscInt,PetscReal,void*),void*,PetscErrorCode (*)(void**)); 128 PETSC_EXTERN PetscErrorCode KSPMonitorCancel(KSP); 129 PETSC_EXTERN PetscErrorCode KSPGetMonitorContext(KSP,void**); 130 PETSC_EXTERN PetscErrorCode KSPGetResidualHistory(KSP,PetscReal*[],PetscInt*); 131 PETSC_EXTERN PetscErrorCode KSPSetResidualHistory(KSP,PetscReal[],PetscInt,PetscBool ); 132 133 PETSC_EXTERN PetscErrorCode KSPBuildSolutionDefault(KSP,Vec,Vec*); 134 PETSC_EXTERN PetscErrorCode KSPBuildResidualDefault(KSP,Vec,Vec,Vec*); 135 PETSC_EXTERN PetscErrorCode KSPDestroyDefault(KSP); 136 PETSC_EXTERN PetscErrorCode KSPSetWorkVecs(KSP,PetscInt); 137 138 PETSC_EXTERN PetscErrorCode PCKSPGetKSP(PC,KSP*); 139 PETSC_EXTERN PetscErrorCode PCBJacobiGetSubKSP(PC,PetscInt*,PetscInt*,KSP*[]); 140 PETSC_EXTERN PetscErrorCode PCASMGetSubKSP(PC,PetscInt*,PetscInt*,KSP*[]); 141 PETSC_EXTERN PetscErrorCode PCGASMGetSubKSP(PC,PetscInt*,PetscInt*,KSP*[]); 142 PETSC_EXTERN PetscErrorCode PCFieldSplitGetSubKSP(PC,PetscInt*,KSP*[]); 143 PETSC_EXTERN PetscErrorCode PCMGGetSmoother(PC,PetscInt,KSP*); 144 PETSC_EXTERN PetscErrorCode PCMGGetSmootherDown(PC,PetscInt,KSP*); 145 PETSC_EXTERN PetscErrorCode PCMGGetSmootherUp(PC,PetscInt,KSP*); 146 PETSC_EXTERN PetscErrorCode PCMGGetCoarseSolve(PC,KSP*); 147 PETSC_EXTERN PetscErrorCode PCGalerkinGetKSP(PC,KSP*); 148 149 PETSC_EXTERN PetscErrorCode KSPBuildSolution(KSP,Vec,Vec*); 150 PETSC_EXTERN PetscErrorCode KSPBuildResidual(KSP,Vec,Vec,Vec*); 151 152 PETSC_EXTERN PetscErrorCode KSPRichardsonSetScale(KSP,PetscReal); 153 PETSC_EXTERN PetscErrorCode KSPRichardsonSetSelfScale(KSP,PetscBool ); 154 PETSC_EXTERN PetscErrorCode KSPChebyshevSetEigenvalues(KSP,PetscReal,PetscReal); 155 PETSC_EXTERN PetscErrorCode KSPChebyshevEstEigSet(KSP,PetscReal,PetscReal,PetscReal,PetscReal); 156 PETSC_EXTERN PetscErrorCode KSPChebyshevEstEigSetUseNoisy(KSP,PetscBool); 157 PETSC_EXTERN PetscErrorCode KSPChebyshevEstEigGetKSP(KSP,KSP*); 158 PETSC_EXTERN PetscErrorCode KSPComputeExtremeSingularValues(KSP,PetscReal*,PetscReal*); 159 PETSC_EXTERN PetscErrorCode KSPComputeEigenvalues(KSP,PetscInt,PetscReal[],PetscReal[],PetscInt*); 160 PETSC_EXTERN PetscErrorCode KSPComputeEigenvaluesExplicitly(KSP,PetscInt,PetscReal[],PetscReal[]); 161 PETSC_EXTERN PetscErrorCode KSPComputeRitz(KSP,PetscBool,PetscBool,PetscInt*,Vec[],PetscReal[],PetscReal[]); 162 163 /*E 164 165 KSPFCDTruncationType - Define how stored directions are used to orthogonalize in flexible conjugate directions (FCD) methods 166 167 KSP_FCD_TRUNC_TYPE_STANDARD uses all (up to mmax) stored directions 168 KSP_FCD_TRUNC_TYPE_NOTAY uses the last max(1,mod(i,mmax)) stored directions at iteration i=0,1.. 169 170 Level: intermediate 171 .seealso : KSPFCG,KSPPIPEFCG,KSPPIPEGCR,KSPFCGSetTruncationType(),KSPFCGGetTruncationType() 172 173 E*/ 174 typedef enum {KSP_FCD_TRUNC_TYPE_STANDARD,KSP_FCD_TRUNC_TYPE_NOTAY} KSPFCDTruncationType; 175 PETSC_EXTERN const char *const KSPFCDTruncationTypes[]; 176 177 PETSC_EXTERN PetscErrorCode KSPFCGSetMmax(KSP,PetscInt); 178 PETSC_EXTERN PetscErrorCode KSPFCGGetMmax(KSP,PetscInt*); 179 PETSC_EXTERN PetscErrorCode KSPFCGSetNprealloc(KSP,PetscInt); 180 PETSC_EXTERN PetscErrorCode KSPFCGGetNprealloc(KSP,PetscInt*); 181 PETSC_EXTERN PetscErrorCode KSPFCGSetTruncationType(KSP,KSPFCDTruncationType); 182 PETSC_EXTERN PetscErrorCode KSPFCGGetTruncationType(KSP,KSPFCDTruncationType*); 183 184 PETSC_EXTERN PetscErrorCode KSPPIPEFCGSetMmax(KSP,PetscInt); 185 PETSC_EXTERN PetscErrorCode KSPPIPEFCGGetMmax(KSP,PetscInt*); 186 PETSC_EXTERN PetscErrorCode KSPPIPEFCGSetNprealloc(KSP,PetscInt); 187 PETSC_EXTERN PetscErrorCode KSPPIPEFCGGetNprealloc(KSP,PetscInt*); 188 PETSC_EXTERN PetscErrorCode KSPPIPEFCGSetTruncationType(KSP,KSPFCDTruncationType); 189 PETSC_EXTERN PetscErrorCode KSPPIPEFCGGetTruncationType(KSP,KSPFCDTruncationType*); 190 191 PETSC_EXTERN PetscErrorCode KSPPIPEGCRSetMmax(KSP,PetscInt); 192 PETSC_EXTERN PetscErrorCode KSPPIPEGCRGetMmax(KSP,PetscInt*); 193 PETSC_EXTERN PetscErrorCode KSPPIPEGCRSetNprealloc(KSP,PetscInt); 194 PETSC_EXTERN PetscErrorCode KSPPIPEGCRGetNprealloc(KSP,PetscInt*); 195 PETSC_EXTERN PetscErrorCode KSPPIPEGCRSetTruncationType(KSP,KSPFCDTruncationType); 196 PETSC_EXTERN PetscErrorCode KSPPIPEGCRGetTruncationType(KSP,KSPFCDTruncationType*); 197 PETSC_EXTERN PetscErrorCode KSPPIPEGCRSetUnrollW(KSP,PetscBool); 198 PETSC_EXTERN PetscErrorCode KSPPIPEGCRGetUnrollW(KSP,PetscBool*); 199 200 PETSC_EXTERN PetscErrorCode KSPGMRESSetRestart(KSP, PetscInt); 201 PETSC_EXTERN PetscErrorCode KSPGMRESGetRestart(KSP, PetscInt*); 202 PETSC_EXTERN PetscErrorCode KSPGMRESSetHapTol(KSP,PetscReal); 203 204 PETSC_EXTERN PetscErrorCode KSPGMRESSetPreAllocateVectors(KSP); 205 PETSC_EXTERN PetscErrorCode KSPGMRESSetOrthogonalization(KSP,PetscErrorCode (*)(KSP,PetscInt)); 206 PETSC_EXTERN PetscErrorCode KSPGMRESGetOrthogonalization(KSP,PetscErrorCode (**)(KSP,PetscInt)); 207 PETSC_EXTERN PetscErrorCode KSPGMRESModifiedGramSchmidtOrthogonalization(KSP,PetscInt); 208 PETSC_EXTERN PetscErrorCode KSPGMRESClassicalGramSchmidtOrthogonalization(KSP,PetscInt); 209 210 PETSC_EXTERN PetscErrorCode KSPLGMRESSetAugDim(KSP,PetscInt); 211 PETSC_EXTERN PetscErrorCode KSPLGMRESSetConstant(KSP); 212 213 PETSC_EXTERN PetscErrorCode KSPPIPEFGMRESSetShift(KSP,PetscScalar); 214 215 PETSC_EXTERN PetscErrorCode KSPGCRSetRestart(KSP,PetscInt); 216 PETSC_EXTERN PetscErrorCode KSPGCRGetRestart(KSP,PetscInt*); 217 PETSC_EXTERN PetscErrorCode KSPGCRSetModifyPC(KSP,PetscErrorCode (*)(KSP,PetscInt,PetscReal,void*),void*,PetscErrorCode(*)(void*)); 218 219 PETSC_EXTERN PetscErrorCode KSPFETIDPGetInnerBDDC(KSP,PC*); 220 PETSC_EXTERN PetscErrorCode KSPFETIDPSetInnerBDDC(KSP,PC); 221 PETSC_EXTERN PetscErrorCode KSPFETIDPGetInnerKSP(KSP,KSP*); 222 PETSC_EXTERN PetscErrorCode KSPFETIDPSetPressureOperator(KSP,Mat); 223 /*E 224 KSPGMRESCGSRefinementType - How the classical (unmodified) Gram-Schmidt is performed. 225 226 Level: advanced 227 228 .seealso: KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESSetOrthogonalization(), KSPGMRESGetOrthogonalization(), 229 KSPGMRESSetCGSRefinementType(), KSPGMRESGetCGSRefinementType(), KSPGMRESModifiedGramSchmidtOrthogonalization() 230 231 E*/ 232 typedef enum {KSP_GMRES_CGS_REFINE_NEVER, KSP_GMRES_CGS_REFINE_IFNEEDED, KSP_GMRES_CGS_REFINE_ALWAYS} KSPGMRESCGSRefinementType; 233 PETSC_EXTERN const char *const KSPGMRESCGSRefinementTypes[]; 234 /*MC 235 KSP_GMRES_CGS_REFINE_NEVER - Do the classical (unmodified) Gram-Schmidt process 236 237 Level: advanced 238 239 Note: Possible unstable, but the fastest to compute 240 241 .seealso: KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESSetOrthogonalization(), KSPGMRESGetOrthogonalization(), 242 KSPGMRESSetCGSRefinementType(), KSPGMRESGetCGSRefinementType(), KSP_GMRES_CGS_REFINE_IFNEEDED, KSP_GMRES_CGS_REFINE_ALWAYS, 243 KSPGMRESModifiedGramSchmidtOrthogonalization() 244 M*/ 245 246 /*MC 247 KSP_GMRES_CGS_REFINE_IFNEEDED - Do the classical (unmodified) Gram-Schmidt process and one step of 248 iterative refinement if an estimate of the orthogonality of the resulting vectors indicates 249 poor orthogonality. 250 251 Level: advanced 252 253 Note: This is slower than KSP_GMRES_CGS_REFINE_NEVER because it requires an extra norm computation to 254 estimate the orthogonality but is more stable. 255 256 .seealso: KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESSetOrthogonalization(), KSPGMRESGetOrthogonalization(), 257 KSPGMRESSetCGSRefinementType(), KSPGMRESGetCGSRefinementType(), KSP_GMRES_CGS_REFINE_NEVER, KSP_GMRES_CGS_REFINE_ALWAYS, 258 KSPGMRESModifiedGramSchmidtOrthogonalization() 259 M*/ 260 261 /*MC 262 KSP_GMRES_CGS_REFINE_NEVER - Do two steps of the classical (unmodified) Gram-Schmidt process. 263 264 Level: advanced 265 266 Note: This is roughly twice the cost of KSP_GMRES_CGS_REFINE_NEVER because it performs the process twice 267 but it saves the extra norm calculation needed by KSP_GMRES_CGS_REFINE_IFNEEDED. 268 269 You should only use this if you absolutely know that the iterative refinement is needed. 270 271 .seealso: KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESSetOrthogonalization(), KSPGMRESGetOrthogonalization(), 272 KSPGMRESSetCGSRefinementType(), KSPGMRESGetCGSRefinementType(), KSP_GMRES_CGS_REFINE_IFNEEDED, KSP_GMRES_CGS_REFINE_ALWAYS, 273 KSPGMRESModifiedGramSchmidtOrthogonalization() 274 M*/ 275 276 PETSC_EXTERN PetscErrorCode KSPGMRESSetCGSRefinementType(KSP,KSPGMRESCGSRefinementType); 277 PETSC_EXTERN PetscErrorCode KSPGMRESGetCGSRefinementType(KSP,KSPGMRESCGSRefinementType*); 278 279 PETSC_EXTERN PetscErrorCode KSPFGMRESModifyPCNoChange(KSP,PetscInt,PetscInt,PetscReal,void*); 280 PETSC_EXTERN PetscErrorCode KSPFGMRESModifyPCKSP(KSP,PetscInt,PetscInt,PetscReal,void*); 281 PETSC_EXTERN PetscErrorCode KSPFGMRESSetModifyPC(KSP,PetscErrorCode (*)(KSP,PetscInt,PetscInt,PetscReal,void*),void*,PetscErrorCode(*)(void*)); 282 283 PETSC_EXTERN PetscErrorCode KSPQCGSetTrustRegionRadius(KSP,PetscReal); 284 PETSC_EXTERN PetscErrorCode KSPQCGGetQuadratic(KSP,PetscReal*); 285 PETSC_EXTERN PetscErrorCode KSPQCGGetTrialStepNorm(KSP,PetscReal*); 286 287 PETSC_EXTERN PetscErrorCode KSPBCGSLSetXRes(KSP,PetscReal); 288 PETSC_EXTERN PetscErrorCode KSPBCGSLSetPol(KSP,PetscBool ); 289 PETSC_EXTERN PetscErrorCode KSPBCGSLSetEll(KSP,PetscInt); 290 PETSC_EXTERN PetscErrorCode KSPBCGSLSetUsePseudoinverse(KSP,PetscBool); 291 292 PETSC_EXTERN PetscErrorCode KSPSetFromOptions(KSP); 293 PETSC_EXTERN PetscErrorCode KSPAddOptionsChecker(PetscErrorCode (*)(KSP)); 294 295 PETSC_EXTERN PetscErrorCode KSPMonitorSingularValue(KSP,PetscInt,PetscReal,PetscViewerAndFormat*); 296 PETSC_EXTERN PetscErrorCode KSPMonitorDefault(KSP,PetscInt,PetscReal,PetscViewerAndFormat*); 297 PETSC_EXTERN PetscErrorCode KSPLSQRMonitorDefault(KSP,PetscInt,PetscReal,PetscViewerAndFormat*); 298 PETSC_EXTERN PetscErrorCode KSPMonitorRange(KSP,PetscInt,PetscReal,PetscViewerAndFormat*); 299 PETSC_EXTERN PetscErrorCode KSPMonitorDynamicTolerance(KSP ksp,PetscInt its,PetscReal fnorm,void *dummy); 300 PETSC_EXTERN PetscErrorCode KSPMonitorDynamicToleranceDestroy(void **dummy); 301 PETSC_EXTERN PetscErrorCode KSPMonitorTrueResidualNorm(KSP,PetscInt,PetscReal,PetscViewerAndFormat*); 302 PETSC_EXTERN PetscErrorCode KSPMonitorTrueResidualMaxNorm(KSP,PetscInt,PetscReal,PetscViewerAndFormat*); 303 PETSC_EXTERN PetscErrorCode KSPMonitorDefaultShort(KSP,PetscInt,PetscReal,PetscViewerAndFormat*); 304 PETSC_EXTERN PetscErrorCode KSPMonitorSolution(KSP,PetscInt,PetscReal,PetscViewerAndFormat*); 305 PETSC_EXTERN PetscErrorCode KSPMonitorSAWs(KSP,PetscInt,PetscReal,void*); 306 PETSC_EXTERN PetscErrorCode KSPMonitorSAWsCreate(KSP,void**); 307 PETSC_EXTERN PetscErrorCode KSPMonitorSAWsDestroy(void**); 308 PETSC_EXTERN PetscErrorCode KSPGMRESMonitorKrylov(KSP,PetscInt,PetscReal,void*); 309 PETSC_EXTERN PetscErrorCode KSPMonitorSetFromOptions(KSP,const char[],const char[],const char [],PetscErrorCode (*)(KSP,PetscInt,PetscReal,PetscViewerAndFormat*)); 310 311 PETSC_EXTERN PetscErrorCode KSPUnwindPreconditioner(KSP,Vec,Vec); 312 PETSC_EXTERN PetscErrorCode KSPInitialResidual(KSP,Vec,Vec,Vec,Vec,Vec); 313 314 PETSC_EXTERN PetscErrorCode KSPSetOperators(KSP,Mat,Mat); 315 PETSC_EXTERN PetscErrorCode KSPGetOperators(KSP,Mat*,Mat*); 316 PETSC_EXTERN PetscErrorCode KSPGetOperatorsSet(KSP,PetscBool*,PetscBool*); 317 PETSC_EXTERN PetscErrorCode KSPSetOptionsPrefix(KSP,const char[]); 318 PETSC_EXTERN PetscErrorCode KSPAppendOptionsPrefix(KSP,const char[]); 319 PETSC_EXTERN PetscErrorCode KSPGetOptionsPrefix(KSP,const char*[]); 320 PETSC_EXTERN PetscErrorCode KSPSetTabLevel(KSP,PetscInt); 321 PETSC_EXTERN PetscErrorCode KSPGetTabLevel(KSP,PetscInt*); 322 323 PETSC_EXTERN PetscErrorCode KSPSetDiagonalScale(KSP,PetscBool ); 324 PETSC_EXTERN PetscErrorCode KSPGetDiagonalScale(KSP,PetscBool*); 325 PETSC_EXTERN PetscErrorCode KSPSetDiagonalScaleFix(KSP,PetscBool ); 326 PETSC_EXTERN PetscErrorCode KSPGetDiagonalScaleFix(KSP,PetscBool*); 327 328 PETSC_EXTERN PetscErrorCode KSPView(KSP,PetscViewer); 329 PETSC_EXTERN PetscErrorCode KSPLoad(KSP,PetscViewer); 330 PETSC_STATIC_INLINE PetscErrorCode KSPViewFromOptions(KSP A,PetscObject obj,const char name[]) {return PetscObjectViewFromOptions((PetscObject)A,obj,name);} 331 PETSC_EXTERN PetscErrorCode KSPReasonView(KSP,PetscViewer); 332 PETSC_EXTERN PetscErrorCode KSPReasonViewFromOptions(KSP); 333 334 #define KSP_FILE_CLASSID 1211223 335 336 PETSC_EXTERN PetscErrorCode KSPLSQRSetStandardErrorVec(KSP,Vec); 337 PETSC_EXTERN PetscErrorCode KSPLSQRGetStandardErrorVec(KSP,Vec*); 338 PETSC_EXTERN PetscErrorCode KSPLSQRGetArnorm(KSP,PetscReal*,PetscReal*,PetscReal*); 339 340 PETSC_EXTERN PetscErrorCode PCRedundantGetKSP(PC,KSP*); 341 PETSC_EXTERN PetscErrorCode PCRedistributeGetKSP(PC,KSP*); 342 PETSC_EXTERN PetscErrorCode PCTelescopeGetKSP(PC,KSP*); 343 344 /*E 345 KSPNormType - Norm that is passed in the Krylov convergence 346 test routines. 347 348 Level: advanced 349 350 Each solver only supports a subset of these and some may support different ones 351 depending on left or right preconditioning, see KSPSetPCSide() 352 353 Notes: 354 this must match petsc/finclude/petscksp.h 355 356 .seealso: KSPSolve(), KSPGetConvergedReason(), KSPSetNormType(), 357 KSPSetConvergenceTest(), KSPSetPCSide() 358 E*/ 359 typedef enum {KSP_NORM_DEFAULT = -1,KSP_NORM_NONE = 0,KSP_NORM_PRECONDITIONED = 1,KSP_NORM_UNPRECONDITIONED = 2,KSP_NORM_NATURAL = 3} KSPNormType; 360 #define KSP_NORM_MAX (KSP_NORM_NATURAL + 1) 361 PETSC_EXTERN const char *const*const KSPNormTypes; 362 363 /*MC 364 KSP_NORM_NONE - Do not compute a norm during the Krylov process. This will 365 possibly save some computation but means the convergence test cannot 366 be based on a norm of a residual etc. 367 368 Level: advanced 369 370 Note: Some Krylov methods need to compute a residual norm (such as GMRES) and then this option is ignored 371 372 .seealso: KSPNormType, KSPSetNormType(), KSP_NORM_PRECONDITIONED, KSP_NORM_UNPRECONDITIONED, KSP_NORM_NATURAL 373 M*/ 374 375 /*MC 376 KSP_NORM_PRECONDITIONED - Compute the norm of the preconditioned residual B*(b - A*x), if left preconditioning, and pass that to the 377 convergence test routine. 378 379 Level: advanced 380 381 .seealso: KSPNormType, KSPSetNormType(), KSP_NORM_NONE, KSP_NORM_UNPRECONDITIONED, KSP_NORM_NATURAL, KSPSetConvergenceTest() 382 M*/ 383 384 /*MC 385 KSP_NORM_UNPRECONDITIONED - Compute the norm of the true residual (b - A*x) and pass that to the 386 convergence test routine. 387 388 Level: advanced 389 390 .seealso: KSPNormType, KSPSetNormType(), KSP_NORM_NONE, KSP_NORM_PRECONDITIONED, KSP_NORM_NATURAL, KSPSetConvergenceTest() 391 M*/ 392 393 /*MC 394 KSP_NORM_NATURAL - Compute the 'natural norm' of residual sqrt((b - A*x)*B*(b - A*x)) and pass that to the 395 convergence test routine. This is only supported by KSPCG, KSPCR, KSPCGNE, KSPCGS, KSPFCG, KSPPIPEFCG, KSPPIPEGCR 396 397 Level: advanced 398 399 .seealso: KSPNormType, KSPSetNormType(), KSP_NORM_NONE, KSP_NORM_PRECONDITIONED, KSP_NORM_UNPRECONDITIONED, KSPSetConvergenceTest() 400 M*/ 401 402 PETSC_EXTERN PetscErrorCode KSPSetNormType(KSP,KSPNormType); 403 PETSC_EXTERN PetscErrorCode KSPGetNormType(KSP,KSPNormType*); 404 PETSC_EXTERN PetscErrorCode KSPSetSupportedNorm(KSP ksp,KSPNormType,PCSide,PetscInt); 405 PETSC_EXTERN PetscErrorCode KSPSetCheckNormIteration(KSP,PetscInt); 406 PETSC_EXTERN PetscErrorCode KSPSetLagNorm(KSP,PetscBool); 407 408 /*E 409 KSPConvergedReason - reason a Krylov method was said to have converged or diverged 410 411 Level: beginner 412 413 Notes: 414 See KSPGetConvergedReason() for explanation of each value 415 416 Developer Notes: 417 this must match petsc/finclude/petscksp.h 418 419 The string versions of these are KSPConvergedReasons; if you change 420 any of the values here also change them that array of names. 421 422 .seealso: KSPSolve(), KSPGetConvergedReason(), KSPSetTolerances() 423 E*/ 424 typedef enum {/* converged */ 425 KSP_CONVERGED_RTOL_NORMAL = 1, 426 KSP_CONVERGED_ATOL_NORMAL = 9, 427 KSP_CONVERGED_RTOL = 2, 428 KSP_CONVERGED_ATOL = 3, 429 KSP_CONVERGED_ITS = 4, 430 KSP_CONVERGED_CG_NEG_CURVE = 5, 431 KSP_CONVERGED_CG_CONSTRAINED = 6, 432 KSP_CONVERGED_STEP_LENGTH = 7, 433 KSP_CONVERGED_HAPPY_BREAKDOWN = 8, 434 /* diverged */ 435 KSP_DIVERGED_NULL = -2, 436 KSP_DIVERGED_ITS = -3, 437 KSP_DIVERGED_DTOL = -4, 438 KSP_DIVERGED_BREAKDOWN = -5, 439 KSP_DIVERGED_BREAKDOWN_BICG = -6, 440 KSP_DIVERGED_NONSYMMETRIC = -7, 441 KSP_DIVERGED_INDEFINITE_PC = -8, 442 KSP_DIVERGED_NANORINF = -9, 443 KSP_DIVERGED_INDEFINITE_MAT = -10, 444 KSP_DIVERGED_PCSETUP_FAILED = -11, 445 446 KSP_CONVERGED_ITERATING = 0} KSPConvergedReason; 447 PETSC_EXTERN const char *const*KSPConvergedReasons; 448 449 /*MC 450 KSP_CONVERGED_RTOL - norm(r) <= rtol*norm(b) or rtol*norm(b - A*x_0) if KSPConvergedDefaultSetUIRNorm() was called 451 452 Level: beginner 453 454 See KSPNormType and KSPSetNormType() for possible norms that may be used. By default 455 for left preconditioning it is the 2-norm of the preconditioned residual, and the 456 2-norm of the residual for right preconditioning 457 458 See also KSP_CONVERGED_ATOL which may apply before this tolerance. 459 460 .seealso: KSP_CONVERGED_ATOL, KSP_DIVERGED_DTOL, KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances() 461 462 M*/ 463 464 /*MC 465 KSP_CONVERGED_ATOL - norm(r) <= atol 466 467 Level: beginner 468 469 See KSPNormType and KSPSetNormType() for possible norms that may be used. By default 470 for left preconditioning it is the 2-norm of the preconditioned residual, and the 471 2-norm of the residual for right preconditioning 472 473 See also KSP_CONVERGED_RTOL which may apply before this tolerance. 474 475 .seealso: KSP_CONVERGED_RTOL, KSP_DIVERGED_DTOL, KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances() 476 477 M*/ 478 479 /*MC 480 KSP_DIVERGED_DTOL - norm(r) >= dtol*norm(b) 481 482 Level: beginner 483 484 See KSPNormType and KSPSetNormType() for possible norms that may be used. By default 485 for left preconditioning it is the 2-norm of the preconditioned residual, and the 486 2-norm of the residual for right preconditioning 487 488 Level: beginner 489 490 .seealso: KSP_CONVERGED_ATOL, KSP_DIVERGED_RTOL, KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances() 491 492 M*/ 493 494 /*MC 495 KSP_DIVERGED_ITS - Ran out of iterations before any convergence criteria was 496 reached 497 498 Level: beginner 499 500 .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances() 501 502 M*/ 503 504 /*MC 505 KSP_CONVERGED_ITS - Used by the KSPPREONLY solver after the single iteration of 506 the preconditioner is applied. Also used when the KSPConvergedSkip() convergence 507 test routine is set in KSP. 508 509 Level: beginner 510 511 .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances() 512 513 M*/ 514 515 /*MC 516 KSP_DIVERGED_BREAKDOWN - A breakdown in the Krylov method was detected so the 517 method could not continue to enlarge the Krylov space. Could be due to a singlular matrix or 518 preconditioner. 519 520 Level: beginner 521 522 .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances() 523 524 M*/ 525 526 /*MC 527 KSP_DIVERGED_BREAKDOWN_BICG - A breakdown in the KSPBICG method was detected so the 528 method could not continue to enlarge the Krylov space. 529 530 Level: beginner 531 532 .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances() 533 534 M*/ 535 536 /*MC 537 KSP_DIVERGED_NONSYMMETRIC - It appears the operator or preconditioner is not 538 symmetric and this Krylov method (KSPCG, KSPMINRES, KSPCR) requires symmetry 539 540 Level: beginner 541 542 .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances() 543 544 M*/ 545 546 /*MC 547 KSP_DIVERGED_INDEFINITE_PC - It appears the preconditioner is indefinite (has both 548 positive and negative eigenvalues) and this Krylov method (KSPCG) requires it to 549 be positive definite 550 551 Level: beginner 552 553 Notes: 554 This can happen with the PCICC preconditioner, use -pc_factor_shift_positive_definite to force 555 the PCICC preconditioner to generate a positive definite preconditioner 556 557 .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances() 558 559 M*/ 560 561 /*MC 562 KSP_DIVERGED_PCSETUP_FAILED - It was not possible to build the requested preconditioner. This is usually due to a 563 zero pivot in a factorization. It can also result from a failure in a subpreconditioner inside a nested preconditioner 564 such as PCFIELDSPLIT. 565 566 Level: beginner 567 568 Notes: 569 Run with -ksp_error_if_not_converged to stop the program when the error is detected and print an error message with details. 570 571 572 .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances() 573 574 M*/ 575 576 /*MC 577 KSP_CONVERGED_ITERATING - This flag is returned if you call KSPGetConvergedReason() 578 while the KSPSolve() is still running. 579 580 Level: beginner 581 582 .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances() 583 584 M*/ 585 586 PETSC_EXTERN PetscErrorCode KSPSetConvergenceTest(KSP,PetscErrorCode (*)(KSP,PetscInt,PetscReal,KSPConvergedReason*,void*),void*,PetscErrorCode (*)(void*)); 587 PETSC_EXTERN PetscErrorCode KSPGetConvergenceContext(KSP,void**); 588 PETSC_EXTERN PetscErrorCode KSPConvergedDefault(KSP,PetscInt,PetscReal,KSPConvergedReason*,void*); 589 PETSC_EXTERN PetscErrorCode KSPConvergedLSQR(KSP,PetscInt,PetscReal,KSPConvergedReason*,void*); 590 PETSC_EXTERN PetscErrorCode KSPConvergedDefaultDestroy(void*); 591 PETSC_EXTERN PetscErrorCode KSPConvergedDefaultCreate(void**); 592 PETSC_EXTERN PetscErrorCode KSPConvergedDefaultSetUIRNorm(KSP); 593 PETSC_EXTERN PetscErrorCode KSPConvergedDefaultSetUMIRNorm(KSP); 594 PETSC_EXTERN PetscErrorCode KSPConvergedSkip(KSP,PetscInt,PetscReal,KSPConvergedReason*,void*); 595 PETSC_EXTERN PetscErrorCode KSPGetConvergedReason(KSP,KSPConvergedReason*); 596 597 PETSC_DEPRECATED("Use KSPConvergedDefault()") PETSC_STATIC_INLINE void KSPDefaultConverged(void) { /* never called */ } 598 #define KSPDefaultConverged (KSPDefaultConverged, KSPConvergedDefault) 599 PETSC_DEPRECATED("Use KSPConvergedDefaultDestroy()") PETSC_STATIC_INLINE void KSPDefaultConvergedDestroy(void) { /* never called */ } 600 #define KSPDefaultConvergedDestroy (KSPDefaultConvergedDestroy, KSPConvergedDefaultDestroy) 601 PETSC_DEPRECATED("Use KSPConvergedDefaultCreate()") PETSC_STATIC_INLINE void KSPDefaultConvergedCreate(void) { /* never called */ } 602 #define KSPDefaultConvergedCreate (KSPDefaultConvergedCreate, KSPConvergedDefaultCreate) 603 PETSC_DEPRECATED("Use KSPConvergedDefaultSetUIRNorm()") PETSC_STATIC_INLINE void KSPDefaultConvergedSetUIRNorm(void) { /* never called */ } 604 #define KSPDefaultConvergedSetUIRNorm (KSPDefaultConvergedSetUIRNorm, KSPConvergedDefaultSetUIRNorm) 605 PETSC_DEPRECATED("Use KSPConvergedDefaultSetUMIRNorm()") PETSC_STATIC_INLINE void KSPDefaultConvergedSetUMIRNorm(void) { /* never called */ } 606 #define KSPDefaultConvergedSetUMIRNorm (KSPDefaultConvergedSetUMIRNorm, KSPConvergedDefaultSetUMIRNorm) 607 PETSC_DEPRECATED("Use KSPConvergedSkip()") PETSC_STATIC_INLINE void KSPSkipConverged(void) { /* never called */ } 608 #define KSPSkipConverged (KSPSkipConverged, KSPConvergedSkip) 609 610 PETSC_EXTERN PetscErrorCode KSPComputeExplicitOperator(KSP,Mat*); 611 612 /*E 613 KSPCGType - Determines what type of CG to use 614 615 Level: beginner 616 617 .seealso: KSPCGSetType() 618 E*/ 619 typedef enum {KSP_CG_SYMMETRIC=0,KSP_CG_HERMITIAN=1} KSPCGType; 620 PETSC_EXTERN const char *const KSPCGTypes[]; 621 622 PETSC_EXTERN PetscErrorCode KSPCGSetType(KSP,KSPCGType); 623 PETSC_EXTERN PetscErrorCode KSPCGUseSingleReduction(KSP,PetscBool ); 624 625 PETSC_EXTERN PetscErrorCode KSPCGSetRadius(KSP,PetscReal); 626 PETSC_EXTERN PetscErrorCode KSPCGGetNormD(KSP,PetscReal*); 627 PETSC_EXTERN PetscErrorCode KSPCGGetObjFcn(KSP,PetscReal*); 628 629 PETSC_EXTERN PetscErrorCode KSPCGGLTRGetMinEig(KSP,PetscReal*); 630 PETSC_EXTERN PetscErrorCode KSPCGGLTRGetLambda(KSP,PetscReal*); 631 632 PETSC_EXTERN PetscErrorCode KSPPythonSetType(KSP,const char[]); 633 634 PETSC_EXTERN PetscErrorCode PCPreSolve(PC,KSP); 635 PETSC_EXTERN PetscErrorCode PCPostSolve(PC,KSP); 636 637 #include <petscdrawtypes.h> 638 PETSC_EXTERN PetscErrorCode KSPMonitorLGResidualNormCreate(MPI_Comm,const char[],const char[],int,int,int,int,PetscDrawLG*); 639 PETSC_EXTERN PetscErrorCode KSPMonitorLGResidualNorm(KSP,PetscInt,PetscReal,void*); 640 PETSC_EXTERN PetscErrorCode KSPMonitorLGTrueResidualNormCreate(MPI_Comm,const char[],const char[],int,int,int,int,PetscDrawLG*); 641 PETSC_EXTERN PetscErrorCode KSPMonitorLGTrueResidualNorm(KSP,PetscInt,PetscReal,void*); 642 PETSC_EXTERN PetscErrorCode KSPMonitorLGRange(KSP,PetscInt,PetscReal,void*); 643 644 PETSC_EXTERN PetscErrorCode PCShellSetPreSolve(PC,PetscErrorCode (*)(PC,KSP,Vec,Vec)); 645 PETSC_EXTERN PetscErrorCode PCShellSetPostSolve(PC,PetscErrorCode (*)(PC,KSP,Vec,Vec)); 646 647 /*S 648 KSPGuess - Abstract PETSc object that manages all initial guess methods in Krylov methods. 649 650 Level: beginner 651 652 Concepts: Krylov methods 653 654 .seealso: KSPCreate(), KSPSetGuessType(), KSPGuessType 655 S*/ 656 typedef struct _p_KSPGuess* KSPGuess; 657 /*J 658 KSPGuessType - String with the name of a PETSc initial guess for Krylov methods. 659 660 Level: beginner 661 662 .seealso: KSPGuess 663 J*/ 664 typedef const char* KSPGuessType; 665 #define KSPGUESSFISCHER "fischer" 666 #define KSPGUESSPOD "pod" 667 PETSC_EXTERN PetscErrorCode KSPSetGuess(KSP,KSPGuess); 668 PETSC_EXTERN PetscErrorCode KSPGetGuess(KSP,KSPGuess*); 669 PETSC_EXTERN PetscErrorCode KSPGuessView(KSPGuess,PetscViewer); 670 PETSC_EXTERN PetscErrorCode KSPGuessDestroy(KSPGuess*); 671 PETSC_EXTERN PetscErrorCode KSPGuessCreate(MPI_Comm,KSPGuess*); 672 PETSC_EXTERN PetscErrorCode KSPGuessSetType(KSPGuess,KSPGuessType); 673 PETSC_EXTERN PetscErrorCode KSPGuessGetType(KSPGuess,KSPGuessType*); 674 PETSC_EXTERN PetscErrorCode KSPGuessSetUp(KSPGuess); 675 PETSC_EXTERN PetscErrorCode KSPGuessUpdate(KSPGuess,Vec,Vec); 676 PETSC_EXTERN PetscErrorCode KSPGuessFormGuess(KSPGuess,Vec,Vec); 677 PETSC_EXTERN PetscErrorCode KSPGuessSetFromOptions(KSPGuess); 678 PETSC_EXTERN PetscErrorCode KSPGuessFischerSetModel(KSPGuess,PetscInt,PetscInt); 679 PETSC_EXTERN PetscErrorCode KSPSetUseFischerGuess(KSP,PetscInt,PetscInt); 680 PETSC_EXTERN PetscErrorCode KSPSetInitialGuessKnoll(KSP,PetscBool); 681 PETSC_EXTERN PetscErrorCode KSPGetInitialGuessKnoll(KSP,PetscBool*); 682 683 /*E 684 MatSchurComplementAinvType - Determines how to approximate the inverse of the (0,0) block in Schur complement preconditioning matrix assembly routines 685 686 Level: intermediate 687 688 .seealso: MatSchurComplementGetAinvType(), MatSchurComplementSetAinvType(), MatSchurComplementGetPmat(), MatGetSchurComplement(), MatCreateSchurComplementPmat() 689 E*/ 690 typedef enum {MAT_SCHUR_COMPLEMENT_AINV_DIAG, MAT_SCHUR_COMPLEMENT_AINV_LUMP, MAT_SCHUR_COMPLEMENT_AINV_BLOCK_DIAG} MatSchurComplementAinvType; 691 PETSC_EXTERN const char *const MatSchurComplementAinvTypes[]; 692 693 PETSC_EXTERN PetscErrorCode MatCreateSchurComplement(Mat,Mat,Mat,Mat,Mat,Mat*); 694 PETSC_EXTERN PetscErrorCode MatSchurComplementGetKSP(Mat,KSP*); 695 PETSC_EXTERN PetscErrorCode MatSchurComplementSetKSP(Mat,KSP); 696 PETSC_EXTERN PetscErrorCode MatSchurComplementSetSubMatrices(Mat,Mat,Mat,Mat,Mat,Mat); 697 PETSC_EXTERN PetscErrorCode MatSchurComplementUpdateSubMatrices(Mat,Mat,Mat,Mat,Mat,Mat); 698 PETSC_EXTERN PetscErrorCode MatSchurComplementGetSubMatrices(Mat,Mat*,Mat*,Mat*,Mat*,Mat*); 699 PETSC_EXTERN PetscErrorCode MatSchurComplementSetAinvType(Mat,MatSchurComplementAinvType); 700 PETSC_EXTERN PetscErrorCode MatSchurComplementGetAinvType(Mat,MatSchurComplementAinvType*); 701 PETSC_EXTERN PetscErrorCode MatSchurComplementGetPmat(Mat,MatReuse,Mat*); 702 PETSC_EXTERN PetscErrorCode MatSchurComplementComputeExplicitOperator(Mat,Mat*); 703 PETSC_EXTERN PetscErrorCode MatGetSchurComplement(Mat,IS,IS,IS,IS,MatReuse,Mat*,MatSchurComplementAinvType,MatReuse,Mat*); 704 PETSC_EXTERN PetscErrorCode MatCreateSchurComplementPmat(Mat,Mat,Mat,Mat,MatSchurComplementAinvType,MatReuse,Mat*); 705 706 PETSC_EXTERN PetscErrorCode KSPSetDM(KSP,DM); 707 PETSC_EXTERN PetscErrorCode KSPSetDMActive(KSP,PetscBool ); 708 PETSC_EXTERN PetscErrorCode KSPGetDM(KSP,DM*); 709 PETSC_EXTERN PetscErrorCode KSPSetApplicationContext(KSP,void*); 710 PETSC_EXTERN PetscErrorCode KSPGetApplicationContext(KSP,void*); 711 PETSC_EXTERN PetscErrorCode KSPSetComputeRHS(KSP,PetscErrorCode (*func)(KSP,Vec,void*),void*); 712 PETSC_EXTERN PetscErrorCode KSPSetComputeOperators(KSP,PetscErrorCode(*)(KSP,Mat,Mat,void*),void*); 713 PETSC_EXTERN PetscErrorCode KSPSetComputeInitialGuess(KSP,PetscErrorCode(*)(KSP,Vec,void*),void*); 714 PETSC_EXTERN PetscErrorCode DMKSPSetComputeOperators(DM,PetscErrorCode(*)(KSP,Mat,Mat,void*),void*); 715 PETSC_EXTERN PetscErrorCode DMKSPGetComputeOperators(DM,PetscErrorCode(**)(KSP,Mat,Mat,void*),void*); 716 PETSC_EXTERN PetscErrorCode DMKSPSetComputeRHS(DM,PetscErrorCode(*)(KSP,Vec,void*),void*); 717 PETSC_EXTERN PetscErrorCode DMKSPGetComputeRHS(DM,PetscErrorCode(**)(KSP,Vec,void*),void*); 718 PETSC_EXTERN PetscErrorCode DMKSPSetComputeInitialGuess(DM,PetscErrorCode(*)(KSP,Vec,void*),void*); 719 PETSC_EXTERN PetscErrorCode DMKSPGetComputeInitialGuess(DM,PetscErrorCode(**)(KSP,Vec,void*),void*); 720 721 PETSC_EXTERN PetscErrorCode DMGlobalToLocalSolve(DM,Vec,Vec); 722 PETSC_EXTERN PetscErrorCode DMProjectField(DM, PetscReal, Vec, 723 void (**)(PetscInt, PetscInt, PetscInt, 724 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 725 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 726 PetscReal, const PetscReal [], PetscInt, const PetscScalar[], PetscScalar []), InsertMode, Vec); 727 #endif 728