[Octave-bug-tracker] [bug #53719] eigs delivers wrong result in case of complex valued general eigenvalue problem

[Dan Sebald]

Follow-up Comment #9, bug #53719 (project octave):

@Marco: I've looked at the code just a little bit.  I have no concrete
answers, but just observations.

First, I agree with your thinking that the FORTRAN (non-sparse construct) and
C (sparse construct) implementations should behave pretty much the same in
terms of convergences.  So, this is a good test.

For the FORTRAN version, I had to remove the debug header file 


zndrv4_53719.f:141: Error: Can't open included file 'debug.h'


I then run the program and get some loop outputs for the residuals, then the
final summary:


 _NDRV4 
 ====== 
  
  Size of the matrix is           10
  The number of Ritz values requested is            4
  The number of Arnoldi vectors generated (NCV) is           10
  What portion of the spectrum: SM
  The number of converged Ritz values is    538187800
  The number of Implicit Arnoldi update iterations taken is            2
  The number of OP*x is           17
  The convergence criterion is    1.1102230246251565E-016


Everything looks reasonable there except the number of converged Ritz values
reported is unrealistic.  It should report 0.  log2(538187800) is 29.004, so
it's not like this could be -1 in disguise or anything.

The other thing that has me wondering is the IDO value reported keeps
repeating with the pattern of 2,2,1,2,2,1,2,2,1,2,2,1...


 IDO =           -1
 IDO =            2
 IDO =            1
 IDO =            2
 IDO =            2
 IDO =            1
 IDO =            2
 IDO =            2
 IDO =            1
 IDO =            2
 IDO =            2
 IDO =            1
 IDO =            2
 IDO =            2
 IDO =            1
 IDO =            2
 IDO =            2
 IDO =            1
 IDO =            2
 IDO =            2
 IDO =            1
 IDO =            2
 IDO =            2
 IDO =            1
 IDO =            2
 IDO =            2
 IDO =            1
 IDO =            2
 IDO =            2
 IDO =            1
 IDO =            2
 IDO =            2
 IDO =            2
 IDO =            1
 IDO =            2
 IDO =            2
 IDO =            1
 IDO =            2
 IDO =            2
 IDO =            1
 IDO =            2
 IDO =            2
 IDO =            1
 IDO =            2
 IDO =            2
 IDO =            1
 IDO =            2
 IDO =            2
 IDO =            1
 IDO =            2
 IDO =            2


The above does look to be about 17 repeated patterns (i.e., the number of
computed OP*x).  However, I'm wondering why in your code for the IDO=-1 and
IDO=1 cases there is the division by something resembling the eigenvalues or
Ritz values (is it supposed to be (A-sigma*M)^-1?).  Won't that throw the
estimates out of sorts?

>From documentation

http://www.caam.rice.edu/software/ARPACK/UG/node138.html

when IDO is 1, the routine wants


c          IDO =  1: compute  Y = OP * Z  and Z = B * X where


and what the code is computing is

 
y <-- OP*x = inv[A-sigma*M]*M*x |


Otherwise, all the FORTRAN vs. C indexing looks about right.  The mx() routine
looks to match what the eigs() example is set up as.

    _______________________________________________________

Reply to this item at:

  <http://savannah.gnu.org/bugs/?53719>

_______________________________________________
  Message sent via/by Savannah
  http://savannah.gnu.org/