Actual source code: mem.c
1: /*
2: EPS routines related to memory management.
3: */
4: #include src/eps/epsimpl.h
8: /*
9: EPSAllocateSolution - Allocate memory storage for common variables such
10: as eigenvalues and eigenvectors.
11: */
12: PetscErrorCode EPSAllocateSolution(EPS eps)
13: {
15:
18: if (eps->allocated_ncv != eps->ncv) {
19: if (eps->allocated_ncv > 0) {
20: PetscFree(eps->eigr);
21: PetscFree(eps->eigi);
22: PetscFree(eps->errest);
23: PetscFree(eps->errest_left);
24: VecDestroyVecs(eps->V,eps->allocated_ncv);
25: VecDestroyVecs(eps->AV,eps->allocated_ncv);
26: if (eps->solverclass == EPS_TWO_SIDE) {
27: VecDestroyVecs(eps->W,eps->allocated_ncv);
28: VecDestroyVecs(eps->AW,eps->allocated_ncv);
29: }
30: }
31: PetscMalloc(eps->ncv*sizeof(PetscScalar),&eps->eigr);
32: PetscMalloc(eps->ncv*sizeof(PetscScalar),&eps->eigi);
33: PetscMalloc(eps->ncv*sizeof(PetscReal),&eps->errest);
34: PetscMalloc(eps->ncv*sizeof(PetscReal),&eps->errest_left);
35: VecDuplicateVecs(eps->vec_initial,eps->ncv,&eps->V);
36: VecDuplicateVecs(eps->vec_initial,eps->ncv,&eps->AV);
37: if (eps->solverclass == EPS_TWO_SIDE) {
38: VecDuplicateVecs(eps->vec_initial,eps->ncv,&eps->W);
39: VecDuplicateVecs(eps->vec_initial,eps->ncv,&eps->AW);
40: }
41: eps->allocated_ncv = eps->ncv;
42: }
43: return(0);
44: }
48: /*
49: EPSFreeSolution - Free memory storage. This routine is related to
50: EPSAllocateSolution().
51: */
52: PetscErrorCode EPSFreeSolution(EPS eps)
53: {
55:
58: if (eps->allocated_ncv > 0) {
59: PetscFree(eps->eigr);
60: PetscFree(eps->eigi);
61: PetscFree(eps->errest);
62: PetscFree(eps->errest_left);
63: VecDestroyVecs(eps->V,eps->allocated_ncv);
64: VecDestroyVecs(eps->AV,eps->allocated_ncv);
65: if (eps->solverclass == EPS_TWO_SIDE) {
66: VecDestroyVecs(eps->W,eps->allocated_ncv);
67: VecDestroyVecs(eps->AW,eps->allocated_ncv);
68: }
69: eps->allocated_ncv = 0;
70: }
71: return(0);
72: }
76: /*
77: EPSAllocateSolutionContiguous - Allocate memory storage for common
78: variables such as eigenvalues and eigenvectors. In this version, all
79: vectors in V (and AV) share a contiguous chunk of memory. This is
80: necessary for external packages such as Arpack.
81: */
82: PetscErrorCode EPSAllocateSolutionContiguous(EPS eps)
83: {
85: int i;
86: PetscInt nloc;
87: PetscScalar *pV,*pW;
90: if (eps->allocated_ncv != eps->ncv) {
91: if (eps->allocated_ncv > 0) {
92: PetscFree(eps->eigr);
93: PetscFree(eps->eigi);
94: PetscFree(eps->errest);
95: PetscFree(eps->errest_left);
96: VecGetArray(eps->V[0],&pV);
97: for (i=0;i<eps->allocated_ncv;i++) {
98: VecDestroy(eps->V[i]);
99: }
100: PetscFree(pV);
101: PetscFree(eps->V);
102: VecGetArray(eps->AV[0],&pV);
103: for (i=0;i<eps->allocated_ncv;i++) {
104: VecDestroy(eps->AV[i]);
105: }
106: PetscFree(pV);
107: PetscFree(eps->AV);
108: if (eps->solverclass == EPS_TWO_SIDE) {
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: VecGetArray(eps->AW[0],&pW);
116: for (i=0;i<eps->allocated_ncv;i++) {
117: VecDestroy(eps->AW[i]);
118: }
119: PetscFree(pW);
120: PetscFree(eps->AW);
121: }
122: }
123: PetscMalloc(eps->ncv*sizeof(PetscScalar),&eps->eigr);
124: PetscMalloc(eps->ncv*sizeof(PetscScalar),&eps->eigi);
125: PetscMalloc(eps->ncv*sizeof(PetscReal),&eps->errest);
126: PetscMalloc(eps->ncv*sizeof(PetscReal),&eps->errest_left);
127: VecGetLocalSize(eps->vec_initial,&nloc);
128: PetscMalloc(eps->ncv*sizeof(Vec),&eps->V);
129: PetscMalloc(eps->ncv*nloc*sizeof(PetscScalar),&pV);
130: for (i=0;i<eps->ncv;i++) {
131: VecCreateMPIWithArray(eps->comm,nloc,PETSC_DECIDE,pV+i*nloc,&eps->V[i]);
132: }
133: PetscMalloc(eps->ncv*sizeof(Vec),&eps->AV);
134: PetscMalloc(eps->ncv*nloc*sizeof(PetscScalar),&pV);
135: for (i=0;i<eps->ncv;i++) {
136: VecCreateMPIWithArray(eps->comm,nloc,PETSC_DECIDE,pV+i*nloc,&eps->AV[i]);
137: }
138: if (eps->solverclass == EPS_TWO_SIDE) {
139: PetscMalloc(eps->ncv*sizeof(Vec),&eps->W);
140: PetscMalloc(eps->ncv*nloc*sizeof(PetscScalar),&pW);
141: for (i=0;i<eps->ncv;i++) {
142: VecCreateMPIWithArray(eps->comm,nloc,PETSC_DECIDE,pW+i*nloc,&eps->W[i]);
143: }
144: PetscMalloc(eps->ncv*sizeof(Vec),&eps->AW);
145: PetscMalloc(eps->ncv*nloc*sizeof(PetscScalar),&pW);
146: for (i=0;i<eps->ncv;i++) {
147: VecCreateMPIWithArray(eps->comm,nloc,PETSC_DECIDE,pW+i*nloc,&eps->AW[i]);
148: }
149: }
150: eps->allocated_ncv = eps->ncv;
151: }
152: return(0);
153: }
157: /*
158: EPSFreeSolution - Free memory storage. This routine is related to
159: EPSAllocateSolutionContiguous().
160: */
161: PetscErrorCode EPSFreeSolutionContiguous(EPS eps)
162: {
164: int i;
165: PetscScalar *pV,*pW;
166:
169: if (eps->allocated_ncv > 0) {
170: PetscFree(eps->eigr);
171: PetscFree(eps->eigi);
172: PetscFree(eps->errest);
173: PetscFree(eps->errest_left);
174: VecGetArray(eps->V[0],&pV);
175: for (i=0;i<eps->allocated_ncv;i++) {
176: VecDestroy(eps->V[i]);
177: }
178: PetscFree(pV);
179: PetscFree(eps->V);
180: VecGetArray(eps->AV[0],&pV);
181: for (i=0;i<eps->allocated_ncv;i++) {
182: VecDestroy(eps->AV[i]);
183: }
184: PetscFree(pV);
185: PetscFree(eps->AV);
186: if (eps->solverclass == EPS_TWO_SIDE) {
187: VecGetArray(eps->W[0],&pW);
188: for (i=0;i<eps->allocated_ncv;i++) {
189: VecDestroy(eps->W[i]);
190: }
191: PetscFree(pW);
192: PetscFree(eps->W);
193: VecGetArray(eps->AW[0],&pW);
194: for (i=0;i<eps->allocated_ncv;i++) {
195: VecDestroy(eps->AW[i]);
196: }
197: PetscFree(pW);
198: PetscFree(eps->AW);
199: }
200: eps->allocated_ncv = 0;
201: }
202: return(0);
203: }