function B=partialtrace(m,n,P)
%function B=partialtrace(m,n,P)
%  evaluate the partial trace to check for trace preserving constraint
%INPUT m,n where m*n is order of P
%OUTPUT:  B   partial trace
B=zeros(n);
for tt=1:n,
	for ss=1:tt,
		temp=trace(P((ss-1)*m+1:ss*m,(tt-1)*m+1:tt*m));
		B(ss,tt)=temp;
		if ss<tt,
		        B(tt,ss)=temp';
		end
	end
end
