#2634 zgeev does not operate on real matrices

[Dominik Wondrousch]

Within the lapack package, zgeev is the complex analogon to dgeev, which computes eigenvalues and ~vectors for real matrices.

Because real numbers should be a special type of complex number (i.e. by accident all imaginary parts are 0 ), zgeev should be able to deal with real valued matrices as well.

On the contrary, zgeev throws the error message "Error in IF [or a callee]: Wrong type error".

(%i58) S:matrix([1,0],[0,1]);
                                 [ 1  0 ]
(%o58)                           [      ]
                                 [ 0  1 ]
(%i59) dgeev(S,true,false);
                                 [ 1.0  0.0 ]
(%o59)              [[1.0, 1.0], [          ], false]
                                 [ 0.0  1.0 ]
(%i60) zgeev(S,true,false);
Maxima encountered a Lisp error:

 Error in IF [or a callee]: Wrong type error

Automatically continuing.
To enable the Lisp debugger set *debugger-hook* to nil.
(%i61) zgeev(%i*S,true,false);
                                  [ 1.0  0.0 ]
(%o61)         [[1.0 %i, 1.0 %i], [          ], false]
                                  [ 0.0  1.0 ]

Tested with Maxima 5.30.0 using GCL 2.6.8 under Windows 7 (standard installer)
Always reproducible.

Also, zgeev is NOT documented in the help file, but dgeev is. Found it's existance via google search in the feature requests.

[Raymond Toy]

I will update the documentation for zgeev. I will also look into adding zheev; that seems like a useful addition to zgeev.

I think zgeev should only accept complex-valued data. I suspect that dgeev would be faster and probably more accurate than using zgeev. If you disagree, please bring this up on the mailing list.

[Dominik Wondrousch]

I you don't mind, I'd like to disagree. It is always possible that by chance of computation or unknown user input a matrix, which is normally considered complex, turns out to be real instead. Taking your suggestion into account would result in somewhat tedious scripting: For each time one would instinctively call zgeev, one would have to check for calling dgeev or zgeev, respectively, depending wether all elements of imagpart are zero or not.

It may be true that dgeev should be prefered over zgeev for real matrices, but in case one does not know the nature of the argument it is natural to call for the most general case which is zgeev because all real numbers are complex numbers as well.

As a compromise, Maxima's zgeev could do a check for real matrices and call dgeev instead, but this may result in an unexpected deviation from LAPACK's standard behaviour.

[Raymond Toy]

Well, it was never really intended that dgeev and zgeev be the actual user interface to LAPACK. This interface was intentionally designed to be a relative bare interface to the LAPACK routines so you could just look at the LAPACK documentation and guess what to do.

The actual interface would have been something at a higher-level appropriate for maxima, possibly hiding, or selecting, the appropriate routine to use.

So far, no one has stepped up to design such an interface. :-(

I'm not completely opposed to making zgeev smarter, but I would really prefer that to be done somewhere else. Obviously, the pain of zgeev hasn't motivated anyone to make a better interface.

[Dominik Wondrousch]

I am sorry for sounding to be a bit nagging, but I still do not get the point.

When programming in C++ or FORTRAN and using LAPACK's routines, the programmer would definitely know how to treat the generated data - as real or as complex, so it is clear which to call: dgeev or zgeev.

In Maxima, the evaluation of a matrix may result in a complex or a real type. As in my case,

S: genmatrix(lambda([i,j],integrate(demoivre(conjugate(B[i])*B[j]),x,0,1)),length(B))

results in a complex hermitian matrix S with a user-provided set of functions B. But for Maxima, S could turn out to be real, given a special case of B. While programming in C++ or FORTRAN, one could simply treat S as being complex anyway because the compiler will not care for a+%i*b being in fact real for b=0. But Maxima will care. If the computation of S[i,j] is real for all i,j then S will be treated as real matrix and can't be used with zgeev anymore. And this incompatibility only depends on the (user-provided and hence unknown) choice of function array B.

I am totally fine with zgeev to be a bare interface to LAPACK. But if done so, an automatic conversion from real to complex should take place (although dgeev may give better results. If you insist, print a warning.) This is the equivalent of the C++ or FORTRAN case of storing real values in complex matrices by setting the imaginary part to zero and using (knowingly) the zgeev function. The knowledge, that the provided Maxima matrix is in fact real, is an artifact of how Maxima treats its equations, and should therefore be ignored in case explicit complex evaluation is requested.

[Robert Dodier]

Actually I think it is absolutely required that zgeev accept real numbers as input, since there is no way, at the user level, to construct a complex number which has zero imaginary part -- Maxima will simplify away any term 0*%i.

[Robert Dodier]

@Dominik: thanks for the suggestion about zheev, but in general I'm pretty sure it's best to separate different topics into different messages.

[Raymond Toy]

I still think it should be done on a different level. And while the user can't simply add 0.0*%i to force a complex value, the user can check if imagpart() is non-zero for some matrix element.

Are we going to make zgeev magically determine the matrix is hermitian and call zheev instead? Certainly possible, but I really think this should all be done at a different level instead of at the lapack routine level.

[Robert Dodier]

I've moved the feature request for zheev to feature requests as #108.

[Robert Dodier]

I've pushed commit 60be33 so zgeev always creates a complex matrix even if all elements have zero imaginary part. Closing this report as fixed.

[Dominik Wondrousch]

Thanks a lot. You are awesome !