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

[Gong Ding]

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? 

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 );
  }