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-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/epsimpl.h
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: PetscInt i;
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->W) {
54: VecGetArray(eps->W[0],&pW);
55: for (i=0;i<eps->allocated_ncv;i++) {
56: VecDestroy(eps->W[i]);
57: }
58: PetscFree(pW);
59: PetscFree(eps->W);
60: }
61: }
62: PetscMalloc(eps->ncv*sizeof(PetscScalar),&eps->eigr);
63: PetscMalloc(eps->ncv*sizeof(PetscScalar),&eps->eigi);
64: PetscMalloc(eps->ncv*sizeof(PetscReal),&eps->errest);
65: PetscMalloc(eps->ncv*sizeof(PetscReal),&eps->errest_left);
66: PetscMalloc(eps->ncv*sizeof(Vec),&eps->V);
67: PetscMalloc(eps->ncv*eps->nloc*sizeof(PetscScalar),&pV);
68: for (i=0;i<eps->ncv;i++) {
69: VecCreateMPIWithArray(((PetscObject)eps)->comm,eps->nloc,PETSC_DECIDE,pV+i*eps->nloc,&eps->V[i]);
70: }
71: if (eps->leftvecs) {
72: PetscMalloc(eps->ncv*sizeof(Vec),&eps->W);
73: PetscMalloc(eps->ncv*eps->nloc*sizeof(PetscScalar),&pW);
74: for (i=0;i<eps->ncv;i++) {
75: VecCreateMPIWithArray(((PetscObject)eps)->comm,eps->nloc,PETSC_DECIDE,pW+i*eps->nloc,&eps->W[i]);
76: }
77: }
78: eps->allocated_ncv = eps->ncv;
79: }
80: return(0);
81: }
85: /*
86: EPSFreeSolution - Free memory storage. This routine is related to
87: EPSAllocateSolution().
88: */
89: PetscErrorCode EPSFreeSolution(EPS eps)
90: {
92: PetscInt i;
93: PetscScalar *pV,*pW;
94:
97: if (eps->allocated_ncv > 0) {
98: PetscFree(eps->eigr);
99: PetscFree(eps->eigi);
100: PetscFree(eps->errest);
101: PetscFree(eps->errest_left);
102: VecGetArray(eps->V[0],&pV);
103: for (i=0;i<eps->allocated_ncv;i++) {
104: VecDestroy(eps->V[i]);
105: }
106: PetscFree(pV);
107: PetscFree(eps->V);
108: if (eps->W) {
109: VecGetArray(eps->W[0],&pW);
110: for (i=0;i<eps->allocated_ncv;i++) {
111: VecDestroy(eps->W[i]);
112: }
113: PetscFree(pW);
114: PetscFree(eps->W);
115: }
116: eps->allocated_ncv = 0;
117: }
118: return(0);
119: }