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-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/eps/epsimpl.h
17: /*
18: EPSAllocateSolution - Allocate memory storage for common variables such
19: as eigenvalues and eigenvectors.
20: */
21: PetscErrorCode EPSAllocateSolution(EPS eps)
22: {
24:
27: if (eps->allocated_ncv != eps->ncv) {
28: if (eps->allocated_ncv > 0) {
29: PetscFree(eps->eigr);
30: PetscFree(eps->eigi);
31: PetscFree(eps->errest);
32: PetscFree(eps->errest_left);
33: VecDestroyVecs(eps->V,eps->allocated_ncv);
34: VecDestroyVecs(eps->AV,eps->allocated_ncv);
35: if (eps->solverclass == EPS_TWO_SIDE) {
36: VecDestroyVecs(eps->W,eps->allocated_ncv);
37: VecDestroyVecs(eps->AW,eps->allocated_ncv);
38: }
39: }
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: VecDuplicateVecs(eps->vec_initial,eps->ncv,&eps->V);
45: VecDuplicateVecs(eps->vec_initial,eps->ncv,&eps->AV);
46: if (eps->solverclass == EPS_TWO_SIDE) {
47: VecDuplicateVecs(eps->vec_initial,eps->ncv,&eps->W);
48: VecDuplicateVecs(eps->vec_initial,eps->ncv,&eps->AW);
49: }
50: eps->allocated_ncv = eps->ncv;
51: }
52: return(0);
53: }
57: /*
58: EPSFreeSolution - Free memory storage. This routine is related to
59: EPSAllocateSolution().
60: */
61: PetscErrorCode EPSFreeSolution(EPS eps)
62: {
64:
67: if (eps->allocated_ncv > 0) {
68: PetscFree(eps->eigr);
69: PetscFree(eps->eigi);
70: PetscFree(eps->errest);
71: PetscFree(eps->errest_left);
72: VecDestroyVecs(eps->V,eps->allocated_ncv);
73: VecDestroyVecs(eps->AV,eps->allocated_ncv);
74: if (eps->solverclass == EPS_TWO_SIDE) {
75: VecDestroyVecs(eps->W,eps->allocated_ncv);
76: VecDestroyVecs(eps->AW,eps->allocated_ncv);
77: }
78: eps->allocated_ncv = 0;
79: }
80: return(0);
81: }
85: /*
86: EPSAllocateSolutionContiguous - Allocate memory storage for common
87: variables such as eigenvalues and eigenvectors. In this version, all
88: vectors in V (and AV) share a contiguous chunk of memory. This is
89: necessary for external packages such as Arpack.
90: */
91: PetscErrorCode EPSAllocateSolutionContiguous(EPS eps)
92: {
94: int i;
95: PetscInt nloc;
96: PetscScalar *pV,*pW;
99: if (eps->allocated_ncv != eps->ncv) {
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: VecGetArray(eps->AV[0],&pV);
112: for (i=0;i<eps->allocated_ncv;i++) {
113: VecDestroy(eps->AV[i]);
114: }
115: PetscFree(pV);
116: PetscFree(eps->AV);
117: if (eps->solverclass == EPS_TWO_SIDE) {
118: VecGetArray(eps->W[0],&pW);
119: for (i=0;i<eps->allocated_ncv;i++) {
120: VecDestroy(eps->W[i]);
121: }
122: PetscFree(pW);
123: PetscFree(eps->W);
124: VecGetArray(eps->AW[0],&pW);
125: for (i=0;i<eps->allocated_ncv;i++) {
126: VecDestroy(eps->AW[i]);
127: }
128: PetscFree(pW);
129: PetscFree(eps->AW);
130: }
131: }
132: PetscMalloc(eps->ncv*sizeof(PetscScalar),&eps->eigr);
133: PetscMalloc(eps->ncv*sizeof(PetscScalar),&eps->eigi);
134: PetscMalloc(eps->ncv*sizeof(PetscReal),&eps->errest);
135: PetscMalloc(eps->ncv*sizeof(PetscReal),&eps->errest_left);
136: VecGetLocalSize(eps->vec_initial,&nloc);
137: PetscMalloc(eps->ncv*sizeof(Vec),&eps->V);
138: PetscMalloc(eps->ncv*nloc*sizeof(PetscScalar),&pV);
139: for (i=0;i<eps->ncv;i++) {
140: VecCreateMPIWithArray(eps->comm,nloc,PETSC_DECIDE,pV+i*nloc,&eps->V[i]);
141: }
142: PetscMalloc(eps->ncv*sizeof(Vec),&eps->AV);
143: PetscMalloc(eps->ncv*nloc*sizeof(PetscScalar),&pV);
144: for (i=0;i<eps->ncv;i++) {
145: VecCreateMPIWithArray(eps->comm,nloc,PETSC_DECIDE,pV+i*nloc,&eps->AV[i]);
146: }
147: if (eps->solverclass == EPS_TWO_SIDE) {
148: PetscMalloc(eps->ncv*sizeof(Vec),&eps->W);
149: PetscMalloc(eps->ncv*nloc*sizeof(PetscScalar),&pW);
150: for (i=0;i<eps->ncv;i++) {
151: VecCreateMPIWithArray(eps->comm,nloc,PETSC_DECIDE,pW+i*nloc,&eps->W[i]);
152: }
153: PetscMalloc(eps->ncv*sizeof(Vec),&eps->AW);
154: PetscMalloc(eps->ncv*nloc*sizeof(PetscScalar),&pW);
155: for (i=0;i<eps->ncv;i++) {
156: VecCreateMPIWithArray(eps->comm,nloc,PETSC_DECIDE,pW+i*nloc,&eps->AW[i]);
157: }
158: }
159: eps->allocated_ncv = eps->ncv;
160: }
161: return(0);
162: }
166: /*
167: EPSFreeSolution - Free memory storage. This routine is related to
168: EPSAllocateSolutionContiguous().
169: */
170: PetscErrorCode EPSFreeSolutionContiguous(EPS eps)
171: {
173: int i;
174: PetscScalar *pV,*pW;
175:
178: if (eps->allocated_ncv > 0) {
179: PetscFree(eps->eigr);
180: PetscFree(eps->eigi);
181: PetscFree(eps->errest);
182: PetscFree(eps->errest_left);
183: VecGetArray(eps->V[0],&pV);
184: for (i=0;i<eps->allocated_ncv;i++) {
185: VecDestroy(eps->V[i]);
186: }
187: PetscFree(pV);
188: PetscFree(eps->V);
189: VecGetArray(eps->AV[0],&pV);
190: for (i=0;i<eps->allocated_ncv;i++) {
191: VecDestroy(eps->AV[i]);
192: }
193: PetscFree(pV);
194: PetscFree(eps->AV);
195: if (eps->solverclass == EPS_TWO_SIDE) {
196: VecGetArray(eps->W[0],&pW);
197: for (i=0;i<eps->allocated_ncv;i++) {
198: VecDestroy(eps->W[i]);
199: }
200: PetscFree(pW);
201: PetscFree(eps->W);
202: VecGetArray(eps->AW[0],&pW);
203: for (i=0;i<eps->allocated_ncv;i++) {
204: VecDestroy(eps->AW[i]);
205: }
206: PetscFree(pW);
207: PetscFree(eps->AW);
208: }
209: eps->allocated_ncv = 0;
210: }
211: return(0);
212: }