Actual source code: ipform.c

  1: /*
  2:      Routines for setting the bilinear form

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

  8:    This file is part of SLEPc.
  9:       
 10:    SLEPc is free software: you can redistribute it and/or modify it under  the
 11:    terms of version 3 of the GNU Lesser General Public License as published by
 12:    the Free Software Foundation.

 14:    SLEPc  is  distributed in the hope that it will be useful, but WITHOUT  ANY 
 15:    WARRANTY;  without even the implied warranty of MERCHANTABILITY or  FITNESS 
 16:    FOR  A  PARTICULAR PURPOSE. See the GNU Lesser General Public  License  for 
 17:    more details.

 19:    You  should have received a copy of the GNU Lesser General  Public  License
 20:    along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
 21:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 22: */

 24:  #include private/ipimpl.h

 28: /*@
 29:    IPSetBilinearForm - Specifies the bilinear form to be used for
 30:    inner products.

 32:    Collective on IP

 34:    Input Parameters:
 35: +  ip    - the inner product context
 36: .  mat   - the matrix of the bilinear form (may be PETSC_NULL)
 37: -  form  - the type of bilinear form

 39:    Note:
 40:    This function is called by EPSSetProblemType() and usually need not be
 41:    called by the user.

 43:    Level: developer

 45: .seealso: IPGetBilinearForm(), IPInnerProduct(), IPNorm(), EPSSetProblemType(),
 46:           IPBilinearForm
 47: @*/
 48: PetscErrorCode IPSetBilinearForm(IP ip,Mat mat,IPBilinearForm form)
 49: {

 54:   if (mat) {
 56:     PetscObjectReference((PetscObject)mat);
 57:   }
 58:   if (ip->matrix) {
 59:     MatDestroy(ip->matrix);
 60:     VecDestroy(ip->Bx);
 61:   }
 62:   ip->matrix = mat;
 63:   ip->bilinear_form = form;
 64:   ip->xid = ip->xstate = 0;
 65:   if (mat) { MatGetVecs(mat,&ip->Bx,PETSC_NULL); }
 66:   return(0);
 67: }

 71: /*@C
 72:    IPGetBilinearForm - Retrieves the bilinear form to be used for
 73:    inner products.

 75:    Input Parameter:
 76: .  ip    - the inner product context

 78:    Output Parameter:
 79: +  mat   - the matrix of the bilinear form (may be PETSC_NULL)
 80: -  form  - the type of bilinear form

 82:    Level: developer

 84: .seealso: IPSetBilinearForm(), IPInnerProduct(), IPNorm(), EPSSetProblemType(),
 85:           IPBilinearForm
 86: @*/
 87: PetscErrorCode IPGetBilinearForm(IP ip,Mat* mat,IPBilinearForm* form)
 88: {
 91:   if (mat) *mat = ip->matrix;
 92:   if (form) *form = ip->bilinear_form;
 93:   return(0);
 94: }

 98: PetscErrorCode IPApplyMatrix_Private(IP ip,Vec x)
 99: {

103:   if (((PetscObject)x)->id != ip->xid || ((PetscObject)x)->state != ip->xstate) {
104:     PetscLogEventBegin(IP_ApplyMatrix,ip,0,0,0);
105:     MatMult(ip->matrix,x,ip->Bx);
106:     ip->xid = ((PetscObject)x)->id;
107:     ip->xstate = ((PetscObject)x)->state;
108:     PetscLogEventEnd(IP_ApplyMatrix,ip,0,0,0);
109:   }
110:   return(0);
111: }

115: /*@
116:    IPApplyMatrix - Multiplies a vector with the matrix associated to the
117:                    bilinear form.

119:    Collective on IP

121:    Input Parameters:
122: +  ip    - the inner product context
123: -  x     - the vector

125:    Output Parameter:
126: .  y     - the result  

128:    Note:
129:    If the bilinear form has no associated matrix this function copies the vector.

131:    Level: developer

133: .seealso: IPSetBilinearForm(), IPInnerProduct(), IPNorm(), EPSSetProblemType() 
134: @*/
135: PetscErrorCode IPApplyMatrix(IP ip,Vec x,Vec y)
136: {

141:   if (ip->matrix) {
142:     IPApplyMatrix_Private(ip,x);
143:     VecCopy(ip->Bx,y);
144:   } else {
145:     VecCopy(x,y);
146:   }
147:   return(0);
148: }