tnt/tnt_lapack.h

#ifndef TNT_LAPACK_HPP
#define TNT_LAPACK_HPP

extern "C"
{
#include "f2c.h"
#include "clapack.h"
}


namespace TNT
{

  Fortran_Array2D<double> Lapack_LinearSolve( const Fortran_Array2D<double> &A,  
      const Fortran_Array2D<double> &B)
  {
      
      integer m = A.dim1();
      integer n = A.dim2();

      if (m != n)
        return Fortran_Array2D<double>();

      integer nrhs = B.dim2();

      if (n != B.dim1())
        return Fortran_Array2D<double>();


      Fortran_Array2D<double> X(B.copy());
      Fortran_Array2D<double> lu(A.copy());
      Fortran_Array1D<integer>  ipiv(m, integer(0));

      integer  info=0;

      dgesv_( &n, &nrhs, &lu(1,1), &m, &ipiv(1), &X(1,1), &n, &info);

      return X;
  }



}


#endif
// TNT_LAPACK_HPP

---