[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
Re: Octave/C++ matrix Inv() comparison
From: |
c. |
Subject: |
Re: Octave/C++ matrix Inv() comparison |
Date: |
Mon, 8 Jul 2013 19:44:29 +0200 |
On 8 Jul 2013, at 18:35, Ed Meyer <address@hidden> wrote:
>
> >> the numbers agree in the first two places so there may not be programming
> >> errors
> >> but the lapack functions used by octave do equilibration and pivoting to
> >> improve
> >> the solution - why not just use lapack (dgesvx)?
>
>
>
> > thank you for your answer
>
> > it "deblocks" me a little ^^. So what you said is that the C++ routine i
> > use are not as precise > as in Octave ? it makes sense now that you say
> > it...
>
> > So I looked at the subroutine dgesvx as you told me
>
> > I see that there are a lot of parameters
>
> > Could you help me with those parameters ? I don't understand a clue
>
> > All i need to do i take the inverse of a <double** matrix> with <int n> the
> > size of the matrix
>
> > but i don't understand what are all those parameters, even with the
> > explanation of the routine
>
> > mucho thanks for the help )))
>
> > Jeff
>
> If you haven't called fortran routines from C++ before there are a few issues
> to
> be aware of; see e.g.
> http://www.yolinux.com/TUTORIALS/LinuxTutorialMixingFortranAndC.html
> the main things you have to be aware of is how fortran stores matrices, all
> args
> are pointers, and the naming convention (e.g. on linux dgesv is dgesv_).
>
> It would probably be worthwhile trying the simpler routine dgesv which
> doesn't do
> equilibration but has a simpler interface:
>
> dgesv (n, nrhs, a, lda, ipiv, b, ldb, info)
>
> where "b" is an (n,n) identity matrix (I like using stl vectors for temp
> arrays):
>
> vector<double> b(n*n, 0.0);
> for (i=0; i<n; i++)
> b[i*(n+1)] = 1.0;
>
> ipiv is an output array of ints:
>
> vector<int> ipiv(n);
> and you call with pointers for ints:
>
> dgesv_ (&n, &n, &a[0], &n, &ipiv[0], &b[0], &ldb, &info);
>
> Calling fortran from C++ is a pain but it does open up a huge amount
> of code once you figure out how to do it.
>
> --
> Ed Meyer
Otherwise,
If you want a linear algebra library with a more typical C++ interface you may
want to have a look at eigen [1]
c.
[1] http://eigen.tuxfamily.org