The class MultiLeastSQquare solve the least square problem when the response b is multidimensional. More...

#include <STK_lapack_MultiLeastSquare.h>

Inheritance diagram for STK::lapack::MultiLeastSquare< ArrayB, ArrayA >:

Public Types
typedef ILeastSquare< MultiLeastSquare< ArrayB, ArrayA > >	Base

Public Member Functions
	MultiLeastSquare (ArrayB const &b, ArrayA const &a, bool isBref=false, bool isAref=false)
	constructor

template<class ArrayB_ , class ArrayA_ >
	MultiLeastSquare (ExprBase< ArrayB_ > const &b, ExprBase< ArrayA_ > const &a)
	template constructor

virtual	~MultiLeastSquare ()
	Destructor.

Real	rcond () const

CVectorX const &	s () const
	return the array with the singular values of A

void	setRcond (Real rcond)

bool	runImpl ()
	solve the multi-linear least square problem.

template<class Weights >
bool	runImpl (Weights const &weights)
	solve the weighted least square problem.

bool	runImpl ()

Public Member Functions inherited from STK::ILeastSquare< MultiLeastSquare< ArrayB, ArrayA > >
virtual	~ILeastSquare ()
	Destructor.

Integer const &	rank () const

ArrayA const &	a () const

ArrayB const &	b () const

ArrayB const &	x () const

virtual bool	run ()
	Compute the Least-Square solution.

bool	run (VecWeights const &weights)
	Compute the weighted Least-Square solution.

bool	run (Weights const &weights)

Public Member Functions inherited from STK::IRunnerBase
String const &	error () const
	get the last error message.

Public Member Functions inherited from STK::IRecursiveTemplate< Derived >
Derived &	asDerived ()
	static cast : return a reference of this with a cast to the derived class.

Derived const &	asDerived () const
	static cast : return a const reference of this with a cast to the derived class.

Derived *	asPtrDerived ()
	static cast : return a ptr on a `Derived` of this with a cast to the derived class.

Derived const *	asPtrDerived () const
	static cast : return a ptr on a constant `Derived` of this with a cast to the derived class.

Derived *	clone () const
	create a leaf using the copy constructor of the Derived class.

Derived *	clone (bool isRef) const
	create a leaf using the copy constructor of the Derived class and a flag determining if the clone is a reference or not.

Protected Attributes
Real	rcond_
	condition number used for determining the effective rank of A

CVectorX	s_
	Array of the singular values.

Protected Attributes inherited from STK::ILeastSquare< MultiLeastSquare< ArrayB, ArrayA > >
ArrayB	b_
	Array or vector of the left hand side.

ArrayA	a_
	Array of the right hand side.

ArrayB	x_
	Array of the solution (a vector if b is a vector, a matrix otherwise)

Integer	rank_
	rank of matrix A

Protected Attributes inherited from STK::IRunnerBase
String	msg_error_
	String with the last error message.

bool	hasRun_
	`true` if run has been used, `false` otherwise

Private Member Functions
bool	computeLS (CArrayXX &b, CArrayXX &a)
	private method for computing the LS solution

Additional Inherited Members
Protected Types inherited from STK::ILeastSquare< MultiLeastSquare< ArrayB, ArrayA > >
typedef hidden::AlgebraTraits< MultiLeastSquare< ArrayB, ArrayA > >::ArrayB	ArrayB

typedef hidden::AlgebraTraits< MultiLeastSquare< ArrayB, ArrayA > >::ArrayA	ArrayA

Protected Member Functions inherited from STK::ILeastSquare< MultiLeastSquare< ArrayB, ArrayA > >
	ILeastSquare (ArrayB const &b, ArrayA const &a, bool isBref=false, bool isAref=false)
	Default constructor.

	ILeastSquare (ExprBase< OtherArrayB > const &b, ExprBase< OtherArrayA > const &a)
	template constructor

Protected Member Functions inherited from STK::IRunnerBase
	IRunnerBase ()
	default constructor

	IRunnerBase (IRunnerBase const &runner)
	copy constructor

virtual	~IRunnerBase ()
	destructor

virtual void	update ()
	update the runner.

Protected Member Functions inherited from STK::IRecursiveTemplate< Derived >
	IRecursiveTemplate ()
	constructor.

	~IRecursiveTemplate ()
	destructor.

Detailed Description

template<class ArrayB, class ArrayA>
class STK::lapack::MultiLeastSquare< ArrayB, ArrayA >

The class MultiLeastSQquare solve the least square problem when the response b is multidimensional.

The class MultiLeastSquare allows to solve the least-square problem

$\min_{x} \|b - A*x\|^2.$

A is an m-by-n matrix which may be rank-deficient.
B can be a vector or a Matrix.

It computes the minimum-norm solution to a real linear least squares problem: minimize 2-norm(| b - A*x |) using the singular value decomposition (SVD) of A. A is an M-by-N matrix which may be rank-deficient.

See also: STK::ILeastSquare, STK::MultiLeastSquare, STK::lapack::MultiLeastSquare

Definition at line 101 of file STK_lapack_MultiLeastSquare.h.

Member Typedef Documentation

◆ Base

template<class ArrayB , class ArrayA >

typedef ILeastSquare< MultiLeastSquare<ArrayB, ArrayA> > STK::lapack::MultiLeastSquare< ArrayB, ArrayA >::Base

Definition at line 104 of file STK_lapack_MultiLeastSquare.h.

Constructor & Destructor Documentation

◆ MultiLeastSquare() [1/2]

template<class ArrayB , class ArrayA >

STK::lapack::MultiLeastSquare< ArrayB, ArrayA >::MultiLeastSquare	(	ArrayB const &	b,
		ArrayA const &	a,
		bool	isBref = `false`,
		bool	isAref = `false`
	)

inline

constructor

Parameters

b,a	the left hand side and the right hand side of the least square problem.
isBref,isAref	are the left hand side and the right hand side references ?

Definition at line 113 of file STK_lapack_MultiLeastSquare.h.

114 : Base(b, a, isBref, isAref), rcond_(-1) {};

STK::ILeastSquare< MultiLeastSquare< ArrayB, ArrayA > >::b

ArrayB const & b() const

Definition STK_ILeastSquare.h:101

STK::ILeastSquare< MultiLeastSquare< ArrayB, ArrayA > >::a

ArrayA const & a() const

Definition STK_ILeastSquare.h:99

STK::lapack::MultiLeastSquare::Base

ILeastSquare< MultiLeastSquare< ArrayB, ArrayA > > Base

Definition STK_lapack_MultiLeastSquare.h:104

STK::lapack::MultiLeastSquare::rcond_

Real rcond_

condition number used for determining the effective rank of A

Definition STK_lapack_MultiLeastSquare.h:145

◆ MultiLeastSquare() [2/2]

template<class ArrayB , class ArrayA >

template<class ArrayB_ , class ArrayA_ >

STK::lapack::MultiLeastSquare< ArrayB, ArrayA >::MultiLeastSquare	(	ExprBase< ArrayB_ > const &	b,
		ExprBase< ArrayA_ > const &	a
	)

inline

template constructor

Parameters

b,a	the left hand side and the right hand side of the least square problem.

Definition at line 120 of file STK_lapack_MultiLeastSquare.h.

121 : Base(b, a), rcond_(-1) {}

◆ ~MultiLeastSquare()

template<class ArrayB , class ArrayA >

virtual STK::lapack::MultiLeastSquare< ArrayB, ArrayA >::~MultiLeastSquare ( )

inlinevirtual

Destructor.

Definition at line 123 of file STK_lapack_MultiLeastSquare.h.

123{};

Member Function Documentation

◆ computeLS()

template<class ArrayB , class ArrayA >

bool STK::lapack::MultiLeastSquare< ArrayB, ArrayA >::computeLS	(	CArrayXX &	b,
		CArrayXX &	a
	)

private

private method for computing the LS solution

Definition at line 196 of file STK_lapack_MultiLeastSquare.h.

{
  Range arows = a.rows(), acols = a.cols(), brows = b.rows(), bcols = b.cols();
  int m = arows.size(), n= acols.size(), nrhs = bcols.size();
  // resize if necessary
  if (brows.size()<n)
  {
    if (!b.isRef())
    {
      CArrayXX tmp(b);
      b.resize(acols,bcols);
      b.sub(brows, bcols) = tmp;
    }
    else
    { STKRUNTIME_ERROR_NO_ARG(MultiLeastSquare::computeLS,b has not enough rows);}
  }
  int lda = a.sizeRows(), ldb = b.sizeRows();
  // shift to 0
  a.shift(0,0);
  b.shift(0,0);
  s_.resize(m); s_ = 0.;
  Range srange = s_.range();
  s_.shift(0);
  // get sizes
  Real wkopt;
  int lWork = -1, iwkopt;
  int info = gelsd( m, n, nrhs
                  , a.p_data(), lda
                  , b.p_data(), ldb
                  , s_.p_data()
                  , &rcond_, &rank_
                  , &wkopt, lWork, &iwkopt);
  if( info != 0 )
  {
    if (info>0)
    { this->msg_error_ = STKERROR_NO_ARG(lapack::MultiLeastSquare::computeLS,internal error);
      return false;
    }
    this->msg_error_= STKERROR_1ARG(lapack::MultiLeastSquare::computeLS get,-info,error parameter);
    return false;
  }
  // get working sizes
  lWork = (int)wkopt;
  Real* work = new Real[lWork];
  int* iwork = new int[iwkopt];
  // Solve the least square problem
  info = gelsd( m, n, nrhs
              , a.p_data(), lda
              , b.p_data(), ldb
              , s_.p_data()
              , &rcond_, &rank_
              , work, lWork, iwork);
  delete[] iwork;
  delete[] work;
  if( info != 0 )
  {
    if (info>0)
    { this->msg_error_ = STKERROR_NO_ARG(lapack::MultiLeastSquare::computeLS,internal error);
      return false;
    }
    this->msg_error_= STKERROR_1ARG(lapack::MultiLeastSquare::computeLS get,-info,error parameter);
    return false;
  }
  // shift back
  a.shift(arows.begin(), acols.begin());
  b.shift(brows.begin(), bcols.begin());
  s_.shift(srange.begin());
  x_ = b.sub(acols,bcols);
  return true;
}

References STK::lapack::MultiLeastSquare< ArrayB, ArrayA >::computeLS(), STK::lapack::gelsd(), STK::ICArray< Derived >::isRef(), STK::ICArray< Derived >::p_data(), STK::ICArray< Derived >::resize(), STK::ICArray< Derived >::shift(), STKERROR_1ARG, STKERROR_NO_ARG, STKRUNTIME_ERROR_NO_ARG, and STK::ICArray< Derived >::sub().

Referenced by STK::lapack::MultiLeastSquare< ArrayB, ArrayA >::computeLS().

◆ rcond()

template<class ArrayB , class ArrayA >

Real STK::lapack::MultiLeastSquare< ArrayB, ArrayA >::rcond ( ) const

inline

Returns: the condition number

Definition at line 125 of file STK_lapack_MultiLeastSquare.h.

125{ return rcond_;}

References STK::lapack::MultiLeastSquare< ArrayB, ArrayA >::rcond_.

Referenced by STK::lapack::MultiLeastSquare< ArrayB, ArrayA >::setRcond().

◆ runImpl() [1/3]

template<class ArrayB , class ArrayA >

bool STK::lapack::MultiLeastSquare< ArrayB, ArrayA >::runImpl ( )

solve the multi-linear least square problem.

Launch gelsd LAPACK routine to perform the decomposition.

Returns: true if no error occur, false otherwise

Definition at line 159 of file STK_lapack_MultiLeastSquare.h.

{
  Range brows = (b_.sizeRows()< a_.sizeCols()) ? a_.cols() : b_.rows();
  // local arrays, b is resized if necessary
  CArrayXX a(a_), b(brows, b_.cols());
  b.sub(b_.rows(), b_.cols()) = b_;
  // start
  return computeLS(b,a);
}

References STK::ICArray< Derived >::sub().

◆ runImpl() [2/3]

bool STK::lapack::MultiLeastSquare< CArrayXX, CArrayXX >::runImpl ( )

inline

Definition at line 190 of file STK_lapack_MultiLeastSquare.h.

191{ return computeLS(b_,a_);}

◆ runImpl() [3/3]

template<class ArrayB , class ArrayA >

template<class Weights >

bool STK::lapack::MultiLeastSquare< ArrayB, ArrayA >::runImpl ( Weights const & weights )

solve the weighted least square problem.

Launch gelsd LAPACK routine to perform the decomposition.

Returns: true if no error occur, false otherwise

Definition at line 174 of file STK_lapack_MultiLeastSquare.h.

{
  STK_STATIC_ASSERT_ONE_DIMENSION_ONLY(Weights);
  Range brows = (b_.sizeRows()< a_.sizeCols()) ? a_.cols() : b_.rows();
  // local arrays, b is resized if necessary
  CArrayXX a(w.sqrt().diagonalize() * a_), b(brows, b_.cols());
  b.sub(b_.rows(), b_.cols()) = w.sqrt().diagonalize() * b_;
  // start
  return computeLS(b,a);
}