Subversion Repositories slepc-dev

Rev

Go to most recent revision | Details | Compare with Previous | Last modification | View Log | RSS feed

Rev Author Line No. Line
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