Actual source code: shift.c

  1: /*
  2:     Shift spectral transformation, applies (A + sigma I) as operator, or 
  3:     inv(B)(A + sigma B) for generalized problems

  5:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  6:       SLEPc - Scalable Library for Eigenvalue Problem Computations
  7:       Copyright (c) 2002-2007, Universidad Politecnica de Valencia, Spain

  9:       This file is part of SLEPc. See the README file for conditions of use
 10:       and additional information.
 11:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 12: */

 14:  #include src/st/stimpl.h

 18: PetscErrorCode STApply_Shift(ST st,Vec x,Vec y)
 19: {

 23:   if (st->B) {
 24:     /* generalized eigenproblem: y = (B^-1 A + sI) x */
 25:     MatMult(st->A,x,st->w);
 26:     STAssociatedKSPSolve(st,st->w,y);
 27:   }
 28:   else {
 29:     /* standard eigenproblem: y = (A + sI) x */
 30:     MatMult(st->A,x,y);
 31:   }
 32:   if (st->sigma != 0.0) {
 33:     VecAXPY(y,st->sigma,x);
 34:   }
 35:   return(0);
 36: }

 40: PetscErrorCode STApplyTranspose_Shift(ST st,Vec x,Vec y)
 41: {

 45:   if (st->B) {
 46:     /* generalized eigenproblem: y = (A^T B^-T + sI) x */
 47:     STAssociatedKSPSolveTranspose(st,x,st->w);
 48:     MatMultTranspose(st->A,st->w,y);
 49:   }
 50:   else {
 51:     /* standard eigenproblem: y = (A^T + sI) x */
 52:     MatMultTranspose(st->A,x,y);
 53:   }
 54:   if (st->sigma != 0.0) {
 55:     VecAXPY(y,st->sigma,x);
 56:   }
 57:   return(0);
 58: }

 62: PetscErrorCode STBackTransform_Shift(ST st,PetscScalar *eigr,PetscScalar *eigi)
 63: {
 65:   if (eigr) *eigr -= st->sigma;
 66:   return(0);
 67: }

 71: PetscErrorCode STSetUp_Shift(ST st)
 72: {

 76:   if (st->B) {
 77:     KSPSetOperators(st->ksp,st->B,st->B,DIFFERENT_NONZERO_PATTERN);
 78:     KSPSetUp(st->ksp);
 79:   }
 80:   return(0);
 81: }

 85: PetscErrorCode STView_Shift(ST st,PetscViewer viewer)
 86: {

 90:   if (st->B) {
 91:     STView_Default(st,viewer);
 92:   }
 93:   return(0);
 94: }

 99: PetscErrorCode STCreate_Shift(ST st)
100: {
102:   st->ops->apply           = STApply_Shift;
103:   st->ops->getbilinearform = STGetBilinearForm_Default;
104:   st->ops->applytrans      = STApplyTranspose_Shift;
105:   st->ops->backtr          = STBackTransform_Shift;
106:   st->ops->setup           = STSetUp_Shift;
107:   st->ops->view            = STView_Shift;
108:   st->checknullspace       = 0;
109:   return(0);
110: }