Commit 3e269c60 authored by Alvarez, Gonzalo's avatar Alvarez, Gonzalo
Browse files

closure operations for matrix

parent 02e6ac9c

drivers/closuresTest.cpp

+12 −1
Original line number Diff line number Diff line
@@ -12,7 +12,7 @@ void testVector() 
	v3 <= v2*1.2;

	std::cout<<v3;



//	// plus

	// plus

	v3 <= v1 + v2;

	std::cout<<v3;

	v3 <= v1 + 1.3*v2;
@@ -44,6 +44,17 @@ void testMatrix() 
	// multiplication by scalar

	m2 = m1*1.2;

	std::cout<<m2;



	// plus

	PsimagLite::Matrix<double> m3;

	m3 = m1 + m2;

	std::cout<<m3;

	m3 = m1 + 1.3*m2;

	std::cout<<m3;

	m3 = 1.3*m2 + m1;

	std::cout<<m3;



	// minus

}



int main(int, char**)

src/Matrix.h

+69 −28
Original line number Diff line number Diff line
@@ -218,6 +218,55 @@ public: 
		MPI::bcast(data_,root,mpiComm);

	}



	// start closure memebers



	template<typename T1>

	void operator=(const std::ClosureOperator<T1,Matrix<T>,std::ClosureOperations::OP_MULT>& c)

	{

		nrow_ = c.v2_.nrow_;

		ncol_ = c.v2_.ncol_;

		this->data_ <= c.v1_*c.v2_.data_;

	}



	void operator=

	(const std::ClosureOperator<Matrix<T>,Matrix<T>,std::ClosureOperations::OP_PLUS>& c)

	{

		nrow_ = c.v1_.nrow_;

		ncol_ = c.v1_.ncol_;

		assert(nrow_ == c.v2_.nrow_);

		assert(ncol_ == c.v2_.ncol_);

		this->data_ <= c.v1_.data_ + c.v2_.data_;

	}



	template<typename T1>

	void operator=

	(const std::ClosureOperator<Matrix<T>,

	 std::ClosureOperator<T1,Matrix<T>,std::ClosureOperations::OP_MULT>,

	 std::ClosureOperations::OP_PLUS>& c)

	{

		nrow_ = c.v1_.nrow_;

		ncol_ = c.v1_.ncol_;

		assert(nrow_ == c.v2_.v2_.nrow_);

		assert(ncol_ == c.v2_.v2_.ncol_);

		this->data_ <= c.v1_.data_ + c.v2_.v1_*c.v2_.v2_.data_;

	}



	template<typename T1>

	void operator=

	(const std::ClosureOperator<

	 std::ClosureOperator<T1,Matrix<T>,std::ClosureOperations::OP_MULT>,

	 Matrix<T>,

	 std::ClosureOperations::OP_PLUS>& c)

	{

		nrow_ = c.v2_.nrow_;

		ncol_ = c.v2_.ncol_;

		assert(nrow_ == c.v1_.v2_.nrow_);

		assert(ncol_ == c.v1_.v2_.ncol_);

		this->data_ <= c.v2_.data_ + c.v1_.v1_*c.v1_.v2_.data_;

	}



	// end closure members



private:



	SizeType nrow_,ncol_;
@@ -444,24 +493,32 @@ bool isZero(const Matrix<T>& m) 


// closures start



template<typename T1,typename T2>

std::ClosureOperator<T1,Matrix<T2>,std::ClosureOperations::OP_MULT> operator*(const T1& val,

                                                                              const Matrix<T2>& a)

{

	return std::ClosureOperator<T1,Matrix<T2>,std::ClosureOperations::OP_MULT>(val,a);

}



template<typename T1,typename T2>

std::ClosureOperator<T1,Matrix<T2>,std::ClosureOperations::OP_MULT> operator*(const Matrix<T2>& a,

                                                                              const T1& val)

{

	return val*a;

}



template<typename T>

Matrix<T> operator+(const Matrix<T>& a,const Matrix<T>& b)

std::ClosureOperator<Matrix<T>,Matrix<T>,std::ClosureOperations::OP_PLUS>

operator+(const Matrix<T>& a,const Matrix<T>& b)

{

	Matrix<T> c(a.n_row(),a.n_col());

	for (SizeType i=0;i<a.n_row();i++)

		for (SizeType j=0;j<a.n_col();j++)

			c(i,j) = a(i,j) + b(i,j);

	return c;

	return std::ClosureOperator<Matrix<T>,Matrix<T>,std::ClosureOperations::OP_PLUS>(a,b);

}



template<typename T>

Matrix<T> operator-(const Matrix<T>& a,const Matrix<T>& b)

std::ClosureOperator<Matrix<T>,Matrix<T>,std::ClosureOperations::OP_MINUS>

operator-(const Matrix<T>& a,const Matrix<T>& b)

{

	Matrix<T> c(a.n_row(),a.n_col());

	for (SizeType i=0;i<a.n_row();i++)

		for (SizeType j=0;j<a.n_col();j++)

			c(i,j) = a(i,j) - b(i,j);

	return c;

	return std::ClosureOperator<Matrix<T>,Matrix<T>,std::ClosureOperations::OP_MINUS>(a,b);

}



template<typename T>
@@ -509,22 +566,6 @@ VectorLikeType operator*(const VectorLikeType& b, 
	return v;

}



template<typename T1,typename T2>

Matrix<T2> operator*(const T1& val,const Matrix<T2>& a)

{

	Matrix<T2> b(a.n_row(),a.n_col());

	for (SizeType i=0;i<a.n_row();i++)

		for (SizeType j=0;j<b.n_col();j++)

			b(i,j) = val*a(i,j);

	return b;

}



template<typename T1,typename T2>

Matrix<T2> operator*(const Matrix<T2>& a,const T1& val)

{

	return val*a;

}



template<typename T>

T norm2(const Matrix<T>& m)

{