算法导论 第4章 源程序 分治策略 最大子数组 Strassen算法

chapter4.h

#ifndef CHAPTER4_H
#define CHAPTER4_H
#include<iostream>
#include<vector>
const int NEGTIVE_INFINITE=-1000000;
struct max_array 
{
	max_array(int left,int right,int sum):
		max_left(left),
		max_right(right),
		max_sum(sum){}
	void display(){
		std::cout<<"max_left="<<max_left<<" max_right="<<max_right
			<<" max_sum="<<max_sum<<std::endl;
	}
	int max_left;
	int max_right;
	int max_sum;
};
//find maximum crossing subarray(P40-41)
//called by max_array find_maximum_subarray (int*a ,int low,int high)
max_array find_max_crossing_subarray
	(int* a,int low,int mid,int high)
{
	int left_sum=NEGTIVE_INFINITE;
	int sum=0;
	int max_left=mid;
	for(int i=mid;i!=low-1;i--)
	{
		sum+=a[i];
		if(sum>left_sum)
		{
			left_sum=sum;
			max_left=i;
		}
	}
	int right_sum=NEGTIVE_INFINITE;
	sum=0;
	int max_right=mid+1;
	for(int i=mid+1;i<high+1;i++)
	{
		sum+=a[i];
		if(sum>right_sum)
		{
			right_sum=sum;
			max_right=i;
		}
	}
	return max_array(max_left,max_right,left_sum+right_sum);
}
// find maximum subarray recursively(P41)
max_array find_maximum_subarray
	(int*a ,int low,int high)
{
	if(low==high)
		return max_array(low,low,a[low]);
	else{
		int mid=(low+high)/2;
		max_array left_subarray=
			find_maximum_subarray(a,low,mid);
		max_array right_subarray=
			find_maximum_subarray(a,mid+1,high);
		max_array cross_subarray=
			find_max_crossing_subarray(a,low,mid,high);
		if(left_subarray.max_sum>=right_subarray.max_sum&&
			left_subarray.max_sum>=cross_subarray.max_sum)
			return left_subarray;
		else if(right_subarray.max_sum>=left_subarray.max_sum&&
			right_subarray.max_sum>=cross_subarray.max_sum)
			return right_subarray;
		else return cross_subarray;
	}
}
//4.1-2 a brute-force way to find maximum subarray
max_array find_max_subarray_brute
	(int* a,int low,int high)
{
	int left=low;
	int right=high;
	int max_sum=NEGTIVE_INFINITE;
	int sum=0;
	for(int i=low;i<=high;i++)
	{
		sum=0;
		for(int j=i;j<=high;j++)
		{
			sum+=a[j];
			if(max_sum<sum)
			{
				max_sum=sum;
				left=i;
				right=j;
			}
		}
	}
	return max_array(left,right,max_sum);
}
//4-1-5 find maximun subarray using linear time
max_array find_max_subarray_linear
	(int*a,int low,int high)
{
	int sum1=a[low];
	int sum2=a[low];
	int sum=a[low];
	int left=low;
	int right=low;
	int left_j=low;
	int right_j=low;
	int max_sum=NEGTIVE_INFINITE;
	int max_sum_j=NEGTIVE_INFINITE;
	for(int j=low;j<high;j++)
	{
		sum=0;
		sum2=NEGTIVE_INFINITE;
		for(int i=j+1;i>=low;i--)//not i<=low
		{
			sum+=a[i];
			if(sum>sum2)
			{
				sum2=sum;
				left_j=i;
			}
		}
		if(sum1<sum2)
		{
			sum1=sum2;
			right_j=j+1;
		}
		if(sum1>max_sum)
		{
			max_sum=sum1;
			left=left_j;
			right=right_j;
		}
	}
	return max_array(left,right,max_sum);
}
//Strassen matrix
class matrix{
private:
	int** p;
	int row;
	int col;
public:
	matrix(const int * a,const int m,const int n);
	int* create_pa(const int numbel);
	void display();
	matrix operator=(matrix& rhs);
	friend matrix operator+
		(const matrix& lhs,const matrix& matrix_rhs);
	friend matrix operator-
		(const matrix& lhs,const matrix& rhs);
	friend matrix operator*
		(const matrix& lhs, const matrix& rhs);
	~matrix();
	matrix split_quarter(int m,int n) const;
	friend matrix zip_quarter
		(const matrix& c11,const matrix& c12,
		 const matrix& c21,const matrix& c22);
	friend matrix strassen
		(const matrix& lhs,const matrix& rhs);
};
matrix::matrix(const int* a,const int m,const int n)
{
	row=m;
	col=n;
	p=new int*[row];
	for(int i=0;i<row;i++)
	{
		int* tmp=new int[col];
		for(int j=0;j<col;j++)
			tmp[j]=a[i*col+j];
		p[i]=tmp;
	}
}
matrix matrix::operator=(matrix& rhs)
{
	row=rhs.row;
	col=rhs.col;
	p=new int*[row];
	for(int i=0;i<row;i++)
	{
		int* tmp=new int[col];
		for(int j=0;j<col;j++)
			tmp[j]=rhs.p[i][j];
		p[i]=tmp;
	}
	return *this;
}
int* matrix::create_pa(const int numbel)
{
	int* p=new int[numbel];
	return p;
}
void matrix::display()
{
	for(int i=0;i<row;i++)
	{
		for(int j=0;j<row;j++)
			std::cout<<p[i][j]<<"\t";
		std::cout<<std::endl;
	}
	std::cout<<std::endl;
}
matrix::~matrix()
{
	for(int i=0;i<row;i++)
		delete[] p[i];
}
matrix operator+(const matrix& lhs,const matrix& rhs)
{
	//比较两个矩阵的size
	if(lhs.row!=rhs.row||lhs.col!=rhs.col)
		return matrix(NULL,0,0);
	int numbel=rhs.row*rhs.col;
	int* c=new int[numbel];
	for(int i=0;i<lhs.row;i++)
		for(int j=0;j<lhs.col;j++)
			c[i*lhs.col+j]=lhs.p[i][j]+rhs.p[i][j];
	return matrix(c,lhs.row,lhs.col);
}
matrix operator-(const matrix& lhs,const matrix& rhs)
{
	//比较两个矩阵的size
	if(lhs.row!=rhs.row||lhs.col!=rhs.col)
		return matrix(NULL,0,0);
	int numbel=rhs.row*rhs.col;
	int* c=new int[numbel];
	for(int i=0;i<lhs.row;i++)
		for(int j=0;j<lhs.col;j++)
			c[i*lhs.col+j]=lhs.p[i][j]-rhs.p[i][j];
	return matrix(c,lhs.row,lhs.col);
}
matrix operator*(const matrix& lhs,const matrix& rhs)
{
	if(lhs.col!=rhs.row)
		return matrix(NULL,0,0);
	int numbel=lhs.row*rhs.col;
	int* c=new int[numbel];
	for(int i=0;i<lhs.row;i++)
		for(int j=0;j<rhs.col;j++)
		{
			c[i*rhs.col+j]=0;
			for(int k=0;k<lhs.col;k++)
				c[i*rhs.col+j]+=lhs.p[i][k]*rhs.p[k][j];
		}
	return matrix(c,lhs.row,rhs.col);
}
matrix matrix::split_quarter(int m,int n) const 
{
	if(m<1||m>3||n<1||n>3)
		return matrix(NULL,0,0);
	int numbel=row*col/4;
	int r=row/2;
	int c=col/2;
	int* pc=new int[numbel];
	for(int i=0;i<r;i++)
		for(int j=0;j<c;j++)
			pc[i*c+j]=p[i+c*(m-1)][j+r*(n-1)];
	return matrix(pc,r,c);
}
matrix zip_quarter
	(const matrix& c11,const matrix& c12,
	const matrix& c21,const matrix& c22)
{
	int r=c11.row;
	int c=c11.col;
	int row=c11.row*2;
	int col=c11.col*2;
	int numbel_quarter=r*c;
	int numbel=row*col;
	int* pc=new int[numbel];
	for(int i=0;i<r;i++)
		for(int j=0;j<c;j++)
			pc[(2*i+0)*c+j]=c11.p[i][j];
	for(int i=0;i<r;i++)
		for(int j=0;j<c;j++)
			pc[(2*i+1)*c+j]=c12.p[i][j];
	for(int i=0;i<r;i++)
		for(int j=0;j<c;j++)
			pc[(2*i+4)*c+j]=c21.p[i][j];
	for(int i=0;i<r;i++)
		for(int j=0;j<c;j++)
			pc[(2*i+5)*c+j]=c22.p[i][j];
	return matrix(pc,row,col);

}
matrix strassen(const matrix& lhs,const matrix& rhs)
{
	int n=lhs.row;
	int* pc=new int[n*n];
	matrix c(pc,n,n);
	if(n==1)
	{	
		c.p[0][0]=lhs.p[0][0]*rhs.p[0][0];
	}
	else
	{
		matrix a11=lhs.split_quarter(1,1);
		matrix a12=lhs.split_quarter(1,2);
		matrix a21=lhs.split_quarter(2,1);
		matrix a22=lhs.split_quarter(2,2);
		matrix b11=rhs.split_quarter(1,1);
		matrix b12=rhs.split_quarter(1,2);
		matrix b21=rhs.split_quarter(2,1);
		matrix b22=rhs.split_quarter(2,2);

		matrix s1=b12-b22;
		matrix s2=a11+a12;
		matrix s3=a21+a22;
		matrix s4=b21-b11;
		matrix s5=a11+a22;
		matrix s6=b11+b22;
		matrix s7=a12-a22;
		matrix s8=b21+b22;
		matrix s9=a11-a21;
		matrix s10=b11+b12;

		matrix p1=strassen(a11,s1);
		matrix p2=strassen(s2,b22);
		matrix p3=strassen(s3,b11);
		matrix p4=strassen(a22,s4);
		matrix p5=strassen(s5,s6);
		matrix p6=strassen(s7,s8);
		matrix p7=strassen(s9,s10);

		matrix c11=p5+p4-p2+p6;
		matrix c12=p1+p2;
		matrix c21=p3+p4;
		matrix c22=p5+p1-p3-p7;
		c=zip_quarter(c11,c12,c21,c22);
	}
	return c;//c.p被销毁了！
}
#endif

chapter4_test.cpp

#include"chapter4.h"
#include<iostream>
void main()
{
	int a[]={13,-3,-25,20,-3,-16,-23,
		18,20,-7,12,-5,-22,15,-4,7};
	int b[]={13,-3,-25,20,-3,-16,-23,
		18,20,-7,12,-5,-22,15,-4,7};
	//int a[]={-1,2,3,2,-1};
	//max_array ma=find_maximum_subarray(a,0,15);
	//max_array ma=find_max_subarray_brute(a,0,15);
	//max_array ma=find_max_subarray_linear(a,0,15);
	//ma.display();
	matrix ma(a,4,4);
	matrix mb(b,4,4);
	//matrix md=ma*mb;
	//md.display();
	matrix mc=strassen(ma,mb);
	mc.display();

}