Actual source code: mem.c

  1: /*
  2:       EPS routines related to memory management.

  4:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  5:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  6:    Copyright (c) 2002-2009, 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/epsimpl.h

 28: /*
 29:   EPSAllocateSolution - Allocate memory storage for common variables such
 30:   as eigenvalues and eigenvectors. In this version, all
 31:   vectors in V (and W) share a contiguous chunk of memory.
 32: */
 33: PetscErrorCode EPSAllocateSolution(EPS eps)
 34: {
 36:   PetscInt       i,nloc;
 37:   PetscScalar    *pV,*pW;
 38: 
 41:   if (eps->allocated_ncv != eps->ncv) {
 42:     if (eps->allocated_ncv > 0) {
 43:       PetscFree(eps->eigr);
 44:       PetscFree(eps->eigi);
 45:       PetscFree(eps->errest);
 46:       PetscFree(eps->errest_left);
 47:       VecGetArray(eps->V[0],&pV);
 48:       for (i=0;i<eps->allocated_ncv;i++) {
 49:         VecDestroy(eps->V[i]);
 50:       }
 51:       PetscFree(pV);
 52:       PetscFree(eps->V);
 53:       if (eps->AV) { VecDestroyVecs(eps->AV,eps->allocated_ncv); }
 54:       if (eps->solverclass == EPS_TWO_SIDE) {
 55:         VecGetArray(eps->W[0],&pW);
 56:         for (i=0;i<eps->allocated_ncv;i++) {
 57:           VecDestroy(eps->W[i]);
 58:         }
 59:         PetscFree(pW);
 60:         PetscFree(eps->W);
 61:       }
 62:     }
 63:     PetscMalloc(eps->ncv*sizeof(PetscScalar),&eps->eigr);
 64:     PetscMalloc(eps->ncv*sizeof(PetscScalar),&eps->eigi);
 65:     PetscMalloc(eps->ncv*sizeof(PetscReal),&eps->errest);
 66:     PetscMalloc(eps->ncv*sizeof(PetscReal),&eps->errest_left);
 67:     VecGetLocalSize(eps->vec_initial,&nloc);
 68:     PetscMalloc(eps->ncv*sizeof(Vec),&eps->V);
 69:     PetscMalloc(eps->ncv*nloc*sizeof(PetscScalar),&pV);
 70:     for (i=0;i<eps->ncv;i++) {
 71:       VecCreateMPIWithArray(((PetscObject)eps)->comm,nloc,PETSC_DECIDE,pV+i*nloc,&eps->V[i]);
 72:     }
 73:     eps->AV = PETSC_NULL;
 74:     if (eps->solverclass == EPS_TWO_SIDE) {
 75:       PetscMalloc(eps->ncv*sizeof(Vec),&eps->W);
 76:       PetscMalloc(eps->ncv*nloc*sizeof(PetscScalar),&pW);
 77:       for (i=0;i<eps->ncv;i++) {
 78:         VecCreateMPIWithArray(((PetscObject)eps)->comm,nloc,PETSC_DECIDE,pW+i*nloc,&eps->W[i]);
 79:       }
 80:     }
 81:     eps->allocated_ncv = eps->ncv;
 82:   }
 83:   return(0);
 84: }

 88: /*
 89:   EPSFreeSolution - Free memory storage. This routine is related to 
 90:   EPSAllocateSolution().
 91: */
 92: PetscErrorCode EPSFreeSolution(EPS eps)
 93: {
 95:   PetscInt       i;
 96:   PetscScalar    *pV,*pW;
 97: 
100:   if (eps->allocated_ncv > 0) {
101:     PetscFree(eps->eigr);
102:     PetscFree(eps->eigi);
103:     PetscFree(eps->errest);
104:     PetscFree(eps->errest_left);
105:     VecGetArray(eps->V[0],&pV);
106:     for (i=0;i<eps->allocated_ncv;i++) {
107:       VecDestroy(eps->V[i]);
108:     }
109:     PetscFree(pV);
110:     PetscFree(eps->V);
111:     if (eps->AV) { VecDestroyVecs(eps->AV,eps->allocated_ncv); }
112:     if (eps->solverclass == EPS_TWO_SIDE) {
113:       VecGetArray(eps->W[0],&pW);
114:       for (i=0;i<eps->allocated_ncv;i++) {
115:         VecDestroy(eps->W[i]);
116:       }
117:       PetscFree(pW);
118:       PetscFree(eps->W);
119:     }
120:     eps->allocated_ncv = 0;
121:   }
122:   return(0);
123: }