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-2011, Universitat 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> /*I "slepceps.h" I*/
28: /*
29: EPSAllocateSolution - Allocate memory storage for common variables such
30: as eigenvalues and eigenvectors. All vectors in V (and W) share a
31: contiguous chunk of memory.
32: */
33: PetscErrorCode EPSAllocateSolution(EPS eps)
34: {
36:
38: if (eps->allocated_ncv != eps->ncv) {
39: EPSFreeSolution(eps);
40: PetscMalloc(eps->ncv*sizeof(PetscScalar),&eps->eigr);
41: PetscMalloc(eps->ncv*sizeof(PetscScalar),&eps->eigi);
42: PetscMalloc(eps->ncv*sizeof(PetscReal),&eps->errest);
43: PetscMalloc(eps->ncv*sizeof(PetscReal),&eps->errest_left);
44: PetscMalloc(eps->ncv*sizeof(PetscInt),&eps->perm);
45: VecDuplicateVecs(eps->t,eps->ncv,&eps->V);
46: if (eps->leftvecs) {
47: VecDuplicateVecs(eps->t,eps->ncv,&eps->W);
48: }
49: eps->allocated_ncv = eps->ncv;
50: }
51: return(0);
52: }
56: /*
57: EPSFreeSolution - Free memory storage. This routine is related to
58: EPSAllocateSolution().
59: */
60: PetscErrorCode EPSFreeSolution(EPS eps)
61: {
63:
65: if (eps->allocated_ncv > 0) {
66: PetscFree(eps->eigr);
67: PetscFree(eps->eigi);
68: PetscFree(eps->errest);
69: PetscFree(eps->errest_left);
70: PetscFree(eps->perm);
71: VecDestroyVecs(eps->allocated_ncv,&eps->V);
72: VecDestroyVecs(eps->allocated_ncv,&eps->W);
73: eps->allocated_ncv = 0;
74: }
75: return(0);
76: }