| Line 50... |
Line 50... |
{
|
{
|
PetscErrorCode ierr;
|
PetscErrorCode ierr;
|
|
|
PetscFunctionBegin;
|
PetscFunctionBegin;
|
if (eps->ncv) { /* ncv set */
|
if (eps->ncv) { /* ncv set */
|
if (eps->ncv<eps->nev) SETERRQ(1,"The value of ncv must be at least nev");
|
if (eps->ncv<eps->nev) SETERRQ(((PetscObject)eps)->comm,1,"The value of ncv must be at least nev");
|
}
|
}
|
else if (eps->mpd) { /* mpd set */
|
else if (eps->mpd) { /* mpd set */
|
eps->ncv = PetscMin(eps->n,eps->nev+eps->mpd);
|
eps->ncv = PetscMin(eps->n,eps->nev+eps->mpd);
|
}
|
}
|
else { /* neither set: defaults depend on nev being small or large */
|
else { /* neither set: defaults depend on nev being small or large */
|
if (eps->nev<500) eps->ncv = PetscMin(eps->n,PetscMax(2*eps->nev,eps->nev+15));
|
if (eps->nev<500) eps->ncv = PetscMin(eps->n,PetscMax(2*eps->nev,eps->nev+15));
|
else { eps->mpd = 500; eps->ncv = PetscMin(eps->n,eps->nev+eps->mpd); }
|
else { eps->mpd = 500; eps->ncv = PetscMin(eps->n,eps->nev+eps->mpd); }
|
}
|
}
|
if (!eps->mpd) eps->mpd = eps->ncv;
|
if (!eps->mpd) eps->mpd = eps->ncv;
|
if (eps->ncv>eps->nev+eps->mpd) SETERRQ(1,"The value of ncv must not be larger than nev+mpd");
|
if (eps->ncv>eps->nev+eps->mpd) SETERRQ(((PetscObject)eps)->comm,1,"The value of ncv must not be larger than nev+mpd");
|
if (!eps->max_it) eps->max_it = PetscMax(100,2*eps->n/eps->ncv);
|
if (!eps->max_it) eps->max_it = PetscMax(100,2*eps->n/eps->ncv);
|
if (!eps->which) eps->which = EPS_LARGEST_MAGNITUDE;
|
if (!eps->which) eps->which = EPS_LARGEST_MAGNITUDE;
|
if (eps->ishermitian && (eps->which==EPS_LARGEST_IMAGINARY || eps->which==EPS_SMALLEST_IMAGINARY))
|
if (eps->ishermitian && (eps->which==EPS_LARGEST_IMAGINARY || eps->which==EPS_SMALLEST_IMAGINARY))
|
SETERRQ(1,"Wrong value of eps->which");
|
SETERRQ(((PetscObject)eps)->comm,1,"Wrong value of eps->which");
|
|
|
if (!eps->extraction) {
|
if (!eps->extraction) {
|
ierr = EPSSetExtraction(eps,EPS_RITZ);CHKERRQ(ierr);
|
ierr = EPSSetExtraction(eps,EPS_RITZ);CHKERRQ(ierr);
|
}
|
}
|
|
|
| Line 83... |
Line 83... |
} else {
|
} else {
|
ierr = EPSDefaultGetWork(eps,1);CHKERRQ(ierr);
|
ierr = EPSDefaultGetWork(eps,1);CHKERRQ(ierr);
|
}
|
}
|
|
|
/* dispatch solve method */
|
/* dispatch solve method */
|
if (eps->leftvecs) SETERRQ(PETSC_ERR_SUP,"Left vectors not supported in this solver");
|
if (eps->leftvecs) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_SUP,"Left vectors not supported in this solver");
|
eps->ops->solve = EPSSolve_ARNOLDI;
|
eps->ops->solve = EPSSolve_ARNOLDI;
|
PetscFunctionReturn(0);
|
PetscFunctionReturn(0);
|
}
|
}
|
|
|
#undef __FUNCT__
|
#undef __FUNCT__
|
| Line 300... |
Line 300... |
v is the resulting vector
|
v is the resulting vector
|
*/
|
*/
|
PetscErrorCode EPSUpdateVectors(EPS eps,PetscInt n_,Vec *U,PetscInt s,PetscInt e,PetscScalar *Q,PetscInt ldq,PetscScalar *H,PetscInt ldh_)
|
PetscErrorCode EPSUpdateVectors(EPS eps,PetscInt n_,Vec *U,PetscInt s,PetscInt e,PetscScalar *Q,PetscInt ldq,PetscScalar *H,PetscInt ldh_)
|
{
|
{
|
#if defined(PETSC_MISSING_LAPACK_GESVD)
|
#if defined(PETSC_MISSING_LAPACK_GESVD)
|
SETERRQ(PETSC_ERR_SUP,"GESVD - Lapack routine is unavailable.");
|
SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_SUP,"GESVD - Lapack routine is unavailable.");
|
#else
|
#else
|
PetscErrorCode ierr;
|
PetscErrorCode ierr;
|
PetscTruth isrefined;
|
PetscTruth isrefined;
|
PetscInt i,j,k;
|
PetscInt i,j,k;
|
PetscBLASInt n1,lwork,idummy=1,info,n=n_,ldh=ldh_;
|
PetscBLASInt n1,lwork,idummy=1,info,n=n_,ldh=ldh_;
|
| Line 336... |
Line 336... |
#if !defined(PETSC_USE_COMPLEX)
|
#if !defined(PETSC_USE_COMPLEX)
|
LAPACKgesvd_("N","O",&n1,&n,B,&n1,sigma,&sdummy,&idummy,&sdummy,&idummy,work,&lwork,&info);
|
LAPACKgesvd_("N","O",&n1,&n,B,&n1,sigma,&sdummy,&idummy,&sdummy,&idummy,work,&lwork,&info);
|
#else
|
#else
|
LAPACKgesvd_("N","O",&n1,&n,B,&n1,sigma,&sdummy,&idummy,&sdummy,&idummy,work,&lwork,sigma+n,&info);
|
LAPACKgesvd_("N","O",&n1,&n,B,&n1,sigma,&sdummy,&idummy,&sdummy,&idummy,work,&lwork,sigma+n,&info);
|
#endif
|
#endif
|
if (info) SETERRQ1(PETSC_ERR_LIB,"Error in Lapack xGESVD %d",info);
|
if (info) SETERRQ1(((PetscObject)eps)->comm,PETSC_ERR_LIB,"Error in Lapack xGESVD %d",info);
|
/* the smallest singular value is the new error estimate */
|
/* the smallest singular value is the new error estimate */
|
eps->errest[k] = sigma[n-1];
|
eps->errest[k] = sigma[n-1];
|
/* update vector with right singular vector associated to smallest singular value */
|
/* update vector with right singular vector associated to smallest singular value */
|
for (i=0;i<n;i++)
|
for (i=0;i<n;i++)
|
Q[k*ldq+i] = B[n-1+i*n1];
|
Q[k*ldq+i] = B[n-1+i*n1];
|
| Line 583... |
Line 583... |
EPS_ARNOLDI *arnoldi = (EPS_ARNOLDI *)eps->data;
|
EPS_ARNOLDI *arnoldi = (EPS_ARNOLDI *)eps->data;
|
|
|
PetscFunctionBegin;
|
PetscFunctionBegin;
|
ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);CHKERRQ(ierr);
|
ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);CHKERRQ(ierr);
|
if (!isascii) {
|
if (!isascii) {
|
SETERRQ1(1,"Viewer type %s not supported for EPSARNOLDI",((PetscObject)viewer)->type_name);
|
SETERRQ1(((PetscObject)eps)->comm,1,"Viewer type %s not supported for EPSARNOLDI",((PetscObject)viewer)->type_name);
|
}
|
}
|
if (arnoldi->delayed) {
|
if (arnoldi->delayed) {
|
ierr = PetscViewerASCIIPrintf(viewer,"using delayed reorthogonalization\n");CHKERRQ(ierr);
|
ierr = PetscViewerASCIIPrintf(viewer,"using delayed reorthogonalization\n");CHKERRQ(ierr);
|
}
|
}
|
PetscFunctionReturn(0);
|
PetscFunctionReturn(0);
|