# ARPACK gives different answers from Matlab and NAG

I'm playing with ARPACK. I looked into the examples they provide, zndrv4.f illustrating the usage of the routine znaupd, in the directory of ARPACK/EXAMPLES/COMPLEX/. I also came cross NAG Fortran Library. In NAG, there are some linear problem solvers F12***. The F12*** routines in NAG are equivalent to the znaupd in ARPACK. So I want to check if they will yield the same results.

I first looked into the example provided in the user guide of F12ARF at http://www.nag.co.uk/numeric/fl/nagdoc_fl22/pdf/F12/f12arf.pdf for example. In the end, it yields the results of

509.9390
380.9092
659.1558
271.9412


around the shift=500. I solved the same generalized eigenvalue problem in Matlab. Matlab gave me the same results.

But when I used znaupd from ARPACK to solve the same problem, I obtained different answers. The 4 eigenvalues are now

Eigenvalue               Residual
501.65650188259684       2.01604D-01
480.15153312181440       2.10853D-01
526.52596256924164       1.91850D-01
461.99019999608828       2.19336D-01


The routines in NAG and ARPACK both use SHIFTED INVERSE mode and solve the same generalized problem. I'm not sure what was wrong, but you can see the residual is too big (usually it should be around 1e-15). I attached my scripts of the ARPACK zndrv4.f (it's basically the same as the example file provided by ARPACK, I just need to change a little bit the matrices around line 174 to be the same as that in NAG.) and the Matlab file zndrv4.m in the compressed file below.

https://www.dropbox.com/s/mcq2jkmrgk2839u/zndrv4.zip?dl=0

Then later I found that actually the M matrix in F12ARF is normalized (entries divided by a six). I followed this and got the above wrong results in ARPACK. Now if I didn't divide the entries by six in zndrv4.f, the ARPACK script gives me the correct answer (the same as Matlab).

Eigenvalue               Residual
4.90161D+02              5.20906D-16
5.50553D+02              8.16401D-16
4.33527D+02              2.29794D-15
6.14759D+02              4.93630D-15


I'm confused by this. It seems to say that the ARPACK routine is not robust. Just scaling the matrix with a constant 6 will lead to such a big change in residual? Thank you very much.

• I'm slightly confused by this, and you probably understand this better, but exactly just what is this "normalization by six"? Why six, and normalized to what? Does this just mean "multiplied"? I googled and found this zndrv4.f, which divides by six in line 175..7, and also this file, which does not. The comment reads "symmetric tridiagonal matrix with 4.0 on the diagonal and 1.0 on the off-diagonals", so why divide that by 6? Do you know? – Kirill Dec 5 '14 at 9:59
• Hi, by normalization, I mean divided by 6, so the sum of a row is 1, otherwise is 6. I checked the second file from the package ARPACK_ng, where it seems funny. Here, the entries are divided by a six, if you run the driver as it is, things are ok. But if you change sigma = one at line 156 to sigma = (5.0D+2,0.0D0) and do not divide a six at lines 175-177, you will find that the residue is at 1e-1, which is wrong. I don't understand what goes wrong.... – jengmge Dec 9 '14 at 0:54