tmatrix.h

将C语言的常用程序集移植到VC开发环境的源代码

/*********************************************************************
 Matrix.h: interface&implementation for the TMatrix class.
 Write by hujinshan@2003.6.26 Airforce Engineering University
*********************************************************************/

#if !defined(HJS_TMATRIX_H)
#define HJS_TMATRIX_H

#include <valarray>
using namespace std;

template<class T>
class Slice_iter{
	valarray<T>* v;
	slice s;
	size_t curr;//index of current element

	T& ref(size_t i)const{return (*v)[s.start()+i*s.stride()];}
public:
	Slice_iter(valarray<T>* vv,slice ss):v(vv),s(ss),curr(0){}
	Slice_iter end()
	{
		Slice_iter t=*this;
		t.curr=s.size();//index of last-plus-one element
		return t;
	}
	Slice_iter& operator++(){curr++;return *this;}
	Slice_iter operator++(int){Slice_iter t=*this;curr++;return t;}
	size_t size(){return v->size();}
	T& operator[](size_t i){return ref(curr=i);}//C-style subscript
	T& operator()(size_t i){return ref(curr=i);}//Fortran-style subscript
	T& operator*(){return ref(curr);}//current element
	//...
};

template<class T>
class Cslice_iter{
	valarray<T>* v;
	slice s;
	size_t curr;//index of current element

	T& ref(size_t i)const{return (*v)[s.start()+i*s.stride()];}
public:
	Cslice_iter(valarray<T>* vv,slice ss):v(vv),s(ss),curr(0){}
	Cslice_iter end()
	{
		Cslice_iter t=*this;
		t.curr=s.size();//index of last-plus-one element
		return t;
	}
	Cslice_iter& operator++(){curr++;return *this;}
	Cslice_iter operator++(int){Cslice_iter t=*this;curr++;return t;}
	size_t size(){return v->size();}
	T& operator[](size_t i){return ref(curr=i);}//C-style subscript
	T& operator()(size_t i){return ref(curr=i);}//Fortran-style subscript
	T& operator*(){return ref(curr);}//current element
	//...
};

template<class T=double>
class TMatrix
{
	valarray<T>* v;
	size_t d1, d2;
public:
	TMatrix(size_t x,size_t y);
	TMatrix(const TMatrix&);
	TMatrix<T>& operator=(const TMatrix&);

	size_t size()const{return d1*d2;}
	size_t dim1()const{return d1;}
	size_t dim2()const{return d2;}

	Slice_iter<T> row(size_t i);
	Cslice_iter<T> row(size_t i)const;
	Slice_iter<T> column(size_t i);	
	Cslice_iter<T> column(size_t i)const;

	T& operator()(size_t x,size_t y);//Fortran-style subscripts
	T operator()(size_t x,size_t y)const;

	Slice_iter<T> operator()(size_t i){return row(i);}//C-style subscripts
	Cslice_iter<T> operator()(size_t i)const{return row(i);}

	Slice_iter<T> operator[](size_t i){return row(i);}//C-style subscripts
	Cslice_iter<T> operator[](size_t i)const{return row(i);}

	TMatrix<T>& operator*=(T);
	valarray<T>& array(){return *v;}
	virtual ~TMatrix();

	void ZeroMatrix();
};

template <class T>
TMatrix<T>::TMatrix(size_t x,size_t y)
{
	d1=x;d2=y;
	v=new valarray<T>(d1*d2);
}
template <class T>
TMatrix<T>::TMatrix(const TMatrix& mat)
{
	(*this)=mat;
}

template <class T>
TMatrix<T>::~TMatrix()
{
	delete v;
}

template <class T>
TMatrix<T>& TMatrix<T>::operator=(const TMatrix& mat)
{
	d1=mat.dim1();
	d2=mat.dim2();
	v=new valarray<T>(d1*d2);
	for(int i=0;i<d1;i++)
		for(int j=0;j<d2;j++)
			(*this) [i][j]=mat[i][j];
	return *this;
}

template <class T>
T& TMatrix<T>::operator() (size_t x,size_t y)
{
	return row(x)[y];
}

template <class T> 
static T mul(Slice_iter<T>& v1,const valarray<T>& v2)
{
	T res=0;
	for(int i=0;i<v1.size();i++)
		res+=v1[i]*v2[i];
	return res;
}
template <class T> 
static valarray<T> operator*(TMatrix<T>& m,const valarray<T>& v)
{
	valarray<T> res(m.dim1());
	for(int i=0;i<m.dim1();i++)
		res[i]=mul(m.row(i),v);
	return res;
}

template <class T>
TMatrix<T>& TMatrix<T>::operator*=(T d)
{
	(*v) *=d;
	return *this;
}

template <class T>
inline Slice_iter<T> TMatrix<T>::row(size_t i)
{
	return Slice_iter<T>(v,slice(i,d2,d1));
}

template <class T>
inline Cslice_iter<T> TMatrix<T>::row(size_t i)const
{
	return Cslice_iter<T>(v,slice(i,d2,d1));
}

template <class T>
inline Slice_iter<T> TMatrix<T>::column(size_t i)
{
	return Slice_iter<T>(v,slice(i*d2,d2,1));
}

template <class T>
inline Cslice_iter<T> TMatrix<T>::column(size_t i)const
{
	return Cslice_iter<T>(v,slice(i*d2,d2,1));
}

template <class T>
inline void TMatrix<T>::ZeroMatrix()
{
	for(int i=0;i<size();i++)
		(*v)[i]=0;
}
#endif // !defined(HJS_TMATRIX_H)
//end of file