matrix.cpp

基于线性规划的回归支持向量机源程序

	int len=m1.p->length;
	int ven=vlen(m1.p->f[0]);

    matrix m(len,ven);
	//do the operations
    vector **f=m.p->f;
    vector **f1=m1.p->f;
	for(int i=0; i<len; i++)
		*f[i]=-*f1[i] ;
	return m;
}

//m+=m1

matrix  &  matrix::operator+=(matrix  & m1)
{
    //allocate the new matrix
    if(m1!=(*this))
    {
        cout<<"matrix's dim are not equal";
    }
	int len=MIN(p->length,m1.p->length);

	//do the operations
    vector **f=p->f;
    vector **f1=m1.p->f;
	for(int i=0; i<len; i++)
		*f[i] +=*f1[i] ;
	return *this;
}


//m-=m1

matrix  &  matrix::operator-=(matrix  & m1)
{
    //allocate the new matrix
    if(m1!=(*this))
    {
        cout<<"matrix's dim are not equal";
    }
	int len=MIN(p->length,m1.p->length);

	//do the operations
    vector **f=p->f;
    vector **f1=m1.p->f;
	for(int i=0; i<len; i++)
		*f[i] -=*f1[i] ;
	return *this;
}


//m/=m1

matrix  &  matrix::operator/=(matrix  & m1)
{
    //allocate the new matrix
    if(m1!=(*this))
    {
        cout<<"matrix's dim are not equal";
    }
	int len=MIN(p->length,m1.p->length);

	//do the operations
    vector **f=p->f;
    vector **f1=m1.p->f;
	for(int i=0; i<len; i++)
		*f[i] /=*f1[i] ;
	return *this;
}


//m+=x

matrix  &  matrix::operator+=(double x)
{
    //allocate the new matrix
	int len=p->length;

	//do the operations
    vector **f=p->f;
	for(int i=0; i<len; i++)
		*f[i] +=x ;
	return *this;
}


//m-=x

matrix  &  matrix::operator-=(double x)
{
    //allocate the new matrix
	int len=p->length;

	//do the operations
    vector **f=p->f;
	for(int i=0; i<len; i++)
		*f[i] -=x ;
	return *this;
}

//m*=x

matrix  &  matrix::operator*=(double x)
{
    //allocate the new matrix
	int len=p->length;

	//do the operations
    vector **f=p->f;
	for(int i=0; i<len; i++)
		*f[i] *=x ;
	return *this;
}

//m/=x

matrix  &  matrix::operator/=(double x)
{
    //allocate the new matrix
	int len=p->length;

	//do the operations
    vector **f=p->f;
	for(int i=0; i<len; i++)
		*f[i] /=x ;
	return *this;
}


//m=m1.prod(m2)

matrix matrix::prod(matrix  & m2)
{
    //allocate the new matrix
    if(m2!=(*this))
    {
        cout<<"matrix's dim are not equal";
    }
	int len=MIN(p->length,m2.p->length);
	int ven=MIN(vlen(p->f[0]),vlen(m2.p->f[0]));

    matrix m(len,ven);
	//do the operations
    vector **f=m.p->f;
    vector **f1=p->f;
    vector **f2=m2.p->f;
	for(int i=0; i<len; i++)
		*f[i]=*f1[i] * *f2[i];
	return m;
}


//m=prod(m1,m2)

matrix  prod(matrix  & m1,matrix  & m2)
{
    //allocate the new matrix
    if(m1!=m2)
    {
        cout<<"matrix's dim are not equal";
    }
	int len=MIN(m1.p->length,m2.p->length);
	int ven=MIN(vlen(m1.p->f[0]),vlen(m2.p->f[0]));

    matrix m(len,ven);
	//do the operations
    vector **f=m.p->f;
    vector **f1=m1.p->f;
    vector **f2=m2.p->f;
	for(int i=0; i<len; i++)
		*f[i]=*f1[i] * *f2[i];
	return m;
}


//v=m[i]

vector  &  matrix::operator[](int i)
{
	if((i>=0)&&(i<p->length))
		return *p->f[i];
	else{
      {
        cout<< Merrorstring[EMAT_OUTOFRANGE];
        throw Merrorstring[EMAT_OUTOFRANGE];
    }
//next statement has no effect, the only function is avoid warnning message.
		return *p->f[0];
	}
}
double  &  matrix::operator()(int i,int j)
{
    if( (i>=0)&&(i<p->length)&&(j>=0)&&(j<vlen(p->f[0]) ) )
        return (*p->f[i])[j];
	else{
      {
        cout<< Merrorstring[EMAT_OUTOFRANGE];
        throw Merrorstring[EMAT_OUTOFRANGE];
    }
//next statement has no effect, the only function is avoid warnning message.
        return (*p->f[0])[0];
	}
}


//x=(m1==m2) :if m1==m2 returns 1. otherwise, returns 0.

int operator==(matrix  & m1, matrix  & m2)
{
	if((m1.p->length !=m2.p->length) ||
		(vlen(m1.p->f[0]) != vlen(m2.p->f[0])))
			return 0;
            else
//    int len=m1.p->length;
	//do the operations
//    vector **f1=m1.p->f;
//    vector **f2=m2.p->f;
//    for(int i=0; i<len; i++)
//        if(*f1[i] != *f2[i])
//            return 0;
	return 1;
}

//x=(m1 !=m2) :if m1!=m2 returns 1. otherwise, returns 0.

int operator!=(matrix  & m1, matrix  & m2)
{
	 return !(m2==m1);
}


//cin>>m1

istream  &  operator>>(istream  & in, matrix  & m1)
{
	int nrow=m1.mrow();
	int ncol=m1.mcol();
	for(int i=0; i<nrow; i++){
		for(int j=0; j<ncol; j++)
			in>>m1[i][j];
	}
		cout<<"\n";
	return in;
}

//cout<<m1

ostream  &  operator<<(ostream  & out, matrix  & m1)
{
    int row=mrow(m1);
    int col=mcol(m1);
    for(int i=0; i<row; i++)
    {
        for(int j=0; j<col; j++)
        out<<m1(i,j)<<" ";
        out<<"\n ";
    }
    return out;
}


//returns a special row of a matrix

vector matrix::row(int n)
{
    if((n>p->length)||(n<=0))
      {
        cout<< Merrorstring[EMAT_OUTOFRANGE];
        throw Merrorstring[EMAT_OUTOFRANGE];
    }
	return ::copy(*(p->f[n-1]));
}


vector row(matrix  & m,int n)
{
    if((n>m.mrow())||(n<=0))
      {
        cout<< Merrorstring[EMAT_OUTOFRANGE];
        throw Merrorstring[EMAT_OUTOFRANGE];
    }
	return copy(*(m.p->f[n-1]));
}


//returns a special column of a matrix

vector matrix::col(int n)
{
	int mrows=p->length;
	int mcols=vlen(p->f[0]);
    if((n>mcols)||(n<=0))
      {
        cout<< Merrorstring[EMAT_OUTOFRANGE];
        throw Merrorstring[EMAT_OUTOFRANGE];
    }
    vector x(mrows);
	for(int i=0; i<mrows; i++)
		x[i]=(*this)[i][n-1];
	return x;
}


vector col(matrix  & m, int n)
{
	int mrows=m.mrow();
	int mcols=m.mcol();
    if((n>mcols)||(n<=0))
      {
        cout<< Merrorstring[EMAT_OUTOFRANGE];
        throw Merrorstring[EMAT_OUTOFRANGE];
    }
    vector x(mrows);
	for(int i=0; i<mrows; i++)
		x[i]=m[i][n-1];
	return x;
}


//mresult=x.slice(n,m): returns a sub matrix of x. it is n*m and
//mresult[0][0]=x[0][0].

matrix matrix::slice(int row,int col)
{
	int mrows=p->length;
	int mcols=vlen(p->f[0]);
    if((row>mrows)||(col>mcols)||(row<0)||(col<0))
      {
        cout<< Merrorstring[EMAT_OUTOFRANGE];
        throw Merrorstring[EMAT_OUTOFRANGE];
    }
    matrix x(row,col);
	for(int i=0; i<row; i++)
		for(int j=0; j<col; j++)
			x[i][j]=(*this)[i][j];
	return x;
}


matrix slice(matrix  & m,int row,int col)
{
	int mrows=m.mrow();
	int mcols=m.mcol();
    if((row>mrows)||(col>mcols)||(row<0)||(col<0))
      {
        cout<< Merrorstring[EMAT_OUTOFRANGE];
        throw Merrorstring[EMAT_OUTOFRANGE];
    }
    matrix x(row,col);
	for(int i=0; i<row; i++)
		for(int j=0; j<col; j++)
			x[i][j]=m[i][j];
	return x;
}

//mresult=x.slice(n1,m1,n2,m2): returns a submatrix of x. it is
//n2*m2 and mresult[0][0]=x[n1][m1].

matrix matrix::slice(int row1,int col1,
				int row2,int col2)
{
	int mrows=p->length;
	int mcols=vlen(p->f[0]);
    if((row2<0)||(col2<0))
      {
        cout<<Merrorstring[EMAT_ASSIGNDATAERR];
        throw Merrorstring[EMAT_ASSIGNDATAERR];
    }
    if((row1+row2>mrows)||(col1+col2>mcols)||(row1<0)||(col1<0))
      {
        cout<< Merrorstring[EMAT_OUTOFRANGE];
        throw Merrorstring[EMAT_OUTOFRANGE];
    }
    matrix x(row2,col2);
	for(int i=row1; i<row1+row2; i++)
		for(int j=col1; j<col1+col2; j++)
			x[i-row1][j-col1]=(*this)[i][j];
	return x;
}


matrix slice(matrix  & m,int row1,
				int col1,int row2,int col2)
{
	int mrows=m.mrow();
	int mcols=m.mcol();
    if((row2<0)||(col2<0))
      {
        cout<<Merrorstring[EMAT_ASSIGNDATAERR];
        throw Merrorstring[EMAT_ASSIGNDATAERR];
    }
    if((row1+row2>mrows)||(col1+col2>mcols)||(row1<0)||(col1<0))
      {
        cout<< Merrorstring[EMAT_OUTOFRANGE];
        throw Merrorstring[EMAT_OUTOFRANGE];
    }
    matrix x(row2,col2);
	for(int i=row1; i<row1+row2; i++)
		for(int j=col1; j<col1+col2; j++)
			x[i-row1][j-col1]=m[i][j];
	return x;
}


matrix commat(matrix  & m1,
                    matrix  & m2, int select)
{
	int i,j,m,n,minc,ninc;
	int mrow1=m1.mrow();
	int mcol1=m1.mcol();
	int mrow2=m2.mrow();
	int mcol2=m2.mcol();
	switch(select)
	{
		case 0:
			m=mrow1+mrow2;
			n=mcol1+mcol2;
			minc=mrow1;
			ninc=mcol1;
			break;
		case 1:
			m=MAX(mrow1,mrow2);
			n=mcol1+mcol2;
			minc=0;
			ninc=mcol1;
			break;
		case 2:
			m=mrow1+mrow2;
			n=MAX(mcol1,mcol2);
			minc=mrow1;
			ninc=0;
			break;
		default:
      {
        cout<<Merrorstring[EMAT_ASSIGNDATAERR];
        throw Merrorstring[EMAT_ASSIGNDATAERR];
    }
            ;
	}
    matrix x(m,n);
	for(i=0; i<mrow1; i++)
		for(j=0; j<mcol1; j++)
			x[i][j]=m1[i][j];
	for(i=0; i<mrow2; i++)
		for(j=0; j<mcol2; j++)
			x[i+minc][j+ninc]=m2[i][j];
	return x;
}


matrix matrix::commat(matrix  & m2,
				 int select)
{
		return ::commat(*this,m2,select);
}


matrix commat(matrix  & m1,
                    vector  & m2, int select)
{
	int i,j,m,n;
	int mrow1=m1.mrow();
	int mcol1=m1.mcol();
	int mrow2=m2.vlen();
	switch(select)
	{
		case 1:
			m=MAX(mrow1,mrow2);
			n=mcol1+1;
			break;
		case 2:
			m=mrow1+1;
			n=MAX(mcol1,mrow2);
			break;
		default:
      {
        cout<<Merrorstring[EMAT_ASSIGNDATAERR];
        throw Merrorstring[EMAT_ASSIGNDATAERR];
    }
            ;
	}
    matrix x(m,n);
	for(i=0; i<mrow1; i++)
		for(j=0; j<mcol1; j++)
			x[i][j]=m1[i][j];
	if(select==1){
		for(i=0; i<mrow2; i++)
			x[i][mcol1]=m2[i];
	}
	else{
		for(j=0; j<mrow2; j++)
			x[mrow1][j]=m2[j];
	}

	return x;
}


matrix matrix::commat(vector  & m2,
				 int select)
{
		return ::commat(*this,m2,select);
}



matrix commat(vector  & m1,
                    matrix  & m2, int select)
{
	int i,j,m,n;
	int mrow1=m1.vlen();
	int mrow2=m2.mrow();
	int mcol2=m2.mcol();
	switch(select)
	{
		case 1:
			m=MAX(mrow1,mrow2);
			n=mcol2+1;
			break;
		case 2:
			m=mrow2+1;
			n=MAX(mrow1,mcol2);
			break;
		default:
      {
        cout<<Merrorstring[EMAT_ASSIGNDATAERR];
        throw Merrorstring[EMAT_ASSIGNDATAERR];
    }
            ;
	}
    matrix x(m,n);
	if(select==1){
		for(i=0; i<mrow1; i++)
			x[i][0]=m1[i];
		for(i=0; i<mrow2; i++)
			for(j=0; j<mcol2; j++)
				x[i][j+1]=m2[i][j];
	}
	else{
		for(j=0; j<mrow1; j++)
			x[0][j]=m1[j];
		for(i=0; i<mrow2; i++)
			for(j=0; j<mcol2; j++)
				x[i+1][j]=m2[i][j];
	}

	return x;
}



matrix commat(vector  & m1,
                    vector  & m2, int select)
{
	int i,m,n;
	int mrow1=m1.vlen();
	int mrow2=m2.vlen();
    matrix tmp;
	switch(select)
	{
		case 1:
			m=MAX(mrow1,mrow2);
			n=2;
			break;
		case 2:
			m=2;
			n=MAX(mrow1,mrow2);
			break;
		default:
      {
        cout<<Merrorstring[EMAT_ASSIGNDATAERR];
        throw Merrorstring[EMAT_ASSIGNDATAERR];
    }
            ;
	}
    matrix x(m,n);
	if(select==1){
			for(i=0; i<mrow1; i++)
				x[i][0]=m1[i];
			for(i=0; i<mrow2; i++)
				x[i][1]=m2[i];
	}
	else{
			for(i=0; i<mrow1; i++)
				x[0][i]=m1[i];
			for(i=0; i<mrow2; i++)
				x[1][i]=m2[i];
	}
	return x;
}



matrix commat(matrix  & m1,
                    matrix  & m2, int srow, int scol)
{
	int i,j;
	int mrow1=m1.mrow();
	int mcol1=m1.mcol();
	int mrow2=m2.mrow();
	int mcol2=m2.mcol();
    if((srow<mrow1)&&(scol<mcol1))
      {
        cout<<Merrorstring[EMAT_ASSIGNDATAERR];
        throw Merrorstring[EMAT_ASSIGNDATAERR];
    }
	int mx=MAX(mrow1,srow+mrow2);
	int my=MAX(mcol1,scol+mcol2);
    matrix x(mx,my);
	for(i=0; i<mrow1; i++)
		for(j=0; j<mcol1; j++)
			x[i][j]=m1[i][j];
	for(i=0; i<mrow2; i++)
		for(j=0; j<mcol2; j++)
			x[i+srow][j+scol]=m2[i][j];
	return x;
}



matrix matrix::commat(matrix  & m2,
			 int srow, int scol)
{
	return ::commat(*this,m2,srow,scol);
}



matrix diag(vector  & m2)
{
	int n=vlen(m2);
    matrix x(n,n,0);
	for(int i=0; i<n; i++)
		x[i][i]=m2[i];
	return x;
}


matrix diag(vector  & m2,int k)
{
	int nn=vlen(m2);
	int kk=abs(k);
	int n=nn+kk;
    matrix x(n,n,0);
	if(k>0){
		for(int i=0; i<nn; i++)
			x[i][i+k]=m2[i];
	}
	else
	{
		for(int i=kk; i<n; i++)
			x[i][i-kk]=m2[i-kk];
	}
	return x;
}


vector diag(matrix  & m2)
{
	int n=MIN(m2.mrow(),m2.mcol());
    vector x(n);
	for(int i=0; i<n; i++)
		x[i]=m2[i][i];
	return x;
}


vector matrix::diag()
{
	int n=MIN(mrow(),mcol());
    vector x(n);
	for(int i=0; i<n; i++)
		x[i]=(*this)[i][i];
	return x;
}



vector diag(matrix  & m2,int k)
{
	int m=m2.mrow();
	int n=m2.mcol();
	int kk=abs(k);
	int nn,i;

    if( ((k>0)&&(k>(n-1)))||((k<0)&&(-k>(m-1))) )
      {
        cout<<Merrorstring[EMAT_ASSIGNDATAERR];
        throw Merrorstring[EMAT_ASSIGNDATAERR];
    }
	if(m>=n){
		if(k>=0)
			nn=n-kk;
		else if(n+kk<=m)
				nn=n;
			 else
				nn=m-kk;
	}
	else{
		if(k>=0)
			if(m+kk<=n)
				nn=m;
			else
				nn=n-kk;
		else
			nn=m-kk;
	}
    vector x(nn);
	if(k>0)
		for(i=0; i<nn; i++)
			x[i]=m2[i][i+kk];
	else