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-2007, Universidad Politecnica de Valencia, Spain

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

 13: #include "src/ip/ipimpl.h"      /*I "slepcip.h" I*/

 17: /*@
 18:    IPSetBilinearForm - Specifies the bilinear form to be used for
 19:    inner products.

 21:    Collective on IP

 23:    Input Parameters:
 24: +  ip    - the inner product context
 25: .  mat   - the matrix of the bilinear form (may be PETSC_NULL)
 26: -  form  - the type of bilinear form

 28:    Note:
 29:    This function is called by EPSSetProblemType() and usually need not be
 30:    called by the user.

 32:    Level: developer

 34: .seealso: IPGetBilinearForm(), IPInnerProduct(), IPNorm(), EPSSetProblemType(),
 35:           IPBilinearForm
 36: @*/
 37: PetscErrorCode IPSetBilinearForm(IP ip,Mat mat,IPBilinearForm form)
 38: {

 43:   if (mat) {
 45:     PetscObjectReference((PetscObject)mat);
 46:   }
 47:   if (ip->matrix) {
 48:     MatDestroy(ip->matrix);
 49:     VecDestroy(ip->Bx);
 50:   }
 51:   ip->matrix = mat;
 52:   ip->bilinear_form = form;
 53:   ip->xid = ip->xstate = 0;
 54:   if (mat) { MatGetVecs(mat,&ip->Bx,PETSC_NULL); }
 55:   return(0);
 56: }

 60: /*@C
 61:    IPGetBilinearForm - Retrieves the bilinear form to be used for
 62:    inner products.

 64:    Input Parameter:
 65: .  ip    - the inner product context

 67:    Output Parameter:
 68: +  mat   - the matrix of the bilinear form (may be PETSC_NULL)
 69: -  form  - the type of bilinear form

 71:    Level: developer

 73: .seealso: IPSetBilinearForm(), IPInnerProduct(), IPNorm(), EPSSetProblemType(),
 74:           IPBilinearForm
 75: @*/
 76: PetscErrorCode IPGetBilinearForm(IP ip,Mat* mat,IPBilinearForm* form)
 77: {
 80:   if (mat) *mat = ip->matrix;
 81:   if (form) *form = ip->bilinear_form;
 82:   return(0);
 83: }

 87: PetscErrorCode IPApplyMatrix_Private(IP ip,Vec x)
 88: {

 92:   if (x->id != ip->xid || x->state != ip->xstate) {
 93:     PetscLogEventBegin(IP_ApplyMatrix,ip,0,0,0);
 94:     MatMult(ip->matrix,x,ip->Bx);
 95:     ip->xid = x->id;
 96:     ip->xstate = x->state;
 97:     PetscLogEventEnd(IP_ApplyMatrix,ip,0,0,0);
 98:   }
 99:   return(0);
100: }

104: /*@
105:    IPApplyMatrix - Multiplies a vector with the matrix associated to the
106:                    bilinear form.

108:    Collective on IP

110:    Input Parameters:
111: +  ip    - the inner product context
112: -  x     - the vector

114:    Output Parameter:
115: .  y     - the result  

117:    Note:
118:    If the bilinear form has no associated matrix this function copies the vector.

120:    Level: developer

122: .seealso: IPSetBilinearForm(), IPInnerProduct(), IPNorm(), EPSSetProblemType() 
123: @*/
124: PetscErrorCode IPApplyMatrix(IP ip,Vec x,Vec y)
125: {

130:   if (ip->matrix) {
131:     IPApplyMatrix_Private(ip,x);
132:     VecCopy(ip->Bx,y);
133:   } else {
134:     VecCopy(x,y);
135:   }
136:   return(0);
137: }