| Line 1... |
Line 1... |
/*
|
/*
|
|
|
SLEPc eigensolver: "krylovschur"
|
SLEPc eigensolver: "krylovschur"
|
|
|
Method: Krylov-Schur with spectrum slicing for symmetric eigenproblems
|
Method: Krylov-Schur with spectrum slicing for symmetric eigenproblems
|
|
|
| Line 91... |
Line 91... |
PC pc;
|
PC pc;
|
KSP ksp;
|
KSP ksp;
|
EPS_KRYLOVSCHUR *ctx = (EPS_KRYLOVSCHUR*)eps->data;
|
EPS_KRYLOVSCHUR *ctx = (EPS_KRYLOVSCHUR*)eps->data;
|
SR sr;
|
SR sr;
|
|
|
PetscFunctionBegin;
|
PetscFunctionBegin;
|
sr = ctx->sr;
|
sr = ctx->sr;
|
if (sr->nPend > 0) {
|
if (sr->nPend > 0) {
|
sr->sPrev = sr->sPres;
|
sr->sPrev = sr->sPres;
|
sr->sPres = sr->pending[--sr->nPend];
|
sr->sPres = sr->pending[--sr->nPend];
|
ierr = STSetShift(eps->st,sr->sPres->value);CHKERRQ(ierr);
|
ierr = STSetShift(eps->st,sr->sPres->value);CHKERRQ(ierr);
|
| Line 129... |
Line 129... |
A[k+i*ld] = sr->S[sr->nS+i];
|
A[k+i*ld] = sr->S[sr->nS+i];
|
}
|
}
|
sr->nS = k;
|
sr->nS = k;
|
ierr = DSRestoreArray(eps->ds,DS_MAT_A,&A);CHKERRQ(ierr);
|
ierr = DSRestoreArray(eps->ds,DS_MAT_A,&A);CHKERRQ(ierr);
|
ierr = DSSetDimensions(eps->ds,0,0,0,k);CHKERRQ(ierr);
|
ierr = DSSetDimensions(eps->ds,0,0,0,k);CHKERRQ(ierr);
|
/* Normalize u and append it to V */
|
/* Normalize u and append it to V */
|
ierr = VecAXPBY(eps->V[sr->nS],1.0/sr->beta,0.0,eps->work[0]);CHKERRQ(ierr);
|
ierr = VecAXPBY(eps->V[sr->nS],1.0/sr->beta,0.0,eps->work[0]);CHKERRQ(ierr);
|
}
|
}
|
eps->nconv = 0;
|
eps->nconv = 0;
|
} else sr->sPres = NULL;
|
} else sr->sPres = NULL;
|
PetscFunctionReturn(0);
|
PetscFunctionReturn(0);
|
| Line 174... |
Line 174... |
count0=0;count1=0; /* Found on both sides */
|
count0=0;count1=0; /* Found on both sides */
|
/* filling in values for the monitor */
|
/* filling in values for the monitor */
|
if (eps->numbermonitors >0) {
|
if (eps->numbermonitors >0) {
|
ierr = PetscObjectTypeCompare((PetscObject)eps->st,STCAYLEY,&iscayley);CHKERRQ(ierr);
|
ierr = PetscObjectTypeCompare((PetscObject)eps->st,STCAYLEY,&iscayley);CHKERRQ(ierr);
|
if (iscayley) {
|
if (iscayley) {
|
ierr = STCayleyGetAntishift(eps->st,&nu);CHKERRQ(ierr);
|
ierr = STCayleyGetAntishift(eps->st,&nu);CHKERRQ(ierr);
|
for (i=0;i<sr->indexEig;i++) {
|
for (i=0;i<sr->indexEig;i++) {
|
sr->monit[i]=(nu + sr->eig[i])/(sr->eig[i] - sPres->value);
|
sr->monit[i]=(nu + sr->eig[i])/(sr->eig[i] - sPres->value);
|
}
|
}
|
} else {
|
} else {
|
for (i=0;i<sr->indexEig;i++) {
|
for (i=0;i<sr->indexEig;i++) {
|
| Line 306... |
Line 306... |
complIterating = PETSC_TRUE;
|
complIterating = PETSC_TRUE;
|
if (n0 >sr->nMAXCompl)sch0 = PETSC_FALSE;
|
if (n0 >sr->nMAXCompl)sch0 = PETSC_FALSE;
|
if (n1 >sr->nMAXCompl)sch1 = PETSC_FALSE;
|
if (n1 >sr->nMAXCompl)sch1 = PETSC_FALSE;
|
iterCompl = sr->iterCompl;
|
iterCompl = sr->iterCompl;
|
} else eps->reason = EPS_CONVERGED_TOL;
|
} else eps->reason = EPS_CONVERGED_TOL;
|
}
|
}
|
}
|
}
|
/* Update l */
|
/* Update l */
|
if (eps->reason == EPS_CONVERGED_ITERATING) l = PetscMax(1,(PetscInt)((nv-k)*ctx->keep));
|
if (eps->reason == EPS_CONVERGED_ITERATING) l = PetscMax(1,(PetscInt)((nv-k)*ctx->keep));
|
else l = nv-k;
|
else l = nv-k;
|
if (breakdown) l=0;
|
if (breakdown) l=0;
|
| Line 382... |
Line 382... |
else sPres->nconv[0]++;
|
else sPres->nconv[0]++;
|
}
|
}
|
sPres->comp[0] = (count0 >= sPres->nsch[0])?PETSC_TRUE:PETSC_FALSE;
|
sPres->comp[0] = (count0 >= sPres->nsch[0])?PETSC_TRUE:PETSC_FALSE;
|
sPres->comp[1] = (count1 >= sPres->nsch[1])?PETSC_TRUE:PETSC_FALSE;
|
sPres->comp[1] = (count1 >= sPres->nsch[1])?PETSC_TRUE:PETSC_FALSE;
|
if (count0 > sPres->nsch[0] || count1 > sPres->nsch[1])SETERRQ(PetscObjectComm((PetscObject)eps),1,"Unexpected error in Spectrum Slicing!\nMismatch between number of values found and information from inertia");
|
if (count0 > sPres->nsch[0] || count1 > sPres->nsch[1])SETERRQ(PetscObjectComm((PetscObject)eps),1,"Unexpected error in Spectrum Slicing!\nMismatch between number of values found and information from inertia");
|
ierr = PetscFree(iwork);CHKERRQ(ierr);
|
ierr = PetscFree(iwork);CHKERRQ(ierr);
|
PetscFunctionReturn(0);
|
PetscFunctionReturn(0);
|
}
|
}
|
|
|
/*
|
/*
|
Obtains value of subsequent shift
|
Obtains value of subsequent shift
|
| Line 399... |
Line 399... |
PetscInt i,idxP;
|
PetscInt i,idxP;
|
SR sr;
|
SR sr;
|
shift sPres,s;
|
shift sPres,s;
|
EPS_KRYLOVSCHUR *ctx = (EPS_KRYLOVSCHUR*)eps->data;
|
EPS_KRYLOVSCHUR *ctx = (EPS_KRYLOVSCHUR*)eps->data;
|
|
|
PetscFunctionBegin;
|
PetscFunctionBegin;
|
sr = ctx->sr;
|
sr = ctx->sr;
|
sPres = sr->sPres;
|
sPres = sr->sPres;
|
if (sPres->neighb[side]) {
|
if (sPres->neighb[side]) {
|
/* Completing a previous interval */
|
/* Completing a previous interval */
|
if (!sPres->neighb[side]->neighb[side] && sPres->neighb[side]->nconv[side]==0) { /* One of the ends might be too far from eigenvalues */
|
if (!sPres->neighb[side]->neighb[side] && sPres->neighb[side]->nconv[side]==0) { /* One of the ends might be too far from eigenvalues */
|
| Line 415... |
Line 415... |
if (sPres->neighb[0]) {
|
if (sPres->neighb[0]) {
|
/* Multiplying by 10 the previous distance */
|
/* Multiplying by 10 the previous distance */
|
*newS = sPres->value + 10*(sr->dir)*PetscAbsReal(sPres->value - sPres->neighb[0]->value);
|
*newS = sPres->value + 10*(sr->dir)*PetscAbsReal(sPres->value - sPres->neighb[0]->value);
|
sr->nleap++;
|
sr->nleap++;
|
/* Stops when the interval is open and no values are found in the last 5 shifts (there might be infinite eigenvalues) */
|
/* Stops when the interval is open and no values are found in the last 5 shifts (there might be infinite eigenvalues) */
|
if (!sr->hasEnd && sr->nleap > 5) SETERRQ(PetscObjectComm((PetscObject)eps),1,"Unable to compute the wanted eigenvalues with open interval");
|
if (!sr->hasEnd && sr->nleap > 5) SETERRQ(PetscObjectComm((PetscObject)eps),1,"Unable to compute the wanted eigenvalues with open interval");
|
} else { /* First shift */
|
} else { /* First shift */
|
if (eps->nconv != 0) {
|
if (eps->nconv != 0) {
|
/* Unaccepted values give information for next shift */
|
/* Unaccepted values give information for next shift */
|
idxP=0;/* Number of values left from shift */
|
idxP=0;/* Number of values left from shift */
|
for (i=0;i<eps->nconv;i++) {
|
for (i=0;i<eps->nconv;i++) {
|
| Line 450... |
Line 450... |
} else { /* First shift. Average distance obtained with values in this shift */
|
} else { /* First shift. Average distance obtained with values in this shift */
|
/* first shift might be too far from first wanted eigenvalue (no values found outside the interval)*/
|
/* first shift might be too far from first wanted eigenvalue (no values found outside the interval)*/
|
if ((sr->dir)*(PetscRealPart(sr->eig[0])-sPres->value)>0 && PetscAbsReal((PetscRealPart(sr->eig[sr->indexEig-1]) - PetscRealPart(sr->eig[0]))/PetscRealPart(sr->eig[0])) > PetscSqrtReal(eps->tol)) {
|
if ((sr->dir)*(PetscRealPart(sr->eig[0])-sPres->value)>0 && PetscAbsReal((PetscRealPart(sr->eig[sr->indexEig-1]) - PetscRealPart(sr->eig[0]))/PetscRealPart(sr->eig[0])) > PetscSqrtReal(eps->tol)) {
|
d_prev = PetscAbsReal((PetscRealPart(sr->eig[sr->indexEig-1]) - PetscRealPart(sr->eig[0])))/(sPres->neigs+0.3);
|
d_prev = PetscAbsReal((PetscRealPart(sr->eig[sr->indexEig-1]) - PetscRealPart(sr->eig[0])))/(sPres->neigs+0.3);
|
} else {
|
} else {
|
d_prev = PetscAbsReal(PetscRealPart(sr->eig[sr->indexEig-1]) - sPres->value)/(sPres->neigs+0.3);
|
d_prev = PetscAbsReal(PetscRealPart(sr->eig[sr->indexEig-1]) - sPres->value)/(sPres->neigs+0.3);
|
}
|
}
|
}
|
}
|
/* Average distance is used for next shift by adding it to value on the right or to shift */
|
/* Average distance is used for next shift by adding it to value on the right or to shift */
|
if ((sr->dir)*(PetscRealPart(sr->eig[sPres->index + sPres->neigs -1]) - sPres->value)>0) {
|
if ((sr->dir)*(PetscRealPart(sr->eig[sPres->index + sPres->neigs -1]) - sPres->value)>0) {
|
*newS = PetscRealPart(sr->eig[sPres->index + sPres->neigs -1])+ ((sr->dir)*d_prev*(eps->nev))/2;
|
*newS = PetscRealPart(sr->eig[sPres->index + sPres->neigs -1])+ ((sr->dir)*d_prev*(eps->nev))/2;
|
} else { /* Last accepted value is on the left of shift. Adding to shift */
|
} else { /* Last accepted value is on the left of shift. Adding to shift */
|
*newS = sPres->value + ((sr->dir)*d_prev*(eps->nev))/2;
|
*newS = sPres->value + ((sr->dir)*d_prev*(eps->nev))/2;
|
}
|
}
|
}
|
}
|
/* End of interval can not be surpassed */
|
/* End of interval can not be surpassed */
|
| Line 600... |
Line 600... |
PetscErrorCode ierr;
|
PetscErrorCode ierr;
|
PetscInt i,lds;
|
PetscInt i,lds;
|
PetscReal newS;
|
PetscReal newS;
|
KSP ksp;
|
KSP ksp;
|
PC pc;
|
PC pc;
|
Mat F;
|
Mat F;
|
PetscReal *errest_left;
|
PetscReal *errest_left;
|
Vec t;
|
Vec t;
|
SR sr;
|
SR sr;
|
shift s;
|
shift s;
|
EPS_KRYLOVSCHUR *ctx = (EPS_KRYLOVSCHUR*)eps->data;
|
EPS_KRYLOVSCHUR *ctx = (EPS_KRYLOVSCHUR*)eps->data;
|