Actual source code: borthog.c

  1: /*
  2:     Routines used for the orthogonalization of the Hessenberg matrix.

  4:     Note that for the complex numbers version, the VecDot() and
  5:     VecMDot() arguments within the code MUST remain in the order
  6:     given for correct computation of inner products.
  7: */
 8:  #include src/ksp/ksp/impls/gmres/gmresp.h

 10: /*@C
 11:      KSPGMRESModifiedGramSchmidtOrthogonalization -  This is the basic orthogonalization routine 
 12:                 using modified Gram-Schmidt.

 14:      Collective on KSP

 16:   Input Parameters:
 17: +   ksp - KSP object, must be associated with GMRES, FGMRES, or LGMRES Krylov method
 18: -   its - one less then the current GMRES restart iteration, i.e. the size of the Krylov space

 20:    Options Database Keys:
 21: .  -ksp_gmres_modifiedgramschmidt - Activates KSPGMRESModifiedGramSchmidtOrthogonalization()

 23:    Level: intermediate

 25: .seealso:  KSPGMRESSetOrthogonalization(), KSPGMRESClassicalGramSchmidtOrthogonalization()

 27: @*/
 30: PetscErrorCode KSPGMRESModifiedGramSchmidtOrthogonalization(KSP ksp,PetscInt it)
 31: {
 32:   KSP_GMRES      *gmres = (KSP_GMRES *)(ksp->data);
 34:   PetscInt       j;
 35:   PetscScalar    *hh,*hes,tmp;

 38:   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:     VecDot(VEC_VV(it+1),VEC_VV(j),hh);
 45:     *hes++ = *hh;
 46:     /* vv(it+1) <- vv(it+1) - hh[it+1][j] vv(j) */
 47:     tmp    = - (*hh++);
 48:     VecAXPY(&tmp,VEC_VV(j),VEC_VV(it+1));
 49:   }
 50:   PetscLogEventEnd(KSP_GMRESOrthogonalization,ksp,0,0,0);
 51:   return(0);
 52: }