1 2 /* 3 This file implements PGMRES (a Pipelined Generalized Minimal Residual method) 4 */ 5 6 #include <../src/ksp/ksp/impls/gmres/pgmres/pgmresimpl.h> /*I "petscksp.h" I*/ 7 #define PGMRES_DELTA_DIRECTIONS 10 8 #define PGMRES_DEFAULT_MAXK 30 9 10 static PetscErrorCode KSPPGMRESUpdateHessenberg(KSP,PetscInt,PetscBool*,PetscReal*); 11 static PetscErrorCode KSPPGMRESBuildSoln(PetscScalar*,Vec,Vec,KSP,PetscInt); 12 13 /* 14 15 KSPSetUp_PGMRES - Sets up the workspace needed by pgmres. 16 17 This is called once, usually automatically by KSPSolve() or KSPSetUp(), 18 but can be called directly by KSPSetUp(). 19 20 */ 21 static PetscErrorCode KSPSetUp_PGMRES(KSP ksp) 22 { 23 PetscFunctionBegin; 24 PetscCall(KSPSetUp_GMRES(ksp)); 25 PetscFunctionReturn(0); 26 } 27 28 /* 29 30 KSPPGMRESCycle - Run pgmres, possibly with restart. Return residual 31 history if requested. 32 33 input parameters: 34 . pgmres - structure containing parameters and work areas 35 36 output parameters: 37 . itcount - number of iterations used. If null, ignored. 38 . converged - 0 if not converged 39 40 Notes: 41 On entry, the value in vector VEC_VV(0) should be 42 the initial residual. 43 44 */ 45 static PetscErrorCode KSPPGMRESCycle(PetscInt *itcount,KSP ksp) 46 { 47 KSP_PGMRES *pgmres = (KSP_PGMRES*)(ksp->data); 48 PetscReal res_norm,res,newnorm; 49 PetscInt it = 0,j,k; 50 PetscBool hapend = PETSC_FALSE; 51 52 PetscFunctionBegin; 53 if (itcount) *itcount = 0; 54 PetscCall(VecNormalize(VEC_VV(0),&res_norm)); 55 KSPCheckNorm(ksp,res_norm); 56 res = res_norm; 57 *RS(0) = res_norm; 58 59 /* check for the convergence */ 60 PetscCall(PetscObjectSAWsTakeAccess((PetscObject)ksp)); 61 if (ksp->normtype != KSP_NORM_NONE) ksp->rnorm = res; 62 else ksp->rnorm = 0; 63 PetscCall(PetscObjectSAWsGrantAccess((PetscObject)ksp)); 64 pgmres->it = it-2; 65 PetscCall(KSPLogResidualHistory(ksp,ksp->rnorm)); 66 PetscCall(KSPMonitor(ksp,ksp->its,ksp->rnorm)); 67 if (!res) { 68 ksp->reason = KSP_CONVERGED_ATOL; 69 PetscCall(PetscInfo(ksp,"Converged due to zero residual norm on entry\n")); 70 PetscFunctionReturn(0); 71 } 72 73 PetscCall((*ksp->converged)(ksp,ksp->its,ksp->rnorm,&ksp->reason,ksp->cnvP)); 74 for (; !ksp->reason; it++) { 75 Vec Zcur,Znext; 76 if (pgmres->vv_allocated <= it + VEC_OFFSET + 1) { 77 PetscCall(KSPGMRESGetNewVectors(ksp,it+1)); 78 } 79 /* VEC_VV(it-1) is orthogonal, it will be normalized once the VecNorm arrives. */ 80 Zcur = VEC_VV(it); /* Zcur is not yet orthogonal, but the VecMDot to orthogonalize it has been started. */ 81 Znext = VEC_VV(it+1); /* This iteration will compute Znext, update with a deferred correction once we know how 82 * Zcur relates to the previous vectors, and start the reduction to orthogonalize it. */ 83 84 if (it < pgmres->max_k+1 && ksp->its+1 < PetscMax(2,ksp->max_it)) { /* We don't know whether what we have computed is enough, so apply the matrix. */ 85 PetscCall(KSP_PCApplyBAorAB(ksp,Zcur,Znext,VEC_TEMP_MATOP)); 86 } 87 88 if (it > 1) { /* Complete the pending reduction */ 89 PetscCall(VecNormEnd(VEC_VV(it-1),NORM_2,&newnorm)); 90 *HH(it-1,it-2) = newnorm; 91 } 92 if (it > 0) { /* Finish the reduction computing the latest column of H */ 93 PetscCall(VecMDotEnd(Zcur,it,&(VEC_VV(0)),HH(0,it-1))); 94 } 95 96 if (it > 1) { 97 /* normalize the base vector from two iterations ago, basis is complete up to here */ 98 PetscCall(VecScale(VEC_VV(it-1),1./ *HH(it-1,it-2))); 99 100 PetscCall(KSPPGMRESUpdateHessenberg(ksp,it-2,&hapend,&res)); 101 pgmres->it = it-2; 102 ksp->its++; 103 if (ksp->normtype != KSP_NORM_NONE) ksp->rnorm = res; 104 else ksp->rnorm = 0; 105 106 PetscCall((*ksp->converged)(ksp,ksp->its,ksp->rnorm,&ksp->reason,ksp->cnvP)); 107 if (it < pgmres->max_k+1 || ksp->reason || ksp->its == ksp->max_it) { /* Monitor if we are done or still iterating, but not before a restart. */ 108 PetscCall(KSPLogResidualHistory(ksp,ksp->rnorm)); 109 PetscCall(KSPMonitor(ksp,ksp->its,ksp->rnorm)); 110 } 111 if (ksp->reason) break; 112 /* Catch error in happy breakdown and signal convergence and break from loop */ 113 if (hapend) { 114 PetscCheck(!ksp->errorifnotconverged,PetscObjectComm((PetscObject)ksp),PETSC_ERR_NOT_CONVERGED,"You reached the happy break down, but convergence was not indicated. Residual norm = %g",(double)res); 115 ksp->reason = KSP_DIVERGED_BREAKDOWN; 116 break; 117 } 118 119 if (!(it < pgmres->max_k+1 && ksp->its < ksp->max_it)) break; 120 121 /* The it-2 column of H was not scaled when we computed Zcur, apply correction */ 122 PetscCall(VecScale(Zcur,1./ *HH(it-1,it-2))); 123 /* And Znext computed in this iteration was computed using the under-scaled Zcur */ 124 PetscCall(VecScale(Znext,1./ *HH(it-1,it-2))); 125 126 /* In the previous iteration, we projected an unnormalized Zcur against the Krylov basis, so we need to fix the column of H resulting from that projection. */ 127 for (k=0; k<it; k++) *HH(k,it-1) /= *HH(it-1,it-2); 128 /* When Zcur was projected against the Krylov basis, VV(it-1) was still not normalized, so fix that too. This 129 * column is complete except for HH(it,it-1) which we won't know until the next iteration. */ 130 *HH(it-1,it-1) /= *HH(it-1,it-2); 131 } 132 133 if (it > 0) { 134 PetscScalar *work; 135 if (!pgmres->orthogwork) PetscCall(PetscMalloc1(pgmres->max_k + 2,&pgmres->orthogwork)); 136 work = pgmres->orthogwork; 137 /* Apply correction computed by the VecMDot in the last iteration to Znext. The original form is 138 * 139 * Znext -= sum_{j=0}^{i-1} Z[j+1] * H[j,i-1] 140 * 141 * where 142 * 143 * Z[j] = sum_{k=0}^j V[k] * H[k,j-1] 144 * 145 * substituting 146 * 147 * Znext -= sum_{j=0}^{i-1} sum_{k=0}^{j+1} V[k] * H[k,j] * H[j,i-1] 148 * 149 * rearranging the iteration space from row-column to column-row 150 * 151 * Znext -= sum_{k=0}^i sum_{j=k-1}^{i-1} V[k] * H[k,j] * H[j,i-1] 152 * 153 * Note that column it-1 of HH is correct. For all previous columns, we must look at HES because HH has already 154 * been transformed to upper triangular form. 155 */ 156 for (k=0; k<it+1; k++) { 157 work[k] = 0; 158 for (j=PetscMax(0,k-1); j<it-1; j++) work[k] -= *HES(k,j) * *HH(j,it-1); 159 } 160 PetscCall(VecMAXPY(Znext,it+1,work,&VEC_VV(0))); 161 PetscCall(VecAXPY(Znext,-*HH(it-1,it-1),Zcur)); 162 163 /* Orthogonalize Zcur against existing basis vectors. */ 164 for (k=0; k<it; k++) work[k] = -*HH(k,it-1); 165 PetscCall(VecMAXPY(Zcur,it,work,&VEC_VV(0))); 166 /* Zcur is now orthogonal, and will be referred to as VEC_VV(it) again, though it is still not normalized. */ 167 /* Begin computing the norm of the new vector, will be normalized after the MatMult in the next iteration. */ 168 PetscCall(VecNormBegin(VEC_VV(it),NORM_2,&newnorm)); 169 } 170 171 /* Compute column of H (to the diagonal, but not the subdiagonal) to be able to orthogonalize the newest vector. */ 172 PetscCall(VecMDotBegin(Znext,it+1,&VEC_VV(0),HH(0,it))); 173 174 /* Start an asynchronous split-mode reduction, the result of the MDot and Norm will be collected on the next iteration. */ 175 PetscCall(PetscCommSplitReductionBegin(PetscObjectComm((PetscObject)Znext))); 176 } 177 178 if (itcount) *itcount = it-1; /* Number of iterations actually completed. */ 179 180 /* 181 Down here we have to solve for the "best" coefficients of the Krylov 182 columns, add the solution values together, and possibly unwind the 183 preconditioning from the solution 184 */ 185 /* Form the solution (or the solution so far) */ 186 PetscCall(KSPPGMRESBuildSoln(RS(0),ksp->vec_sol,ksp->vec_sol,ksp,it-2)); 187 PetscFunctionReturn(0); 188 } 189 190 /* 191 KSPSolve_PGMRES - This routine applies the PGMRES method. 192 193 Input Parameter: 194 . ksp - the Krylov space object that was set to use pgmres 195 196 Output Parameter: 197 . outits - number of iterations used 198 199 */ 200 static PetscErrorCode KSPSolve_PGMRES(KSP ksp) 201 { 202 PetscInt its,itcount; 203 KSP_PGMRES *pgmres = (KSP_PGMRES*)ksp->data; 204 PetscBool guess_zero = ksp->guess_zero; 205 206 PetscFunctionBegin; 207 PetscCheck(!ksp->calc_sings || pgmres->Rsvd,PetscObjectComm((PetscObject)ksp),PETSC_ERR_ORDER,"Must call KSPSetComputeSingularValues() before KSPSetUp() is called"); 208 PetscCall(PetscObjectSAWsTakeAccess((PetscObject)ksp)); 209 ksp->its = 0; 210 PetscCall(PetscObjectSAWsGrantAccess((PetscObject)ksp)); 211 212 itcount = 0; 213 ksp->reason = KSP_CONVERGED_ITERATING; 214 while (!ksp->reason) { 215 PetscCall(KSPInitialResidual(ksp,ksp->vec_sol,VEC_TEMP,VEC_TEMP_MATOP,VEC_VV(0),ksp->vec_rhs)); 216 PetscCall(KSPPGMRESCycle(&its,ksp)); 217 itcount += its; 218 if (itcount >= ksp->max_it) { 219 if (!ksp->reason) ksp->reason = KSP_DIVERGED_ITS; 220 break; 221 } 222 ksp->guess_zero = PETSC_FALSE; /* every future call to KSPInitialResidual() will have nonzero guess */ 223 } 224 ksp->guess_zero = guess_zero; /* restore if user provided nonzero initial guess */ 225 PetscFunctionReturn(0); 226 } 227 228 static PetscErrorCode KSPDestroy_PGMRES(KSP ksp) 229 { 230 PetscFunctionBegin; 231 PetscCall(KSPDestroy_GMRES(ksp)); 232 PetscFunctionReturn(0); 233 } 234 235 /* 236 KSPPGMRESBuildSoln - create the solution from the starting vector and the 237 current iterates. 238 239 Input parameters: 240 nrs - work area of size it + 1. 241 vguess - index of initial guess 242 vdest - index of result. Note that vguess may == vdest (replace 243 guess with the solution). 244 it - HH upper triangular part is a block of size (it+1) x (it+1) 245 246 This is an internal routine that knows about the PGMRES internals. 247 */ 248 static PetscErrorCode KSPPGMRESBuildSoln(PetscScalar *nrs,Vec vguess,Vec vdest,KSP ksp,PetscInt it) 249 { 250 PetscScalar tt; 251 PetscInt k,j; 252 KSP_PGMRES *pgmres = (KSP_PGMRES*)(ksp->data); 253 254 PetscFunctionBegin; 255 /* Solve for solution vector that minimizes the residual */ 256 257 if (it < 0) { /* no pgmres steps have been performed */ 258 PetscCall(VecCopy(vguess,vdest)); /* VecCopy() is smart, exits immediately if vguess == vdest */ 259 PetscFunctionReturn(0); 260 } 261 262 /* solve the upper triangular system - RS is the right side and HH is 263 the upper triangular matrix - put soln in nrs */ 264 if (*HH(it,it) != 0.0) nrs[it] = *RS(it) / *HH(it,it); 265 else nrs[it] = 0.0; 266 267 for (k=it-1; k>=0; k--) { 268 tt = *RS(k); 269 for (j=k+1; j<=it; j++) tt -= *HH(k,j) * nrs[j]; 270 nrs[k] = tt / *HH(k,k); 271 } 272 273 /* Accumulate the correction to the solution of the preconditioned problem in TEMP */ 274 PetscCall(VecZeroEntries(VEC_TEMP)); 275 PetscCall(VecMAXPY(VEC_TEMP,it+1,nrs,&VEC_VV(0))); 276 PetscCall(KSPUnwindPreconditioner(ksp,VEC_TEMP,VEC_TEMP_MATOP)); 277 /* add solution to previous solution */ 278 if (vdest == vguess) { 279 PetscCall(VecAXPY(vdest,1.0,VEC_TEMP)); 280 } else { 281 PetscCall(VecWAXPY(vdest,1.0,VEC_TEMP,vguess)); 282 } 283 PetscFunctionReturn(0); 284 } 285 286 /* 287 288 KSPPGMRESUpdateHessenberg - Do the scalar work for the orthogonalization. 289 Return new residual. 290 291 input parameters: 292 293 . ksp - Krylov space object 294 . it - plane rotations are applied to the (it+1)th column of the 295 modified hessenberg (i.e. HH(:,it)) 296 . hapend - PETSC_FALSE not happy breakdown ending. 297 298 output parameters: 299 . res - the new residual 300 301 */ 302 /* 303 . it - column of the Hessenberg that is complete, PGMRES is actually computing two columns ahead of this 304 */ 305 static PetscErrorCode KSPPGMRESUpdateHessenberg(KSP ksp,PetscInt it,PetscBool *hapend,PetscReal *res) 306 { 307 PetscScalar *hh,*cc,*ss,*rs; 308 PetscInt j; 309 PetscReal hapbnd; 310 KSP_PGMRES *pgmres = (KSP_PGMRES*)(ksp->data); 311 312 PetscFunctionBegin; 313 hh = HH(0,it); /* pointer to beginning of column to update */ 314 cc = CC(0); /* beginning of cosine rotations */ 315 ss = SS(0); /* beginning of sine rotations */ 316 rs = RS(0); /* right hand side of least squares system */ 317 318 /* The Hessenberg matrix is now correct through column it, save that form for possible spectral analysis */ 319 for (j=0; j<=it+1; j++) *HES(j,it) = hh[j]; 320 321 /* check for the happy breakdown */ 322 hapbnd = PetscMin(PetscAbsScalar(hh[it+1] / rs[it]),pgmres->haptol); 323 if (PetscAbsScalar(hh[it+1]) < hapbnd) { 324 PetscCall(PetscInfo(ksp,"Detected happy breakdown, current hapbnd = %14.12e H(%" PetscInt_FMT ",%" PetscInt_FMT ") = %14.12e\n",(double)hapbnd,it+1,it,(double)PetscAbsScalar(*HH(it+1,it)))); 325 *hapend = PETSC_TRUE; 326 } 327 328 /* Apply all the previously computed plane rotations to the new column 329 of the Hessenberg matrix */ 330 /* Note: this uses the rotation [conj(c) s ; -s c], c= cos(theta), s= sin(theta), 331 and some refs have [c s ; -conj(s) c] (don't be confused!) */ 332 333 for (j=0; j<it; j++) { 334 PetscScalar hhj = hh[j]; 335 hh[j] = PetscConj(cc[j])*hhj + ss[j]*hh[j+1]; 336 hh[j+1] = -ss[j] *hhj + cc[j]*hh[j+1]; 337 } 338 339 /* 340 compute the new plane rotation, and apply it to: 341 1) the right-hand-side of the Hessenberg system (RS) 342 note: it affects RS(it) and RS(it+1) 343 2) the new column of the Hessenberg matrix 344 note: it affects HH(it,it) which is currently pointed to 345 by hh and HH(it+1, it) (*(hh+1)) 346 thus obtaining the updated value of the residual... 347 */ 348 349 /* compute new plane rotation */ 350 351 if (!*hapend) { 352 PetscReal delta = PetscSqrtReal(PetscSqr(PetscAbsScalar(hh[it])) + PetscSqr(PetscAbsScalar(hh[it+1]))); 353 if (delta == 0.0) { 354 ksp->reason = KSP_DIVERGED_NULL; 355 PetscFunctionReturn(0); 356 } 357 358 cc[it] = hh[it] / delta; /* new cosine value */ 359 ss[it] = hh[it+1] / delta; /* new sine value */ 360 361 hh[it] = PetscConj(cc[it])*hh[it] + ss[it]*hh[it+1]; 362 rs[it+1] = -ss[it]*rs[it]; 363 rs[it] = PetscConj(cc[it])*rs[it]; 364 *res = PetscAbsScalar(rs[it+1]); 365 } else { /* happy breakdown: HH(it+1, it) = 0, therefore we don't need to apply 366 another rotation matrix (so RH doesn't change). The new residual is 367 always the new sine term times the residual from last time (RS(it)), 368 but now the new sine rotation would be zero...so the residual should 369 be zero...so we will multiply "zero" by the last residual. This might 370 not be exactly what we want to do here -could just return "zero". */ 371 372 *res = 0.0; 373 } 374 PetscFunctionReturn(0); 375 } 376 377 /* 378 KSPBuildSolution_PGMRES 379 380 Input Parameter: 381 . ksp - the Krylov space object 382 . ptr- 383 384 Output Parameter: 385 . result - the solution 386 387 Note: this calls KSPPGMRESBuildSoln - the same function that KSPPGMRESCycle 388 calls directly. 389 390 */ 391 PetscErrorCode KSPBuildSolution_PGMRES(KSP ksp,Vec ptr,Vec *result) 392 { 393 KSP_PGMRES *pgmres = (KSP_PGMRES*)ksp->data; 394 395 PetscFunctionBegin; 396 if (!ptr) { 397 if (!pgmres->sol_temp) { 398 PetscCall(VecDuplicate(ksp->vec_sol,&pgmres->sol_temp)); 399 PetscCall(PetscLogObjectParent((PetscObject)ksp,(PetscObject)pgmres->sol_temp)); 400 } 401 ptr = pgmres->sol_temp; 402 } 403 if (!pgmres->nrs) { 404 /* allocate the work area */ 405 PetscCall(PetscMalloc1(pgmres->max_k,&pgmres->nrs)); 406 PetscCall(PetscLogObjectMemory((PetscObject)ksp,pgmres->max_k*sizeof(PetscScalar))); 407 } 408 409 PetscCall(KSPPGMRESBuildSoln(pgmres->nrs,ksp->vec_sol,ptr,ksp,pgmres->it)); 410 if (result) *result = ptr; 411 PetscFunctionReturn(0); 412 } 413 414 PetscErrorCode KSPSetFromOptions_PGMRES(PetscOptionItems *PetscOptionsObject,KSP ksp) 415 { 416 PetscFunctionBegin; 417 PetscCall(KSPSetFromOptions_GMRES(PetscOptionsObject,ksp)); 418 PetscOptionsHeadBegin(PetscOptionsObject,"KSP pipelined GMRES Options"); 419 PetscOptionsHeadEnd(); 420 PetscFunctionReturn(0); 421 } 422 423 PetscErrorCode KSPReset_PGMRES(KSP ksp) 424 { 425 PetscFunctionBegin; 426 PetscCall(KSPReset_GMRES(ksp)); 427 PetscFunctionReturn(0); 428 } 429 430 /*MC 431 KSPPGMRES - Implements the Pipelined Generalized Minimal Residual method. 432 433 Options Database Keys: 434 + -ksp_gmres_restart <restart> - the number of Krylov directions to orthogonalize against 435 . -ksp_gmres_haptol <tol> - sets the tolerance for "happy ending" (exact convergence) 436 . -ksp_gmres_preallocate - preallocate all the Krylov search directions initially (otherwise groups of 437 vectors are allocated as needed) 438 . -ksp_gmres_classicalgramschmidt - use classical (unmodified) Gram-Schmidt to orthogonalize against the Krylov space (fast) (the default) 439 . -ksp_gmres_modifiedgramschmidt - use modified Gram-Schmidt in the orthogonalization (more stable, but slower) 440 . -ksp_gmres_cgs_refinement_type <refine_never,refine_ifneeded,refine_always> - determine if iterative refinement is used to increase the 441 stability of the classical Gram-Schmidt orthogonalization. 442 - -ksp_gmres_krylov_monitor - plot the Krylov space generated 443 444 Level: beginner 445 446 Notes: 447 MPI configuration may be necessary for reductions to make asynchronous progress, which is important for performance of pipelined methods. 448 See the FAQ on the PETSc website for details. 449 450 Reference: 451 Ghysels, Ashby, Meerbergen, Vanroose, Hiding global communication latencies in the GMRES algorithm on massively parallel machines, 2012. 452 453 Developer Notes: 454 This object is subclassed off of KSPGMRES 455 456 .seealso: `KSPCreate()`, `KSPSetType()`, `KSPType`, `KSP`, `KSPGMRES`, `KSPLGMRES`, `KSPPIPECG`, `KSPPIPECR`, 457 `KSPGMRESSetRestart()`, `KSPGMRESSetHapTol()`, `KSPGMRESSetPreAllocateVectors()`, `KSPGMRESSetOrthogonalization()`, `KSPGMRESGetOrthogonalization()`, 458 `KSPGMRESClassicalGramSchmidtOrthogonalization()`, `KSPGMRESModifiedGramSchmidtOrthogonalization()`, 459 `KSPGMRESCGSRefinementType`, `KSPGMRESSetCGSRefinementType()`, `KSPGMRESGetCGSRefinementType()`, `KSPGMRESMonitorKrylov()` 460 M*/ 461 462 PETSC_EXTERN PetscErrorCode KSPCreate_PGMRES(KSP ksp) 463 { 464 KSP_PGMRES *pgmres; 465 466 PetscFunctionBegin; 467 PetscCall(PetscNewLog(ksp,&pgmres)); 468 469 ksp->data = (void*)pgmres; 470 ksp->ops->buildsolution = KSPBuildSolution_PGMRES; 471 ksp->ops->setup = KSPSetUp_PGMRES; 472 ksp->ops->solve = KSPSolve_PGMRES; 473 ksp->ops->reset = KSPReset_PGMRES; 474 ksp->ops->destroy = KSPDestroy_PGMRES; 475 ksp->ops->view = KSPView_GMRES; 476 ksp->ops->setfromoptions = KSPSetFromOptions_PGMRES; 477 ksp->ops->computeextremesingularvalues = KSPComputeExtremeSingularValues_GMRES; 478 ksp->ops->computeeigenvalues = KSPComputeEigenvalues_GMRES; 479 480 PetscCall(KSPSetSupportedNorm(ksp,KSP_NORM_PRECONDITIONED,PC_LEFT,3)); 481 PetscCall(KSPSetSupportedNorm(ksp,KSP_NORM_UNPRECONDITIONED,PC_RIGHT,2)); 482 PetscCall(KSPSetSupportedNorm(ksp,KSP_NORM_NONE,PC_RIGHT,1)); 483 484 PetscCall(PetscObjectComposeFunction((PetscObject)ksp,"KSPGMRESSetPreAllocateVectors_C",KSPGMRESSetPreAllocateVectors_GMRES)); 485 PetscCall(PetscObjectComposeFunction((PetscObject)ksp,"KSPGMRESSetOrthogonalization_C",KSPGMRESSetOrthogonalization_GMRES)); 486 PetscCall(PetscObjectComposeFunction((PetscObject)ksp,"KSPGMRESGetOrthogonalization_C",KSPGMRESGetOrthogonalization_GMRES)); 487 PetscCall(PetscObjectComposeFunction((PetscObject)ksp,"KSPGMRESSetRestart_C",KSPGMRESSetRestart_GMRES)); 488 PetscCall(PetscObjectComposeFunction((PetscObject)ksp,"KSPGMRESGetRestart_C",KSPGMRESGetRestart_GMRES)); 489 PetscCall(PetscObjectComposeFunction((PetscObject)ksp,"KSPGMRESSetCGSRefinementType_C",KSPGMRESSetCGSRefinementType_GMRES)); 490 PetscCall(PetscObjectComposeFunction((PetscObject)ksp,"KSPGMRESGetCGSRefinementType_C",KSPGMRESGetCGSRefinementType_GMRES)); 491 492 pgmres->nextra_vecs = 1; 493 pgmres->haptol = 1.0e-30; 494 pgmres->q_preallocate = 0; 495 pgmres->delta_allocate = PGMRES_DELTA_DIRECTIONS; 496 pgmres->orthog = KSPGMRESClassicalGramSchmidtOrthogonalization; 497 pgmres->nrs = NULL; 498 pgmres->sol_temp = NULL; 499 pgmres->max_k = PGMRES_DEFAULT_MAXK; 500 pgmres->Rsvd = NULL; 501 pgmres->orthogwork = NULL; 502 pgmres->cgstype = KSP_GMRES_CGS_REFINE_NEVER; 503 PetscFunctionReturn(0); 504 } 505