[petsc-users] Slepc eigen value solver gives strange behavior

[Jose E. Roman]

El 05/04/2011, a las 05:07, Gong Ding escribió:

> Hi,
> I use slepc eigen value solver to evaluate the eigen values of jacobian matrix on each nonlinear iteration.
> 
> First, I create the EPS structure and set the operator to jacobian matrix. And then call the EPSSolve in the SNES build 
> jacobian matrix functuion. This method create EPS once, call EPSSolve many times. However, it seems EPSSolve only work at
> first time, and result eigen value never changes (however, relative error becomes larger and larger) in the following solve procedure. 
> 
> Then I create EPS each time in the SNES build jacobian matrix functuion, do EPSSolve and delete EPS at the end of function.
> This method gives eigen value for the jacobian matrix with small relative error ~1e-8. Of course, create and destroy the EPS solver each time
> is not efficient.
> 
> Does something get wrong in the first method? 

You will always solve the same eigenproblem, unless you call EPSSetOperators every time the matrix changes.
When EPSSetOperators is called, the EPS object is reset and therefore EPSSetUp will be called in the next EPSSolve. So basically the first approach will have the same cost as the second approach.

By the way, you are not using STSINVERT correctly. You should set EPS_TARGET_MAGNITUDE (instead of EPS_SMALLEST_MAGNITUDE) together with target=0.0 with EPSSetTarget.

Jose

> 
> The code I used is attached here
> 
>  // create the EPS solver for smallest and largest eigen value
>  EPS            eps_s;
>  EPS            eps_l;
> 
>  FVM_NonlinearSolver & nonlinear_solver = dynamic_cast<FVM_NonlinearSolver &>(_solver);
>  Mat & J = nonlinear_solver.jacobian_matrix();
> 
>  // create eigen value problem solver
>  EPSCreate(PETSC_COMM_WORLD, &eps_s);
>  EPSCreate(PETSC_COMM_WORLD, &eps_l);
>  // Set operator
>  EPSSetOperators(eps_s, J, PETSC_NULL);
>  EPSSetOperators(eps_l, J, PETSC_NULL);
> 
>  // calculate smallest and largest eigen value
>  EPSSetWhichEigenpairs(eps_s, EPS_SMALLEST_MAGNITUDE);
>  EPSSetWhichEigenpairs(eps_l, EPS_LARGEST_MAGNITUDE);
> 
>  // shift-and-invert spectral transformation to enhance convergence of eigenvalues near zero
>  ST st_s;
>  EPSGetST(eps_s, &st_s);
>  STSetType(st_s, STSINVERT);
> 
> 
>  // Set solver parameters at runtime
>  EPSSetFromOptions(eps_s);
>  EPSSetFromOptions(eps_l);
> 
> 
> ///////////////////////////////////////////////////////////////////////////////
> 
>  // this part is called after jacobian matrix assemly
> 
>  PetscScalar kr_s, ki_s;
>  PetscScalar kr_l, ki_l;
>  PetscReal error_s;
>  PetscReal error_l;
>  PetscInt nconv_s;
>  PetscInt nconv_l;
> 
>  // get the smallest eigen value
>  EPSSolve( eps_s );
>  EPSGetConverged( eps_s, &nconv_s );
>  if( nconv_s > 0 )
>  {
>    EPSGetEigenvalue( eps_s, 0, &kr_s, &ki_s );
>    EPSComputeRelativeError( eps_s, 0, &error_s );
>  }
> 
>  // get the largest eigen value
>  EPSSolve( eps_l );
>  EPSGetConverged( eps_l, &nconv_l );
>  if( nconv_l > 0 )
>  {
>    EPSGetEigenvalue( eps_l, 0, &kr_l, &ki_l );
>    EPSComputeRelativeError( eps_l, 0, &error_l );
>  }
> 
>