eigendecomposition is not usable?

[idavydov]

6015d269f11686cecdb1622a6ac1060ca8201862 makes it impossible to use matrix eigendecomposition outside the library. Change from Eigen to eigen makes the function non-exported.

[btracey]

Could you explain what you mean?

var eigen mat64.Eigen
ok := eigen.Factorize(a)
if !ok {
   log.Fatal("bad decomp")
}
evs := eigen.Values(nil, false)
fmt.Println("eigenvalues = ", evs)

[idavydov]

@btracey is it possible to get eigenvectors as well?

[btracey]

Ah, no, not right at the moment. Have you been using that extraction? The plan is to put it in, but I don't know the resolution to https://github.com/gonum/matrix/issues/308 . If you need eigenvectors, we can temporarily make EigenFactors still public until that issue is resolved.

[idavydov]

@btracey yes I've been using eigenvalues and eigenvectors in my code: https://bitbucket.org/Davydov/godon/src/c45db765c6f7b364454665151241cadc7572a259/cmodel/ematrix.go?at=master&fileviewer=file-view-default#ematrix.go-46 Having a way to access them will be very nice, thanks.

[btracey]

PR submitted. Couple of side notes though

1) You are aware we have an Exp function in mat64 already? https://godoc.org/github.com/gonum/matrix/mat64#Dense.Exp 2) Do you know your matrix is symmetric? Are your sure your code works even if complex eigenvalues are generated? 3) It's almost always a bad idea to compute an explicit matrix inverse. The code

res := mat64.NewDense(cols, rows, nil)
res.Mul(m.v, cD)
res.Mul(res, m.iv)

is better as

res := mat64.NewDense(cols, rows, nil)
err := res.Solve(m.v.T(), cD.T())
if err != nil {
    log.Fatal("bad solve")
}
res.Mul(m.v, res)

[idavydov]

@btracey thanks a lot. In my problem I have to compute exp(t*Q) multiple times, where Q is a constant matrix and t is a variable scalar. So performing a single eigendecomposition followed by a number of dot products is much faster. The matrix is non-symmetric. But since my Q matrix corresponds to the time-reversible continuous-time Markov proccess, there should be no complex eigenvalues. And the exponential should be strictly positive. Thanks a lot.

[kortschak]

This use case is actually what I wanted the real power Exp function for. Nice to see gonum being used in this space.

[btracey]

We could add a PowEig function that takes in eigenvalues []float64 and eigenvectors

[kortschak]

SGTM

[idavydov]

@btracey having PowEig function like this would be very helpful indeed.

[idavydov]

Not sure if that's a separate bug report. But I don't think Factorize works. Sample code When I run it, the program hangs. This matrix was successfully factorized previously.

[kortschak]

Is this an issue with the epsilon being used? Would you try replacing the zero in the eigen call with whatever epsilon you were using in the Eigen call previously (if not zero).

[idavydov]

@kortschak epsilon function isn't exported anymore (since commit 6015d269f11686cecdb1622a6ac1060ca8201862). And there's no way to specify epsilon for the Factorize function. Am I missing something?

[kortschak]

Make a change to the eigen.go file and see if it fixes the problem.

[idavydov]

@kortschak yes, replacing 0 in the eigen call with 1e-9 or 1e-6 solves the problem.

[kortschak]

OK. Let's go with that as an interim fix.

[idavydov]

HI, it seems that in some cases eigen still hangs. Here's an example matrix: https://gist.github.com/idavydov/eb4c88b4663a7240c92f It works ok with a smaller epsilon value (1e-16 hangs, 1e-15 works). I propose 1e-14 to be safe, as far as I understand this code will be replaced anyway in future.