mat.cpp

矩阵类,自己编写,通过较多考验,符合matlab编程习惯,可处理无穷维矩阵和向量

	int m=out.row;
	int n=out.col;
	for (int i=1;i<=m;i++)
	{
		for (int j=1;j<=n;j++)
		{
			if ( i< rr)
			{
				if ( j< cc)
				{
					out.set(i,j,r(i,j));
				}
				else
				{
					out.set(i,j,r(i,j+1));
				}
			}
			else
			{
				if (i>=rr)
				{
					if ( j< cc)
					{
						out.set(i,j,r(i+1,j));
					}
					else
					{
						out.set(i,j,r(i+1,j+1));
					}
				}
			}
		}
	}
	return out;
}

//--------------------逆矩阵--------------------------------
mat mat::invmat() const
{
	mat out;
	out.zero(row,col);
	if (col!=row)
	{
		exit(1);
	}
	else
	{
		double d=this->detmat();
		if (d !=0.0)
		{
			double d1=1/d;
			mat rusult=adj();
			int m=rusult.row;
			int n=rusult.col;
			for (int i=1;i<=m;i++)
			{
				for (int j=1;j<=n;j++)
				{
					out.set(i,j,rusult.r(i,j)*d1);			
				}
			}
		}
	}
	return out;
}


//---------单变量数值积分-------------------
double mat::Intf(double dydx_new,double dydx,double y,double int_step) const
{	
	// dydx_new 当前速度
	// dydx     上一时刻速度
	// y        上一时刻位移
	// 输出     当前位移
	int const n=8;
	double z[n+1];
	double h=int_step/n;
	double new_df;
	z[0]=y;	z[1]=z[0]+h*dydx;
	for(int i=1;i<=n-1;i++)
	{
		new_df=dydx+i*((dydx_new-dydx)/n);
		z[i+1]=z[i-1]+2*h*new_df;
	}
	return 0.5*(z[n]+z[n-1]+h*dydx_new);
}

//---------列向量数值积分-------------------
mat mat::intm(const mat& dydx_new, const mat& dydx, const mat& y, double int_step) const
{
	int m =dydx_new.row;
	int m1=dydx.row;
	int m2=y.row;
	int n =dydx_new.col;
	int n1=dydx.col;
	int n2=y.col;
	mat out;out.zero(m,1);
	if ((m==m1)&&(m==m2)&&(n==1)&&(n1==1)&&(n2==1))
	{
		for (int i=1;i<=m;i++)
		{
			out.set(i,1,Intf(dydx_new.r(i,1),dydx.r(i,1),y.r(i,1),int_step));
		}
	}
	return out;
}

//---------开辟线性空间矢量-----------
void mat::linespace(double start, double step, double end)
{
	double m=(end-start)/step;
	int n;
	n=int(ceil(m-0.5));
	zero(1,n+1);
	for (int i=1;i<=n+1;i++)
	{
		this->set(1,i,start+(i-1)*step);
	}
}

//------2个线性空间的结合----------------
void mat::operator||(mat a)
{	
	int m1=a.row;	int n1=a.col;
	int m=row;      int n=col;
	if (m1 != m)
	{
		exit(0);
	}
	else
	{
		mat temp;temp.zero(m,n);
		for (int i=1;i<=m;i++)
		{
			for (int j=1;j<=n;j++)
			{
				temp.set(i,j,r(i,j));
			}
		}
		zero(m1,n+n1);		
		for (int ii=1;ii<=row;ii++)
		{
			for (int jj=1;jj<=col;jj++)
			{
				if (jj<=n)
					this->set(ii,jj,temp.r(ii,jj));
				else
					this->set(ii,jj,a.r(ii,jj-n));
			}
		}
	} 
}

//------2个线性空间的结合----------------
void mat::operator||(double a)
{
	int m=row;      int n=col;
	if (m  >= 1.5 )
	{
		exit(0);
	}
	else
	{
		mat temp;temp.zero(1,n);
		for (int j=1;j<=n;j++)
		{
			temp.set(1,j,r(1,j));
		}
		zero(m,n+1);
		for (int jj=1;jj<=col;jj++)
		{
			if (jj<=n)
				this->set(1,jj,temp.r(1,jj));
			else
				this->set(1,jj,a);
		}
	} 
}

//--------------按照列截断段，得到包含该列的前面所有数据---------
void mat::Jieduan(int a)
{
	mat temp(*this);	
/*	temp.zero(row,a);
	for(int i=1;i<=row;i++)
	{
		for (int j=1;j<=a;j++)
		{
			temp.set(i,j,r(i,j));
		}
	}
*/
	zero(row,a);
	for(int ii=1;ii<=row;ii++)
	{
		for (int jj=1;jj<=a;jj++)
		{
			this->set(ii,jj,temp.r(ii,jj));
		}
	}
}

//----------将数组赋值给矩阵--------------
void mat::setarray(double d[])
{
	int n=int(d[0]);
	if (n<=0.0)
	{
		zero(1,1);
	}
	else
	{
		zero(1,n);
		for (int i=1;i<=n;i++)
		{
			this->set(1,i,d[i]);
		}
	}
}

//----------插值--------------
double mat::inperp(mat Value,double t) const
{
	mat a=Value;
	int m=a.row;
	int n=a.col;
	int in;
	double out;
	if ((m != row) && (n !=col))
	{
		out=0.0;
	}
	else
	{
		if (m == 1)
		{
			if ( t<= r(1,1))
			{
				out= a.r(1,1);
			}
			else
			{
				if ( t>=r(1,col))
				{	
					out=a.r(1,col);
				}
				else
				{
					for (int i=1;i<col;i++)
					{
						if ( (t>r(1,i)) && (t <= r(1,i+1)) )
						{
							in=i;
							break;
						}
					}
					out=(t-r(1,in)) * ( a.r(1,in+1)-a.r(1,in) )/( r(1,in+1) -r(1,in) )+a.r(1,in);
				}
			}
		}
		else
		{
			if (n == 1.0 )
			{
				if ( t<=r(1,1))
				{
					out=a.r(1,1);
				}
				else 
				{
					if ( t>=r(row,1))
					{
						out=a.r(row,1);
					}
					else
					{
						for (int i=1;i<row;i++)
						{
							if ( (t>r(i,1)) && (t <= r(i+1,1)) )
							{
								in=i;
								break;
							}
						}
						out=(t-r(in,1)) * ( a.r(in+1,1)-a.r(in,1) )/( r(in+1,1) - r(in,1) )+a.r(in,1);
					}
				}
			}
			else
			{
				out=0.0;
			}
		}
	}
	return out;
}

int mat::length()
{
	return col*row;
}

//----------求元素中的最大值，并返回所在行和列---------
mat mat::maxmat() const
{
	double out=r(1,1);
	int m=1;
	int n=1;
	int i,j;
	for (i=1;i<=row;i++)
	{
		for(j=1;j<=col;j++)
		{
			if ( out <= r(i,j) )
			{
				out=r(i,j);
				m=i;
				n=j;
			}
		}
	}
	mat put;put.zero(3,1);
	put.set(1,1,out);
	put.set(2,1,m);
	put.set(3,1,n);
	return put;
}

//----------求除去某元素后的的最大值，并返回所在行和列---------
mat mat::max1mat(int m, int n) const
{ 
	double out;
	if( (m==1) && (n==1))
	{
		out=r(1,2);
	}
	else
	{
		out=r(1,1);
	}
	int e=1;
	int f=1;
	for (int i=1;i<=row;i++)
	{
		for(int j=1;j<=col;j++)
		{
			if( (i==m ) && (j==n ) )
			{
				out=out;
			}
			else
			{
				if ( out <= r(i,j) )
				{
					out=r(i,j);
					e=i;
					f=j;
				}
			}
		}
	}
	mat put;put.zero(3,1);
	put.set(1,1,out);
	put.set(2,1,e);
	put.set(3,1,f);
	return put;
}

//-------------按低到高排列----------------------------
mat mat::Sortmat() const
{
	int m=row;
	int n=col;
	mat out;out.zero(m,n);
	mat maxx;	
	mat temp;
	mat maxtemp;
	if ( m == 1)
	{
		maxx=maxmat();
		out.set(1,n,maxx.r(1,1));
		temp=Delc(int(maxx.r(3,1)));
		for (int i=n-1; i>=1;i--)
		{			
			maxtemp=temp.maxmat();
			if (i>=2)
			{
				out.set(1,i,maxtemp.r(1,1));
				temp=temp.Delc(int(maxtemp.r(3,1)));
			}
			else
			{
				out.set(1,1,temp.r(1,1));
			}
		}
	}
	else
	{
		if (n ==1 )
		{
			maxx=maxmat();
			out.set(m,1,maxx.r(1,1));
			temp=Delr(int(maxx.r(2,1)));
			for (int i=m-1; i>=1;i--)
			{
				maxtemp=temp.maxmat();
				if (i>=2)
				{
					out.set(i,1,maxtemp.r(1,1));
					temp=temp.Delr(int(maxtemp.r(2,1)));
				}
				else
				{
					out.set(1,1,temp.r(1,1));
				}
			}
		}
	}
	return out;
}

//-------------删除某一列----------------------------
mat mat::Delc(int n) const
{
	mat out;out.zero(row,col-1);
	if ( n <= col)
	{
		for (int i=1;i<=row;i++)
		{
			for (int j=1;j<=col-1;j++)
			{
				if (j < n)
					out.set(i,j,r(i,j));
				else
					out.set(i,j,r(i,j+1));
			}
		}		
	}
	return out;
}

//-------------删除某一行----------------------------
mat mat::Delr(int n) const
{
	mat out;out.zero(row-1,col);
	if ( n <= row)
	{
		for (int j=1;j<=col;j++)
		{
			for (int i=1;i<=row-1;i++)
			{
				if (i < n)
					out.set(i,j,r(i,j));
				else
					out.set(i,j,r(i+1,j));
			}
		}		
	}
	return out;
}