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: {
14: PetscErrorCode ierr;
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: VecDestroyVecs(eps->V,eps->allocated_ncv);
24: VecDestroyVecs(eps->AV,eps->allocated_ncv);
25: }
26: PetscMalloc(eps->ncv*sizeof(PetscScalar),&eps->eigr);
27: PetscMalloc(eps->ncv*sizeof(PetscScalar),&eps->eigi);
28: PetscMalloc(eps->ncv*sizeof(PetscReal),&eps->errest);
29: VecDuplicateVecs(eps->vec_initial,eps->ncv,&eps->V);
30: VecDuplicateVecs(eps->vec_initial,eps->ncv,&eps->AV);
31: eps->allocated_ncv = eps->ncv;
32: }
33: return(0);
34: }
38: /*
39: EPSFreeSolution - Free memory storage. This routine is related to
40: EPSAllocateSolution().
41: */
42: PetscErrorCode EPSFreeSolution(EPS eps)
43: {
44: PetscErrorCode ierr;
45:
48: if (eps->allocated_ncv > 0) {
49: PetscFree(eps->eigr);
50: PetscFree(eps->eigi);
51: PetscFree(eps->errest);
52: VecDestroyVecs(eps->V,eps->allocated_ncv);
53: VecDestroyVecs(eps->AV,eps->allocated_ncv);
54: eps->allocated_ncv = 0;
55: }
56: return(0);
57: }
61: /*
62: EPSAllocateSolutionContiguous - Allocate memory storage for common
63: variables such as eigenvalues and eigenvectors. In this version, all
64: vectors in V (and AV) share a contiguous chunk of memory. This is
65: necessary for external packages such as Arpack.
66: */
67: PetscErrorCode EPSAllocateSolutionContiguous(EPS eps)
68: {
69: PetscErrorCode ierr;
70: int i,nloc;
71: PetscScalar *pV;
74: if (eps->allocated_ncv != eps->ncv) {
75: if (eps->allocated_ncv > 0) {
76: PetscFree(eps->eigr);
77: PetscFree(eps->eigi);
78: PetscFree(eps->errest);
79: VecGetArray(eps->V[0],&pV);
80: for (i=0;i<eps->allocated_ncv;i++) {
81: VecDestroy(eps->V[i]);
82: }
83: PetscFree(pV);
84: PetscFree(eps->V);
85: VecGetArray(eps->AV[0],&pV);
86: for (i=0;i<eps->allocated_ncv;i++) {
87: VecDestroy(eps->AV[i]);
88: }
89: PetscFree(pV);
90: PetscFree(eps->AV);
91: }
92: PetscMalloc(eps->ncv*sizeof(PetscScalar),&eps->eigr);
93: PetscMalloc(eps->ncv*sizeof(PetscScalar),&eps->eigi);
94: PetscMalloc(eps->ncv*sizeof(PetscReal),&eps->errest);
95: VecGetLocalSize(eps->vec_initial,&nloc);
96: PetscMalloc(eps->ncv*sizeof(Vec),&eps->V);
97: PetscMalloc(eps->ncv*nloc*sizeof(PetscScalar),&pV);
98: for (i=0;i<eps->ncv;i++) {
99: VecCreateMPIWithArray(eps->comm,nloc,PETSC_DECIDE,pV+i*nloc,&eps->V[i]);
100: }
101: PetscMalloc(eps->ncv*sizeof(Vec),&eps->AV);
102: PetscMalloc(eps->ncv*nloc*sizeof(PetscScalar),&pV);
103: for (i=0;i<eps->ncv;i++) {
104: VecCreateMPIWithArray(eps->comm,nloc,PETSC_DECIDE,pV+i*nloc,&eps->AV[i]);
105: }
106: eps->allocated_ncv = eps->ncv;
107: }
108: return(0);
109: }
113: /*
114: EPSFreeSolution - Free memory storage. This routine is related to
115: EPSAllocateSolutionContiguous().
116: */
117: PetscErrorCode EPSFreeSolutionContiguous(EPS eps)
118: {
119: PetscErrorCode ierr;
120: int i;
121: PetscScalar *pV;
122:
125: if (eps->allocated_ncv > 0) {
126: PetscFree(eps->eigr);
127: PetscFree(eps->eigi);
128: PetscFree(eps->errest);
129: VecGetArray(eps->V[0],&pV);
130: for (i=0;i<eps->allocated_ncv;i++) {
131: VecDestroy(eps->V[i]);
132: }
133: PetscFree(pV);
134: PetscFree(eps->V);
135: VecGetArray(eps->AV[0],&pV);
136: for (i=0;i<eps->allocated_ncv;i++) {
137: VecDestroy(eps->AV[i]);
138: }
139: PetscFree(pV);
140: PetscFree(eps->AV);
141: eps->allocated_ncv = 0;
142: }
143: return(0);
144: }