Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Compute trace and determinant #9

Open
andreaferretti opened this issue Jun 13, 2017 · 4 comments
Open

Compute trace and determinant #9

andreaferretti opened this issue Jun 13, 2017 · 4 comments

Comments

@andreaferretti
Copy link
Owner

The computation of eigenvalues can be used for the determinant, or perhaps there is a more specific LAPACK function

@ghost
Copy link

ghost commented Jun 15, 2017

The Faq mentions:

Are there routines in LAPACK to compute determinants?

No. There are no routines in LAPACK to compute determinants. This is discussed in the "Accuracy and Stability" chapter in the LAPACK Users' Guide.

In gonum this is done using the LU factorization of the matrix.

var lu LU
lu.Factorize(a)
lu.LogDet()
det, sign := lu.LogDet()
res := math.Exp(det) * sign

@andreaferretti
Copy link
Owner Author

Yeah, I know it's not very stable. There is no LU factorization here yet, but the Schur form can be used for now (although it is more costly; the determinant will be unstable for large matrices anyway)

@andreaferretti
Copy link
Owner Author

742c60f introduced the two functions. det can be optimized, but as remarked it only makes sense for small matrices anyway, so it not much of an issue for now

@ryandvmartin
Copy link

I'm too new to nim and this to make a PR, but I cobbled together a cholesky decomp. I'm not sure if I used the T right since I think that has to be float64 with the dpotrf call..
anyway here it is:

proc cholesky*[T](a: Matrix[T]): Matrix[T] =
  var
    h = a.clone()
    n = a.ld.cint
    info: cint
    ulo: cstring
  ulo = "L"
  dpotrf(ulo, addr n, h.fp, addr n, addr info)
  if info > 0:
    raise newException(LinearAlgebraError, "ERROR finding the LU decomp")
  for i in 0 ..< n:
    for j in i + 1 ..< n:
      h[i, j] = 0.0
  return h

awesome work on this package!

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

No branches or pull requests

2 participants