#ifndef  IntfLapack_H_SEEN
#define  IntfLapack_H_SEEN

#include "machdefs.h"
#include "tarray.h"

namespace SOPHYA {

template <class T>
class LapackServer { 
public:
  LapackServer();
  virtual ~LapackServer();

  virtual int LinSolve(TArray<T>& a, TArray<T> & b); 
  virtual int SVD(TArray<T>& a, TArray<T> & s); 
  virtual int SVD(TArray<T>& a, TArray<T> & s, TArray<T> & u, TArray<T> & vt); 

  //! Set the workspace size factor for LAPACK routines
  inline void SetWorkSpaceSizeFactor(int f = 2)
  { wspace_size_factor = (f > 1) ? f : 1; }

  //! Returns the workspace size factor
  inline int  GetWorkSpaceSizeFactor() 
  { return wspace_size_factor; }

private:
  int SVDDriver(TArray<T>& a, TArray<T> & s, 
		TArray<T>* up=NULL, TArray<T> * vtp=NULL);

  int wspace_size_factor;
};

/*! \ingroup LinAlg
    \fn LapackLinSolve(TArray<T>&, TArray<T> &)
    \brief Solves the linear system A*X = B using LapackServer. 
*/
template <class T>
inline int LapackLinSolve(TArray<T>& a, TArray<T> & b)
{ LapackServer<T> lps; return( lps.LinSolve(a, b) );  }

/*! \ingroup LinAlg
    \fn LapackSVD(TArray<T>&, TArray<T> &)
    \brief SVD decomposition using LapackServer. 
*/
template <class T>
inline int LapackSVD(TArray<T>& a, TArray<T> & s)
{ LapackServer<T> lps; return( lps.SVD(a, s) ); }


/*! \ingroup LinAlg
    \fn LapackSVD(TArray<T>&, TArray<T> &, TArray<T> &, TArray<T> &)
    \brief SVD decomposition using LapackServer. 
*/
template <class T>
inline int LapackSVD(TArray<T>& a, TArray<T> & s, TArray<T> & u, TArray<T> & vt)
{ LapackServer<T> lps; return( lps.SVD(a, s, u, vt) ); }


} // Fin du namespace

#endif