0

How can I tell Eigensystem that a matrix $M$, which I would like to diagonalize, is a numerical matrix of complex numbers? My idea is that this information could speed up the calculation, since Mathematica would not need to determine the type of the matrix.

m_goldberg
  • 107,779
  • 16
  • 103
  • 257
lagoa
  • 835
  • 1
  • 7
  • 10
  • It's not clear to me what you mean. Are you trying to find a symbolic solution first by using assumptions, or do you want to define a function that takes a numerical argument and returns the Eigensystem, or is M already numerical but written in terms of arbitrary-precision numbers, or something else? Please include a minimum example to clarify what you're trying to do. – Jens Oct 18 '14 at 21:06
  • You are already given M. Let's say M=RandomReal[{0,1},{1000,1000}]. Since you know in this case this matrix is made of real numbers (numeric, no symbols) how do you pass this information to eigensystem to speed up the computation? – lagoa Oct 18 '14 at 21:09
  • For example, myeig = Compile[{{mat, Complex, 2}}, Eigenvalues[mat], {{Eigenvalues[], _Complex, 1}}, RuntimeAttributes -> {Listable}, CompilationTarget -> "C"]; will speed up the computation of eigenvalues of a matrix of complex numbers. But I couldn't generalize this to use Eigensystem. – lagoa Oct 18 '14 at 21:21
  • 1
    @Iagoa: What you are looking for does not exist, to my knowledge. In particular, Eigenvalues is not a compilable function, since it is outsourced to a backend LAPACK/BLAS solver that operates independently of Mathematica. Likewise, you can't tell it to use CompilationTarget -> "C" because the computations are not being done in C, but are done by LAPACK. For more information, see tutorial/SomeNotesOnInternalImplementation and the section "Approximate Numerical Linear Algebra". – DumpsterDoofus Oct 18 '14 at 21:29
  • Beyond that, I personally don't know how the internals work. I have noticed that Mathematica seems to auto-detect Hermiticity, symmetry, and other matrix properties when Eigenvalues is called: for example, calling Eigenvalues on an approximate complex Hermitian matrix returns a list of purely-real floating point numbers, with no machine-precision-size complex parts present. Other than that, I can't say what sort of auto-detection is or isn't occurring. Perhaps someone else with more knowledge can chime in with more info. – DumpsterDoofus Oct 18 '14 at 21:31
  • Related: Is there any faster way than Eigensystem to diagonalize a Hermitian matrix? – ybeltukov Oct 18 '14 at 22:17
  • Your gravatar is a near-twin of mine. I am surprised I never noticed this before since you have been participating in Mathematica.SE for well over a year. However, I only noticed now when I saw my gravatar next to yours after editing this question. – m_goldberg Oct 18 '14 at 22:36

1 Answers1

2

Eigensystem is already informed

A = RandomReal[1, {1000, 1000}];
A += Transpose[A];

Eigensystem[A]; // AbsoluteTiming
(* {1.034878, Null} *)

Eigensystem[A + 0. I]; // AbsoluteTiming
(* {2.645509, Null} *)

Update: type detection timings:

Needs["GeneralUtilities`"];

r = RandomReal[1, {100, 100}];
c = RandomComplex[1 + I, {100, 100}];

MatrixQ[r, Internal`RealValuedNumericQ] // AccurateTiming
MatrixQ[c, Internal`RealValuedNumericQ] // AccurateTiming
(* 8.60352*10^-7 *)
(* 8.61328*10^-7 *)

r = RandomReal[1, {1000, 1000}];
c = RandomComplex[1 + I, {1000, 1000}];

MatrixQ[r, Internal`RealValuedNumericQ] // AccurateTiming
MatrixQ[c, Internal`RealValuedNumericQ] // AccurateTiming
(* 8.4375*10^-7 *)
(* 8.92578*10^-7 *)

0.9 microsecond is fast enough to forget about it. Moreover, the timing doesn't depend on the matrix size because matrices are packed with the proper information about the type

r // Developer`PackedArrayForm
c // Developer`PackedArrayForm

(* PackedArray[Real, <1000, 1000>] *)
(* PackedArray[Complex, <1000, 1000>] *)
ybeltukov
  • 43,673
  • 5
  • 108
  • 212