Actual source code: borthog.c
2: /*
3: Routines used for the orthogonalization of the Hessenberg matrix.
5: Note that for the complex numbers version, the VecDot() and
6: VecMDot() arguments within the code MUST remain in the order
7: given for correct computation of inner products.
8: */
9: #include <../src/ksp/ksp/impls/gmres/gmresimpl.h>
11: /*@C
12: KSPGMRESModifiedGramSchmidtOrthogonalization - This is the basic orthogonalization routine
13: using modified Gram-Schmidt.
15: Collective
17: Input Parameters:
18: + ksp - KSP object, must be associated with `KSPGMRES`, `KSPFGMRES`, or `KSPLGMRES` Krylov method
19: - its - one less then the current GMRES restart iteration, i.e. the size of the Krylov space
21: Options Database Keys:
22: . -ksp_gmres_modifiedgramschmidt - Activates `KSPGMRESModifiedGramSchmidtOrthogonalization()`
24: Level: intermediate
26: Notes:
27: In general this is much slower than `KSPGMRESClassicalGramSchmidtOrthogonalization()` but has better stability properties.
29: .seealso: [](ch_ksp), `KSPGMRESSetOrthogonalization()`, `KSPGMRESClassicalGramSchmidtOrthogonalization()`, `KSPGMRESGetOrthogonalization()`
30: @*/
31: PetscErrorCode KSPGMRESModifiedGramSchmidtOrthogonalization(KSP ksp, PetscInt it)
32: {
33: KSP_GMRES *gmres = (KSP_GMRES *)(ksp->data);
34: PetscInt j;
35: PetscScalar *hh, *hes;
37: PetscFunctionBegin;
38: PetscCall(PetscLogEventBegin(KSP_GMRESOrthogonalization, ksp, 0, 0, 0));
39: /* update Hessenberg matrix and do Gram-Schmidt */
40: hh = HH(0, it);
41: hes = HES(0, it);
42: for (j = 0; j <= it; j++) {
43: /* (vv(it+1), vv(j)) */
44: PetscCall(VecDot(VEC_VV(it + 1), VEC_VV(j), hh));
45: KSPCheckDot(ksp, *hh);
46: if (ksp->reason) break;
47: *hes++ = *hh;
48: /* vv(it+1) <- vv(it+1) - hh[it+1][j] vv(j) */
49: PetscCall(VecAXPY(VEC_VV(it + 1), -(*hh++), VEC_VV(j)));
50: }
51: PetscCall(PetscLogEventEnd(KSP_GMRESOrthogonalization, ksp, 0, 0, 0));
52: PetscFunctionReturn(PETSC_SUCCESS);
53: }