| 1376 |
slepc |
1 |
/*
|
|
|
2 |
This file implements a wrapper to the BLZPACK package
|
| 6 |
dsic.upv.es!jroman |
3 |
|
| 1376 |
slepc |
4 |
- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
|
|
|
5 |
SLEPc - Scalable Library for Eigenvalue Problem Computations
|
|
|
6 |
Copyright (c) 2002-2007, Universidad Politecnica de Valencia, Spain
|
|
|
7 |
|
|
|
8 |
This file is part of SLEPc. See the README file for conditions of use
|
|
|
9 |
and additional information.
|
|
|
10 |
- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
|
| 6 |
dsic.upv.es!jroman |
11 |
*/
|
| 1376 |
slepc |
12 |
|
| 1521 |
slepc |
13 |
#include "src/eps/impls/blzpack/blzpackp.h"
|
| 6 |
dsic.upv.es!jroman |
14 |
|
| 114 |
dsic.upv.es!antodo |
15 |
const char* blzpack_error[33] = {
|
|
|
16 |
"",
|
|
|
17 |
"illegal data, LFLAG ",
|
|
|
18 |
"illegal data, dimension of (U), (V), (X) ",
|
|
|
19 |
"illegal data, leading dimension of (U), (V), (X) ",
|
|
|
20 |
"illegal data, leading dimension of (EIG) ",
|
|
|
21 |
"illegal data, number of required eigenpairs ",
|
|
|
22 |
"illegal data, Lanczos algorithm block size ",
|
|
|
23 |
"illegal data, maximum number of steps ",
|
|
|
24 |
"illegal data, number of starting vectors ",
|
|
|
25 |
"illegal data, number of eigenpairs provided ",
|
|
|
26 |
"illegal data, problem type flag ",
|
|
|
27 |
"illegal data, spectrum slicing flag ",
|
|
|
28 |
"illegal data, eigenvectors purification flag ",
|
|
|
29 |
"illegal data, level of output ",
|
|
|
30 |
"illegal data, output file unit ",
|
|
|
31 |
"illegal data, LCOMM (MPI or PVM) ",
|
|
|
32 |
"illegal data, dimension of ISTOR ",
|
|
|
33 |
"illegal data, convergence threshold ",
|
|
|
34 |
"illegal data, dimension of RSTOR ",
|
|
|
35 |
"illegal data on at least one PE ",
|
|
|
36 |
"ISTOR(3:14) must be equal on all PEs ",
|
|
|
37 |
"RSTOR(1:3) must be equal on all PEs ",
|
|
|
38 |
"not enough space in ISTOR to start eigensolution ",
|
|
|
39 |
"not enough space in RSTOR to start eigensolution ",
|
|
|
40 |
"illegal data, number of negative eigenvalues ",
|
|
|
41 |
"illegal data, entries of V ",
|
|
|
42 |
"illegal data, entries of X ",
|
|
|
43 |
"failure in computational subinterval ",
|
|
|
44 |
"file I/O error, blzpack.__.BQ ",
|
|
|
45 |
"file I/O error, blzpack.__.BX ",
|
|
|
46 |
"file I/O error, blzpack.__.Q ",
|
|
|
47 |
"file I/O error, blzpack.__.X ",
|
|
|
48 |
"parallel interface error "
|
|
|
49 |
};
|
|
|
50 |
|
| 6 |
dsic.upv.es!jroman |
51 |
#undef __FUNCT__
|
|
|
52 |
#define __FUNCT__ "EPSSetUp_BLZPACK"
|
| 476 |
dsic.upv.es!antodo |
53 |
PetscErrorCode EPSSetUp_BLZPACK(EPS eps)
|
| 6 |
dsic.upv.es!jroman |
54 |
{
|
| 476 |
dsic.upv.es!antodo |
55 |
PetscErrorCode ierr;
|
| 1089 |
slepc |
56 |
PetscInt N, n;
|
| 1509 |
slepc |
57 |
PetscBLASInt listor, lrstor, ncuv, k1, k2, k3, k4;
|
| 476 |
dsic.upv.es!antodo |
58 |
EPS_BLZPACK *blz = (EPS_BLZPACK *)eps->data;
|
|
|
59 |
PetscTruth flg;
|
|
|
60 |
KSP ksp;
|
|
|
61 |
PC pc;
|
| 6 |
dsic.upv.es!jroman |
62 |
|
|
|
63 |
PetscFunctionBegin;
|
| 1386 |
slepc |
64 |
ierr = VecGetSize(eps->vec_initial,&N);CHKERRQ(ierr);
|
|
|
65 |
ierr = VecGetLocalSize(eps->vec_initial,&n);CHKERRQ(ierr);
|
| 260 |
dsic.upv.es!antodo |
66 |
if (eps->ncv) {
|
|
|
67 |
if( eps->ncv < PetscMin(eps->nev+10,eps->nev*2) )
|
|
|
68 |
SETERRQ(0,"Warning: BLZpack recommends that ncv be larger than min(nev+10,nev*2)");
|
|
|
69 |
}
|
|
|
70 |
else eps->ncv = PetscMin(eps->nev+10,eps->nev*2);
|
| 1579 |
slepc |
71 |
if (eps->mpd) PetscInfo(eps,"Warning: parameter mpd ignored\n");
|
| 1220 |
slepc |
72 |
if (!eps->max_it) eps->max_it = PetscMax(1000,N);
|
| 260 |
dsic.upv.es!antodo |
73 |
|
| 6 |
dsic.upv.es!jroman |
74 |
if (!eps->ishermitian)
|
|
|
75 |
SETERRQ(PETSC_ERR_SUP,"Requested method is only available for Hermitian problems");
|
| 1229 |
slepc |
76 |
if (blz->slice || eps->isgeneralized) {
|
| 420 |
dsic.upv.es!antodo |
77 |
ierr = PetscTypeCompare((PetscObject)eps->OP,STSINV,&flg);CHKERRQ(ierr);
|
|
|
78 |
if (!flg)
|
| 1229 |
slepc |
79 |
SETERRQ(PETSC_ERR_SUP,"Shift-and-invert ST is needed for generalized problems or spectrum slicing");
|
| 420 |
dsic.upv.es!antodo |
80 |
ierr = STGetKSP(eps->OP,&ksp);CHKERRQ(ierr);
|
|
|
81 |
ierr = PetscTypeCompare((PetscObject)ksp,KSPPREONLY,&flg);CHKERRQ(ierr);
|
|
|
82 |
if (!flg)
|
| 1229 |
slepc |
83 |
SETERRQ(PETSC_ERR_SUP,"Preonly KSP is needed for generalized problems or spectrum slicing");
|
| 420 |
dsic.upv.es!antodo |
84 |
ierr = KSPGetPC(ksp,&pc);CHKERRQ(ierr);
|
|
|
85 |
ierr = PetscTypeCompare((PetscObject)pc,PCCHOLESKY,&flg);CHKERRQ(ierr);
|
|
|
86 |
if (!flg)
|
| 1229 |
slepc |
87 |
SETERRQ(PETSC_ERR_SUP,"Cholesky PC is needed for generalized problems or spectrum slicing");
|
| 420 |
dsic.upv.es!antodo |
88 |
}
|
| 92 |
dsic.upv.es!antodo |
89 |
if (eps->which!=EPS_SMALLEST_REAL)
|
| 6 |
dsic.upv.es!jroman |
90 |
SETERRQ(1,"Wrong value of eps->which");
|
|
|
91 |
|
|
|
92 |
k1 = PetscMin(N,180);
|
| 211 |
dsic.upv.es!antodo |
93 |
k2 = blz->block_size;
|
| 6 |
dsic.upv.es!jroman |
94 |
k4 = PetscMin(eps->ncv,N);
|
|
|
95 |
k3 = 484+k1*(13+k1*2+k2+PetscMax(18,k2+2))+k2*k2*3+k4*2;
|
|
|
96 |
|
|
|
97 |
listor = 123+k1*12;
|
| 1040 |
slepc |
98 |
ierr = PetscFree(blz->istor);CHKERRQ(ierr);
|
| 1509 |
slepc |
99 |
ierr = PetscMalloc((17+listor)*sizeof(PetscBLASInt),&blz->istor);CHKERRQ(ierr);
|
| 6 |
dsic.upv.es!jroman |
100 |
blz->istor[14] = listor;
|
|
|
101 |
|
| 420 |
dsic.upv.es!antodo |
102 |
if (blz->slice) lrstor = n*(k2*4+k1*2+k4)+k3;
|
| 6 |
dsic.upv.es!jroman |
103 |
else lrstor = n*(k2*4+k1)+k3;
|
| 1040 |
slepc |
104 |
ierr = PetscFree(blz->rstor);CHKERRQ(ierr);
|
| 6 |
dsic.upv.es!jroman |
105 |
ierr = PetscMalloc((4+lrstor)*sizeof(PetscReal),&blz->rstor);CHKERRQ(ierr);
|
|
|
106 |
blz->rstor[3] = lrstor;
|
|
|
107 |
|
|
|
108 |
ncuv = PetscMax(3,blz->block_size);
|
| 1040 |
slepc |
109 |
ierr = PetscFree(blz->u);CHKERRQ(ierr);
|
| 6 |
dsic.upv.es!jroman |
110 |
ierr = PetscMalloc(ncuv*n*sizeof(PetscScalar),&blz->u);CHKERRQ(ierr);
|
| 1040 |
slepc |
111 |
ierr = PetscFree(blz->v);CHKERRQ(ierr);
|
| 6 |
dsic.upv.es!jroman |
112 |
ierr = PetscMalloc(ncuv*n*sizeof(PetscScalar),&blz->v);CHKERRQ(ierr);
|
|
|
113 |
|
| 1040 |
slepc |
114 |
ierr = PetscFree(blz->eig);CHKERRQ(ierr);
|
| 6 |
dsic.upv.es!jroman |
115 |
ierr = PetscMalloc(2*eps->ncv*sizeof(PetscReal),&blz->eig);CHKERRQ(ierr);
|
|
|
116 |
|
| 1560 |
slepc |
117 |
if (eps->extraction) {
|
|
|
118 |
ierr = PetscInfo(eps,"Warning: extraction type ignored\n");CHKERRQ(ierr);
|
| 1426 |
slepc |
119 |
}
|
|
|
120 |
|
| 260 |
dsic.upv.es!antodo |
121 |
ierr = EPSAllocateSolutionContiguous(eps);CHKERRQ(ierr);
|
| 1345 |
slepc |
122 |
ierr = EPSDefaultGetWork(eps,1);CHKERRQ(ierr);
|
| 6 |
dsic.upv.es!jroman |
123 |
PetscFunctionReturn(0);
|
|
|
124 |
}
|
|
|
125 |
|
|
|
126 |
#undef __FUNCT__
|
|
|
127 |
#define __FUNCT__ "EPSSolve_BLZPACK"
|
| 476 |
dsic.upv.es!antodo |
128 |
PetscErrorCode EPSSolve_BLZPACK(EPS eps)
|
| 6 |
dsic.upv.es!jroman |
129 |
{
|
| 476 |
dsic.upv.es!antodo |
130 |
PetscErrorCode ierr;
|
| 1089 |
slepc |
131 |
EPS_BLZPACK *blz = (EPS_BLZPACK *)eps->data;
|
|
|
132 |
PetscInt n, nn;
|
| 1509 |
slepc |
133 |
PetscBLASInt i, nneig, lflag, nvopu;
|
| 476 |
dsic.upv.es!antodo |
134 |
Vec x, y;
|
|
|
135 |
PetscScalar sigma,*pV;
|
|
|
136 |
Mat A;
|
|
|
137 |
KSP ksp;
|
|
|
138 |
PC pc;
|
| 6 |
dsic.upv.es!jroman |
139 |
|
|
|
140 |
PetscFunctionBegin;
|
|
|
141 |
|
| 1386 |
slepc |
142 |
ierr = VecGetLocalSize(eps->vec_initial,&n); CHKERRQ(ierr);
|
| 1429 |
slepc |
143 |
ierr = VecCreateMPIWithArray(((PetscObject)eps)->comm,n,PETSC_DECIDE,PETSC_NULL,&x);CHKERRQ(ierr);
|
|
|
144 |
ierr = VecCreateMPIWithArray(((PetscObject)eps)->comm,n,PETSC_DECIDE,PETSC_NULL,&y);CHKERRQ(ierr);
|
| 6 |
dsic.upv.es!jroman |
145 |
ierr = VecGetArray(eps->V[0],&pV);CHKERRQ(ierr);
|
|
|
146 |
|
| 1229 |
slepc |
147 |
if (eps->isgeneralized && !blz->slice) {
|
|
|
148 |
ierr = STGetShift(eps->OP,&sigma);CHKERRQ(ierr); /* shift of origin */
|
|
|
149 |
blz->rstor[0] = sigma; /* lower limit of eigenvalue interval */
|
|
|
150 |
blz->rstor[1] = sigma; /* upper limit of eigenvalue interval */
|
|
|
151 |
} else {
|
|
|
152 |
sigma = 0.0;
|
|
|
153 |
blz->rstor[0] = blz->initial; /* lower limit of eigenvalue interval */
|
|
|
154 |
blz->rstor[1] = blz->final; /* upper limit of eigenvalue interval */
|
|
|
155 |
}
|
| 420 |
dsic.upv.es!antodo |
156 |
nneig = 0; /* no. of eigs less than sigma */
|
|
|
157 |
|
| 6 |
dsic.upv.es!jroman |
158 |
blz->istor[0] = n; /* number of rows of U, V, X*/
|
|
|
159 |
blz->istor[1] = n; /* leading dimension of U, V, X */
|
|
|
160 |
blz->istor[2] = eps->nev; /* number of required eigenpairs */
|
|
|
161 |
blz->istor[3] = eps->ncv; /* number of working eigenpairs */
|
|
|
162 |
blz->istor[4] = blz->block_size; /* number of vectors in a block */
|
| 422 |
dsic.upv.es!antodo |
163 |
blz->istor[5] = blz->nsteps; /* maximun number of steps per run */
|
| 6 |
dsic.upv.es!jroman |
164 |
blz->istor[6] = 1; /* number of starting vectors as input */
|
|
|
165 |
blz->istor[7] = 0; /* number of eigenpairs given as input */
|
| 1229 |
slepc |
166 |
blz->istor[8] = (blz->slice || eps->isgeneralized) ? 1 : 0; /* problem type */
|
| 420 |
dsic.upv.es!antodo |
167 |
blz->istor[9] = blz->slice; /* spectrum slicing */
|
| 1229 |
slepc |
168 |
blz->istor[10] = eps->isgeneralized ? 1 : 0; /* solutions refinement (purify) */
|
| 420 |
dsic.upv.es!antodo |
169 |
blz->istor[11] = 0; /* level of printing */
|
| 6 |
dsic.upv.es!jroman |
170 |
blz->istor[12] = 6; /* file unit for output */
|
| 1429 |
slepc |
171 |
blz->istor[13] = MPI_Comm_c2f(((PetscObject)eps)->comm); /* communicator */
|
| 6 |
dsic.upv.es!jroman |
172 |
|
|
|
173 |
blz->rstor[2] = eps->tol; /* threshold for convergence */
|
|
|
174 |
|
|
|
175 |
lflag = 0; /* reverse communication interface flag */
|
|
|
176 |
|
| 422 |
dsic.upv.es!antodo |
177 |
do {
|
| 6 |
dsic.upv.es!jroman |
178 |
|
|
|
179 |
BLZpack_( blz->istor, blz->rstor, &sigma, &nneig, blz->u, blz->v,
|
|
|
180 |
&lflag, &nvopu, blz->eig, pV );
|
|
|
181 |
|
| 422 |
dsic.upv.es!antodo |
182 |
switch (lflag) {
|
|
|
183 |
case 1:
|
| 6 |
dsic.upv.es!jroman |
184 |
/* compute v = OP u */
|
|
|
185 |
for (i=0;i<nvopu;i++) {
|
|
|
186 |
ierr = VecPlaceArray( x, blz->u+i*n );CHKERRQ(ierr);
|
|
|
187 |
ierr = VecPlaceArray( y, blz->v+i*n );CHKERRQ(ierr);
|
| 1229 |
slepc |
188 |
if (blz->slice || eps->isgeneralized) {
|
|
|
189 |
ierr = STAssociatedKSPSolve( eps->OP, x, y );CHKERRQ(ierr);
|
| 6 |
dsic.upv.es!jroman |
190 |
} else {
|
|
|
191 |
ierr = STApply( eps->OP, x, y ); CHKERRQ(ierr);
|
|
|
192 |
}
|
| 1538 |
slepc |
193 |
ierr = IPOrthogonalize(eps->ip,eps->nds,PETSC_NULL,eps->DS,y,PETSC_NULL,PETSC_NULL,PETSC_NULL,eps->work[0],PETSC_NULL);CHKERRQ(ierr);
|
| 850 |
dsic.upv.es!antodo |
194 |
ierr = VecResetArray(x);CHKERRQ(ierr);
|
|
|
195 |
ierr = VecResetArray(y);CHKERRQ(ierr);
|
| 6 |
dsic.upv.es!jroman |
196 |
}
|
| 422 |
dsic.upv.es!antodo |
197 |
/* monitor */
|
|
|
198 |
eps->nconv = BLZistorr_(blz->istor,"NTEIG",5);
|
|
|
199 |
EPSMonitor(eps,eps->its,eps->nconv,
|
|
|
200 |
blz->rstor+BLZistorr_(blz->istor,"IRITZ",5),
|
|
|
201 |
eps->eigi,
|
|
|
202 |
blz->rstor+BLZistorr_(blz->istor,"IRITZ",5)+BLZistorr_(blz->istor,"JT",2),
|
|
|
203 |
BLZistorr_(blz->istor,"NRITZ",5));
|
| 624 |
dsic.upv.es!antodo |
204 |
eps->its = eps->its + 1;
|
| 422 |
dsic.upv.es!antodo |
205 |
if (eps->its >= eps->max_it || eps->nconv >= eps->nev) lflag = 5;
|
|
|
206 |
break;
|
|
|
207 |
case 2:
|
| 6 |
dsic.upv.es!jroman |
208 |
/* compute v = B u */
|
|
|
209 |
for (i=0;i<nvopu;i++) {
|
|
|
210 |
ierr = VecPlaceArray( x, blz->u+i*n );CHKERRQ(ierr);
|
|
|
211 |
ierr = VecPlaceArray( y, blz->v+i*n );CHKERRQ(ierr);
|
| 1358 |
slepc |
212 |
ierr = IPApplyMatrix(eps->ip, x, y ); CHKERRQ(ierr);
|
| 850 |
dsic.upv.es!antodo |
213 |
ierr = VecResetArray(x);CHKERRQ(ierr);
|
|
|
214 |
ierr = VecResetArray(y);CHKERRQ(ierr);
|
| 6 |
dsic.upv.es!jroman |
215 |
}
|
| 422 |
dsic.upv.es!antodo |
216 |
break;
|
|
|
217 |
case 3:
|
| 6 |
dsic.upv.es!jroman |
218 |
/* update shift */
|
| 1229 |
slepc |
219 |
PetscInfo1(eps,"Factorization update (sigma=%g)\n",sigma);
|
| 420 |
dsic.upv.es!antodo |
220 |
ierr = STSetShift(eps->OP,sigma);CHKERRQ(ierr);
|
|
|
221 |
ierr = STGetKSP(eps->OP,&ksp);CHKERRQ(ierr);
|
|
|
222 |
ierr = KSPGetPC(ksp,&pc);CHKERRQ(ierr);
|
|
|
223 |
ierr = PCGetFactoredMatrix(pc,&A);CHKERRQ(ierr);
|
| 1089 |
slepc |
224 |
ierr = MatGetInertia(A,&nn,PETSC_NULL,PETSC_NULL);CHKERRQ(ierr);
|
|
|
225 |
nneig = nn;
|
| 422 |
dsic.upv.es!antodo |
226 |
break;
|
|
|
227 |
case 4:
|
| 6 |
dsic.upv.es!jroman |
228 |
/* copy the initial vector */
|
|
|
229 |
ierr = VecPlaceArray(x,blz->v);CHKERRQ(ierr);
|
| 1229 |
slepc |
230 |
ierr = EPSGetStartVector(eps,0,x,PETSC_NULL);CHKERRQ(ierr);
|
| 850 |
dsic.upv.es!antodo |
231 |
ierr = VecResetArray(x);CHKERRQ(ierr);
|
| 422 |
dsic.upv.es!antodo |
232 |
break;
|
| 6 |
dsic.upv.es!jroman |
233 |
}
|
| 422 |
dsic.upv.es!antodo |
234 |
|
|
|
235 |
} while (lflag > 0);
|
| 6 |
dsic.upv.es!jroman |
236 |
|
|
|
237 |
ierr = VecRestoreArray( eps->V[0], &pV ); CHKERRQ(ierr);
|
|
|
238 |
|
|
|
239 |
eps->nconv = BLZistorr_(blz->istor,"NTEIG",5);
|
|
|
240 |
eps->reason = EPS_CONVERGED_TOL;
|
|
|
241 |
|
|
|
242 |
for (i=0;i<eps->nconv;i++) {
|
|
|
243 |
eps->eigr[i]=blz->eig[i];
|
|
|
244 |
}
|
|
|
245 |
|
| 114 |
dsic.upv.es!antodo |
246 |
if (lflag!=0) {
|
|
|
247 |
char msg[2048] = "";
|
|
|
248 |
for (i = 0; i < 33; i++) {
|
| 134 |
dsic.upv.es!antodo |
249 |
if (blz->istor[15] & (1 << i)) PetscStrcat(msg, blzpack_error[i]);
|
| 114 |
dsic.upv.es!antodo |
250 |
}
|
|
|
251 |
SETERRQ2(PETSC_ERR_LIB,"Error in BLZPACK (code=%d): '%s'",blz->istor[15], msg);
|
|
|
252 |
}
|
| 6 |
dsic.upv.es!jroman |
253 |
ierr = VecDestroy(x);CHKERRQ(ierr);
|
|
|
254 |
ierr = VecDestroy(y);CHKERRQ(ierr);
|
|
|
255 |
|
|
|
256 |
PetscFunctionReturn(0);
|
|
|
257 |
}
|
|
|
258 |
|
|
|
259 |
#undef __FUNCT__
|
| 420 |
dsic.upv.es!antodo |
260 |
#define __FUNCT__ "EPSBackTransform_BLZPACK"
|
| 476 |
dsic.upv.es!antodo |
261 |
PetscErrorCode EPSBackTransform_BLZPACK(EPS eps)
|
| 420 |
dsic.upv.es!antodo |
262 |
{
|
| 476 |
dsic.upv.es!antodo |
263 |
PetscErrorCode ierr;
|
|
|
264 |
EPS_BLZPACK *blz = (EPS_BLZPACK *)eps->data;
|
| 420 |
dsic.upv.es!antodo |
265 |
|
|
|
266 |
PetscFunctionBegin;
|
| 1229 |
slepc |
267 |
if (!blz->slice && !eps->isgeneralized) {
|
| 442 |
dsic.upv.es!antodo |
268 |
ierr = EPSBackTransform_Default(eps);CHKERRQ(ierr);
|
|
|
269 |
}
|
| 420 |
dsic.upv.es!antodo |
270 |
PetscFunctionReturn(0);
|
|
|
271 |
}
|
|
|
272 |
|
|
|
273 |
#undef __FUNCT__
|
| 6 |
dsic.upv.es!jroman |
274 |
#define __FUNCT__ "EPSDestroy_BLZPACK"
|
| 476 |
dsic.upv.es!antodo |
275 |
PetscErrorCode EPSDestroy_BLZPACK(EPS eps)
|
| 6 |
dsic.upv.es!jroman |
276 |
{
|
| 476 |
dsic.upv.es!antodo |
277 |
PetscErrorCode ierr;
|
|
|
278 |
EPS_BLZPACK *blz = (EPS_BLZPACK *)eps->data;
|
| 6 |
dsic.upv.es!jroman |
279 |
|
|
|
280 |
PetscFunctionBegin;
|
| 25 |
dsic.upv.es!jroman |
281 |
PetscValidHeaderSpecific(eps,EPS_COOKIE,1);
|
| 1040 |
slepc |
282 |
ierr = PetscFree(blz->istor);CHKERRQ(ierr);
|
|
|
283 |
ierr = PetscFree(blz->rstor);CHKERRQ(ierr);
|
|
|
284 |
ierr = PetscFree(blz->u);CHKERRQ(ierr);
|
|
|
285 |
ierr = PetscFree(blz->v);CHKERRQ(ierr);
|
|
|
286 |
ierr = PetscFree(blz->eig);CHKERRQ(ierr);
|
|
|
287 |
ierr = PetscFree(eps->data);CHKERRQ(ierr);
|
| 260 |
dsic.upv.es!antodo |
288 |
ierr = EPSFreeSolutionContiguous(eps);CHKERRQ(ierr);
|
| 6 |
dsic.upv.es!jroman |
289 |
PetscFunctionReturn(0);
|
|
|
290 |
}
|
|
|
291 |
|
|
|
292 |
#undef __FUNCT__
|
|
|
293 |
#define __FUNCT__ "EPSView_BLZPACK"
|
| 476 |
dsic.upv.es!antodo |
294 |
PetscErrorCode EPSView_BLZPACK(EPS eps,PetscViewer viewer)
|
| 6 |
dsic.upv.es!jroman |
295 |
{
|
| 476 |
dsic.upv.es!antodo |
296 |
PetscErrorCode ierr;
|
|
|
297 |
EPS_BLZPACK *blz = (EPS_BLZPACK *) eps->data;
|
|
|
298 |
PetscTruth isascii;
|
| 6 |
dsic.upv.es!jroman |
299 |
|
|
|
300 |
PetscFunctionBegin;
|
|
|
301 |
ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);CHKERRQ(ierr);
|
|
|
302 |
if (!isascii) {
|
|
|
303 |
SETERRQ1(1,"Viewer type %s not supported for EPSBLZPACK",((PetscObject)viewer)->type_name);
|
|
|
304 |
}
|
|
|
305 |
ierr = PetscViewerASCIIPrintf(viewer,"block size of the block-Lanczos algorithm: %d\n",blz->block_size);CHKERRQ(ierr);
|
|
|
306 |
ierr = PetscViewerASCIIPrintf(viewer,"computational interval: [%f,%f]\n",blz->initial,blz->final);CHKERRQ(ierr);
|
|
|
307 |
PetscFunctionReturn(0);
|
|
|
308 |
}
|
|
|
309 |
|
|
|
310 |
#undef __FUNCT__
|
|
|
311 |
#define __FUNCT__ "EPSSetFromOptions_BLZPACK"
|
| 476 |
dsic.upv.es!antodo |
312 |
PetscErrorCode EPSSetFromOptions_BLZPACK(EPS eps)
|
| 6 |
dsic.upv.es!jroman |
313 |
{
|
| 476 |
dsic.upv.es!antodo |
314 |
PetscErrorCode ierr;
|
|
|
315 |
EPS_BLZPACK *blz = (EPS_BLZPACK *)eps->data;
|
| 1089 |
slepc |
316 |
PetscInt bs,n;
|
| 476 |
dsic.upv.es!antodo |
317 |
PetscReal interval[2];
|
|
|
318 |
PetscTruth flg;
|
|
|
319 |
KSP ksp;
|
|
|
320 |
PC pc;
|
| 6 |
dsic.upv.es!jroman |
321 |
|
|
|
322 |
PetscFunctionBegin;
|
|
|
323 |
ierr = PetscOptionsHead("BLZPACK options");CHKERRQ(ierr);
|
|
|
324 |
|
| 444 |
dsic.upv.es!antodo |
325 |
bs = blz->block_size;
|
|
|
326 |
ierr = PetscOptionsInt("-eps_blzpack_block_size","Block size","EPSBlzpackSetBlockSize",bs,&bs,&flg);CHKERRQ(ierr);
|
|
|
327 |
if (flg) {ierr = EPSBlzpackSetBlockSize(eps,bs);CHKERRQ(ierr);}
|
| 6 |
dsic.upv.es!jroman |
328 |
|
| 444 |
dsic.upv.es!antodo |
329 |
n = blz->nsteps;
|
|
|
330 |
ierr = PetscOptionsInt("-eps_blzpack_nsteps","Number of steps","EPSBlzpackSetNSteps",n,&n,&flg);CHKERRQ(ierr);
|
|
|
331 |
if (flg) {ierr = EPSBlzpackSetNSteps(eps,n);CHKERRQ(ierr);}
|
| 422 |
dsic.upv.es!antodo |
332 |
|
| 444 |
dsic.upv.es!antodo |
333 |
interval[0] = blz->initial;
|
|
|
334 |
interval[1] = blz->final;
|
|
|
335 |
n = 2;
|
|
|
336 |
ierr = PetscOptionsRealArray("-eps_blzpack_interval","Computational interval","EPSBlzpackSetInterval",interval,&n,&flg);CHKERRQ(ierr);
|
|
|
337 |
if (flg) {
|
|
|
338 |
if (n==1) interval[1]=interval[0];
|
|
|
339 |
ierr = EPSBlzpackSetInterval(eps,interval[0],interval[1]);CHKERRQ(ierr);
|
|
|
340 |
}
|
| 6 |
dsic.upv.es!jroman |
341 |
|
| 1229 |
slepc |
342 |
if (blz->slice || eps->isgeneralized) {
|
| 444 |
dsic.upv.es!antodo |
343 |
ierr = STSetType(eps->OP,STSINV);CHKERRQ(ierr);
|
|
|
344 |
ierr = STGetKSP(eps->OP,&ksp);CHKERRQ(ierr);
|
|
|
345 |
ierr = KSPSetType(ksp,KSPPREONLY);CHKERRQ(ierr);
|
|
|
346 |
ierr = KSPGetPC(ksp,&pc);CHKERRQ(ierr);
|
|
|
347 |
ierr = PCSetType(pc,PCCHOLESKY);CHKERRQ(ierr);
|
|
|
348 |
}
|
|
|
349 |
|
| 6 |
dsic.upv.es!jroman |
350 |
ierr = PetscOptionsTail();CHKERRQ(ierr);
|
|
|
351 |
PetscFunctionReturn(0);
|
|
|
352 |
}
|
|
|
353 |
|
|
|
354 |
EXTERN_C_BEGIN
|
|
|
355 |
#undef __FUNCT__
|
|
|
356 |
#define __FUNCT__ "EPSBlzpackSetBlockSize_BLZPACK"
|
| 1509 |
slepc |
357 |
PetscErrorCode EPSBlzpackSetBlockSize_BLZPACK(EPS eps,PetscInt bs)
|
| 6 |
dsic.upv.es!jroman |
358 |
{
|
| 476 |
dsic.upv.es!antodo |
359 |
EPS_BLZPACK *blz = (EPS_BLZPACK *) eps->data;;
|
| 6 |
dsic.upv.es!jroman |
360 |
|
|
|
361 |
PetscFunctionBegin;
|
| 211 |
dsic.upv.es!antodo |
362 |
if (bs == PETSC_DEFAULT) blz->block_size = 3;
|
|
|
363 |
else if (bs <= 0) {
|
|
|
364 |
SETERRQ(1, "Incorrect block size");
|
|
|
365 |
} else blz->block_size = bs;
|
| 6 |
dsic.upv.es!jroman |
366 |
PetscFunctionReturn(0);
|
|
|
367 |
}
|
|
|
368 |
EXTERN_C_END
|
|
|
369 |
|
|
|
370 |
#undef __FUNCT__
|
|
|
371 |
#define __FUNCT__ "EPSBlzpackSetBlockSize"
|
|
|
372 |
/*@
|
|
|
373 |
EPSBlzpackSetBlockSize - Sets the block size for the BLZPACK package.
|
|
|
374 |
|
|
|
375 |
Collective on EPS
|
|
|
376 |
|
|
|
377 |
Input Parameters:
|
|
|
378 |
+ eps - the eigenproblem solver context
|
|
|
379 |
- bs - block size
|
|
|
380 |
|
|
|
381 |
Options Database Key:
|
|
|
382 |
. -eps_blzpack_block_size - Sets the value of the block size
|
|
|
383 |
|
|
|
384 |
Level: advanced
|
|
|
385 |
|
| 553 |
dsic.upv.es!jroman |
386 |
.seealso: EPSBlzpackSetInterval()
|
| 6 |
dsic.upv.es!jroman |
387 |
@*/
|
| 1509 |
slepc |
388 |
PetscErrorCode EPSBlzpackSetBlockSize(EPS eps,PetscInt bs)
|
| 6 |
dsic.upv.es!jroman |
389 |
{
|
| 1509 |
slepc |
390 |
PetscErrorCode ierr, (*f)(EPS,PetscInt);
|
| 6 |
dsic.upv.es!jroman |
391 |
|
|
|
392 |
PetscFunctionBegin;
|
| 25 |
dsic.upv.es!jroman |
393 |
PetscValidHeaderSpecific(eps,EPS_COOKIE,1);
|
| 6 |
dsic.upv.es!jroman |
394 |
ierr = PetscObjectQueryFunction((PetscObject)eps,"EPSBlzpackSetBlockSize_C",(void (**)())&f);CHKERRQ(ierr);
|
|
|
395 |
if (f) {
|
|
|
396 |
ierr = (*f)(eps,bs);CHKERRQ(ierr);
|
|
|
397 |
}
|
|
|
398 |
PetscFunctionReturn(0);
|
|
|
399 |
}
|
|
|
400 |
|
|
|
401 |
EXTERN_C_BEGIN
|
|
|
402 |
#undef __FUNCT__
|
|
|
403 |
#define __FUNCT__ "EPSBlzpackSetInterval_BLZPACK"
|
| 476 |
dsic.upv.es!antodo |
404 |
PetscErrorCode EPSBlzpackSetInterval_BLZPACK(EPS eps,PetscReal initial,PetscReal final)
|
| 6 |
dsic.upv.es!jroman |
405 |
{
|
| 1229 |
slepc |
406 |
PetscErrorCode ierr;
|
| 476 |
dsic.upv.es!antodo |
407 |
EPS_BLZPACK *blz = (EPS_BLZPACK *) eps->data;;
|
| 6 |
dsic.upv.es!jroman |
408 |
|
|
|
409 |
PetscFunctionBegin;
|
|
|
410 |
blz->initial = initial;
|
|
|
411 |
blz->final = final;
|
| 420 |
dsic.upv.es!antodo |
412 |
blz->slice = 1;
|
| 1229 |
slepc |
413 |
ierr = STSetShift(eps->OP,initial);CHKERRQ(ierr);
|
| 6 |
dsic.upv.es!jroman |
414 |
PetscFunctionReturn(0);
|
|
|
415 |
}
|
|
|
416 |
EXTERN_C_END
|
|
|
417 |
|
|
|
418 |
#undef __FUNCT__
|
|
|
419 |
#define __FUNCT__ "EPSBlzpackSetInterval"
|
|
|
420 |
/*@
|
|
|
421 |
EPSBlzpackSetInterval - Sets the computational interval for the BLZPACK
|
|
|
422 |
package.
|
|
|
423 |
|
|
|
424 |
Collective on EPS
|
|
|
425 |
|
|
|
426 |
Input Parameters:
|
|
|
427 |
+ eps - the eigenproblem solver context
|
|
|
428 |
. initial - lower bound of the interval
|
|
|
429 |
- final - upper bound of the interval
|
|
|
430 |
|
|
|
431 |
Options Database Key:
|
|
|
432 |
. -eps_blzpack_interval - Sets the bounds of the interval (two values
|
|
|
433 |
separated by commas)
|
|
|
434 |
|
|
|
435 |
Note:
|
|
|
436 |
The following possibilities are accepted (see Blzpack user's guide for
|
|
|
437 |
details).
|
|
|
438 |
initial>final: start seeking for eigenpairs in the upper bound
|
|
|
439 |
initial<final: start in the lower bound
|
|
|
440 |
initial=final: run around a single value (no interval)
|
|
|
441 |
|
|
|
442 |
Level: advanced
|
|
|
443 |
|
|
|
444 |
.seealso: EPSBlzpackSetBlockSize()
|
|
|
445 |
@*/
|
| 476 |
dsic.upv.es!antodo |
446 |
PetscErrorCode EPSBlzpackSetInterval(EPS eps,PetscReal initial,PetscReal final)
|
| 6 |
dsic.upv.es!jroman |
447 |
{
|
| 476 |
dsic.upv.es!antodo |
448 |
PetscErrorCode ierr, (*f)(EPS,PetscReal,PetscReal);
|
| 6 |
dsic.upv.es!jroman |
449 |
|
|
|
450 |
PetscFunctionBegin;
|
| 25 |
dsic.upv.es!jroman |
451 |
PetscValidHeaderSpecific(eps,EPS_COOKIE,1);
|
| 6 |
dsic.upv.es!jroman |
452 |
ierr = PetscObjectQueryFunction((PetscObject)eps,"EPSBlzpackSetInterval_C",(void (**)())&f);CHKERRQ(ierr);
|
|
|
453 |
if (f) {
|
|
|
454 |
ierr = (*f)(eps,initial,final);CHKERRQ(ierr);
|
|
|
455 |
}
|
|
|
456 |
PetscFunctionReturn(0);
|
|
|
457 |
}
|
|
|
458 |
|
|
|
459 |
EXTERN_C_BEGIN
|
|
|
460 |
#undef __FUNCT__
|
| 422 |
dsic.upv.es!antodo |
461 |
#define __FUNCT__ "EPSBlzpackSetNSteps_BLZPACK"
|
| 1509 |
slepc |
462 |
PetscErrorCode EPSBlzpackSetNSteps_BLZPACK(EPS eps,PetscInt nsteps)
|
| 422 |
dsic.upv.es!antodo |
463 |
{
|
| 476 |
dsic.upv.es!antodo |
464 |
EPS_BLZPACK *blz = (EPS_BLZPACK *) eps->data;
|
| 422 |
dsic.upv.es!antodo |
465 |
|
|
|
466 |
PetscFunctionBegin;
|
|
|
467 |
blz->nsteps = nsteps == PETSC_DEFAULT ? 0 : nsteps;
|
|
|
468 |
PetscFunctionReturn(0);
|
|
|
469 |
}
|
|
|
470 |
EXTERN_C_END
|
|
|
471 |
|
|
|
472 |
#undef __FUNCT__
|
|
|
473 |
#define __FUNCT__ "EPSBlzpackSetNSteps"
|
|
|
474 |
/*@
|
|
|
475 |
EPSBlzpackSetNSteps - Sets the maximum number of steps per run for the BLZPACK
|
|
|
476 |
package.
|
|
|
477 |
|
|
|
478 |
Collective on EPS
|
|
|
479 |
|
|
|
480 |
Input Parameters:
|
|
|
481 |
+ eps - the eigenproblem solver context
|
|
|
482 |
- nsteps - maximum number of steps
|
|
|
483 |
|
|
|
484 |
Options Database Key:
|
| 424 |
dsic.upv.es!antodo |
485 |
. -eps_blzpack_nsteps - Sets the maximum number of steps per run
|
| 422 |
dsic.upv.es!antodo |
486 |
|
|
|
487 |
Level: advanced
|
|
|
488 |
|
|
|
489 |
@*/
|
| 1509 |
slepc |
490 |
PetscErrorCode EPSBlzpackSetNSteps(EPS eps,PetscInt nsteps)
|
| 422 |
dsic.upv.es!antodo |
491 |
{
|
| 1509 |
slepc |
492 |
PetscErrorCode ierr, (*f)(EPS,PetscInt);
|
| 422 |
dsic.upv.es!antodo |
493 |
|
|
|
494 |
PetscFunctionBegin;
|
|
|
495 |
PetscValidHeaderSpecific(eps,EPS_COOKIE,1);
|
|
|
496 |
ierr = PetscObjectQueryFunction((PetscObject)eps,"EPSBlzpackSetNSteps_C",(void (**)())&f);CHKERRQ(ierr);
|
|
|
497 |
if (f) {
|
|
|
498 |
ierr = (*f)(eps,nsteps);CHKERRQ(ierr);
|
|
|
499 |
}
|
|
|
500 |
PetscFunctionReturn(0);
|
|
|
501 |
}
|
|
|
502 |
|
|
|
503 |
EXTERN_C_BEGIN
|
|
|
504 |
#undef __FUNCT__
|
| 6 |
dsic.upv.es!jroman |
505 |
#define __FUNCT__ "EPSCreate_BLZPACK"
|
| 476 |
dsic.upv.es!antodo |
506 |
PetscErrorCode EPSCreate_BLZPACK(EPS eps)
|
| 6 |
dsic.upv.es!jroman |
507 |
{
|
| 476 |
dsic.upv.es!antodo |
508 |
PetscErrorCode ierr;
|
|
|
509 |
EPS_BLZPACK *blzpack;
|
| 6 |
dsic.upv.es!jroman |
510 |
|
|
|
511 |
PetscFunctionBegin;
|
|
|
512 |
ierr = PetscNew(EPS_BLZPACK,&blzpack);CHKERRQ(ierr);
|
|
|
513 |
PetscLogObjectMemory(eps,sizeof(EPS_BLZPACK));
|
|
|
514 |
eps->data = (void *) blzpack;
|
| 503 |
dsic.upv.es!antodo |
515 |
eps->ops->solve = EPSSolve_BLZPACK;
|
| 6 |
dsic.upv.es!jroman |
516 |
eps->ops->setup = EPSSetUp_BLZPACK;
|
| 503 |
dsic.upv.es!antodo |
517 |
eps->ops->setfromoptions = EPSSetFromOptions_BLZPACK;
|
| 6 |
dsic.upv.es!jroman |
518 |
eps->ops->destroy = EPSDestroy_BLZPACK;
|
|
|
519 |
eps->ops->view = EPSView_BLZPACK;
|
| 420 |
dsic.upv.es!antodo |
520 |
eps->ops->backtransform = EPSBackTransform_BLZPACK;
|
| 503 |
dsic.upv.es!antodo |
521 |
eps->ops->computevectors = EPSComputeVectors_Default;
|
| 6 |
dsic.upv.es!jroman |
522 |
|
| 211 |
dsic.upv.es!antodo |
523 |
blzpack->block_size = 3;
|
| 6 |
dsic.upv.es!jroman |
524 |
blzpack->initial = 0.0;
|
|
|
525 |
blzpack->final = 0.0;
|
| 420 |
dsic.upv.es!antodo |
526 |
blzpack->slice = 0;
|
| 422 |
dsic.upv.es!antodo |
527 |
blzpack->nsteps = 0;
|
|
|
528 |
|
|
|
529 |
ierr = PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSBlzpackSetBlockSize_C","EPSBlzpackSetBlockSize_BLZPACK",EPSBlzpackSetBlockSize_BLZPACK);CHKERRQ(ierr);
|
| 6 |
dsic.upv.es!jroman |
530 |
ierr = PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSBlzpackSetInterval_C","EPSBlzpackSetInterval_BLZPACK",EPSBlzpackSetInterval_BLZPACK);CHKERRQ(ierr);
|
| 422 |
dsic.upv.es!antodo |
531 |
ierr = PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSBlzpackSetNSteps_C","EPSBlzpackSetNSteps_BLZPACK",EPSBlzpackSetNSteps_BLZPACK);CHKERRQ(ierr);
|
| 6 |
dsic.upv.es!jroman |
532 |
|
| 104 |
dsic.upv.es!antodo |
533 |
eps->which = EPS_SMALLEST_REAL;
|
|
|
534 |
|
| 6 |
dsic.upv.es!jroman |
535 |
PetscFunctionReturn(0);
|
|
|
536 |
}
|
|
|
537 |
EXTERN_C_END
|