| Line 180... |
Line 180... |
ierr = MatMult(A,v,e);CHKERRQ(ierr);
|
ierr = MatMult(A,v,e);CHKERRQ(ierr);
|
ierr = VecDot(v,e,&alpha2);CHKERRQ(ierr);
|
ierr = VecDot(v,e,&alpha2);CHKERRQ(ierr);
|
alpha2 = alpha2 / (beta1 * beta1);
|
alpha2 = alpha2 / (beta1 * beta1);
|
|
|
/* choose the eigenvalue of [rho beta1; beta1 alpha2] closest to rho */
|
/* choose the eigenvalue of [rho beta1; beta1 alpha2] closest to rho */
|
|
ierr = PetscFPTrapPush(PETSC_FP_TRAP_OFF);CHKERRQ(ierr);
|
LAPACKlaev2_(&rho,&beta1,&alpha2,&rt1,&rt2,&cs1,&sn1);
|
LAPACKlaev2_(&rho,&beta1,&alpha2,&rt1,&rt2,&cs1,&sn1);
|
|
ierr = PetscFPTrapPop();CHKERRQ(ierr);
|
if (PetscAbsScalar(rt1-rho) < PetscAbsScalar(rt2-rho)) rho = rt1;
|
if (PetscAbsScalar(rt1-rho) < PetscAbsScalar(rt2-rho)) rho = rt1;
|
else rho = rt2;
|
else rho = rt2;
|
}
|
}
|
/* update operator according to new shift */
|
/* update operator according to new shift */
|
PetscPushErrorHandler(PetscIgnoreErrorHandler,PETSC_NULL);
|
PetscPushErrorHandler(PetscIgnoreErrorHandler,PETSC_NULL);
|
| Line 325... |
Line 327... |
ierr = MatMult(A,v,e);CHKERRQ(ierr);
|
ierr = MatMult(A,v,e);CHKERRQ(ierr);
|
ierr = VecDot(v,e,&alpha2);CHKERRQ(ierr);
|
ierr = VecDot(v,e,&alpha2);CHKERRQ(ierr);
|
alpha2 = alpha2 / (beta1 * beta1);
|
alpha2 = alpha2 / (beta1 * beta1);
|
|
|
/* choose the eigenvalue of [rho beta1; beta1 alpha2] closest to rho */
|
/* choose the eigenvalue of [rho beta1; beta1 alpha2] closest to rho */
|
|
ierr = PetscFPTrapPush(PETSC_FP_TRAP_OFF);CHKERRQ(ierr);
|
LAPACKlaev2_(&rho,&beta1,&alpha2,&rt1,&rt2,&cs1,&sn1);
|
LAPACKlaev2_(&rho,&beta1,&alpha2,&rt1,&rt2,&cs1,&sn1);
|
|
ierr = PetscFPTrapPop();CHKERRQ(ierr);
|
if (PetscAbsScalar(rt1-rho) < PetscAbsScalar(rt2-rho)) rho = rt1;
|
if (PetscAbsScalar(rt1-rho) < PetscAbsScalar(rt2-rho)) rho = rt1;
|
else rho = rt2;
|
else rho = rt2;
|
}
|
}
|
/* update operator according to new shift */
|
/* update operator according to new shift */
|
PetscPushErrorHandler(PetscIgnoreErrorHandler,PETSC_NULL);
|
PetscPushErrorHandler(PetscIgnoreErrorHandler,PETSC_NULL);
|