#include
#include

#define?N?10
#define?L?50
double?a[N][N];

double?Q[N][N]R[N][N]RQ[N][N];
double?lamda[N]Re[N]Im[N];
double?M[N][N]QK[N][N];
double?biaozhi;
int?mcnt;

void?initial（void）;				//初始化矩陣A
void?nishangsanjiaohua（void）;	//擬上三角化矩陣A
void?QR_resolve（void）;			//QR分解擬上三角陣A（n-1）
void?eig（void）;					//求矩陣的特征向量
void?QR_resolveM（void）;			//對矩陣Mk進行QR分解
void?solve（void）;				//解一元二次方程
void?jiexiangliang（void）;		//求對應于實數特征值的特征向量



void?main（void）
{
	double?b[N][N];
	int?ijl;

	initial（）;			//初始化矩陣A
	
//顯示矩陣A
	printf（“初始矩陣A顯示如下\n“）;
	for（i=0;i<=N-1;i++）
	{
		for（j=0;j<=N-1;j++）
		{
			printf（“a[%d][%d]?=?%.11e		“i+1j+1a[i][j]）;
		}
		printf（“\n“）;
	}

	nishangsanjiaohua（）;		//對矩陣A進行擬上三角化得到A（n-1）

//顯示擬上三角化后的矩陣A（n-1）
	printf（“擬上三角化后矩陣A（n-1）顯示如下\n“）;
	for（i=0;i<=N-1;i++）
	{
		for（j=0;j<=N-1;j++）
		{
			printf（“a[%d][%d]?=?%.11e		“i+1j+1a[i][j]）;
		}
		printf（“\n“）;
	}

	for（i=0;i<=N-1;i++）
		for（j=0;j<=N-1;j++）
			b[i][j]=a[i][j];

	QR_resolve（）;		//對A（n-1）進行QR分解，得矩陣Q，R，RQ


//顯示矩陣Q
	printf（“對A（n-1）進行QR分解后矩陣Q顯示如下\n“）;	
	for（i=0;i<=N-1;i++）
	{
		for（j=0;j<=N-1;j++）
			printf（“Q[%d][%d]?=?%.11e		“i+1j+1Q[i][j]）;
		printf（“\n“）;
	}

//顯示矩陣R
	printf（“對A（n-1）進行QR分解后矩陣R顯示如下\n“）;	
	for（i=0;i<=N-1;i++）
	{
		for（j=0;j<=N-1;j++）
		{
			R[i][j]=a[i][j];
			printf（“R[%d][%d]?=?%.11e		“i+1j+1R[i][j]）;
		}
		printf（“\n“）;
	}

//顯示矩陣RQ
	printf（“對A（n-1）進行QR分解后矩陣RQ顯示如下\n“）;
	for（i=0;i<=N-1;i++）
		for（j=0;j<=N-1;j++）
		{
			RQ[i][j]=0;
			for（l=0;l<=N-1;l++）
				RQ[i][j]+=R[i][l]*Q[l][j];
		}

	for（i=0;i<=N-1;i++）
		for（j=0;j<=N-1;j++）
		{
			printf（“RQ[%d][%d]?=?%.11e		“i+1j+1RQ[i][j]）;
		}

	for（i=0;i<=N-1;i++）
		for（j=0;j<=N-1;j++）
			a[i][j]=b[i][j];

//求出所有的特征值
	eig（）;
//顯示特征值
	printf（“所有特征值顯示如下:\n“）;
	for（i=0;i<=N-1;i++）
		printf（“λ%d（%.11e??%.11e）\n“i+1Re[i]Im[i]）;

	initial（）;		//初始化矩陣A，求A的特征向量

//求對應于實特征值的特征向量
	jiexiangliang（）;

}



//初始化矩陣A
void?initial（void）
{
	int?ij;
	for（i=0;i<=N-1;i++）
		for（j=0;j<=N-1;j++）
		{
			if（i==j）
				a[i][j]=1.5*cos（（i+1）+1.2*（j+1））;
			else
				a[i][j]=sin（0.5*（i+1）+0.2*（j+1））;
		}
}

//對矩陣A進行擬上三角化得到A（n-1）
void?nishangsanjiaohua（void）
{
		
	int?ijr;
	double?sumcdht;
	double?u[N]v[N]p[N]q[N]w[N];

	for（r=0;r<=N-3;r++）
	{
		for（i=r+2;i<=N-1;i++）
			if（fabs（a[i][r]）>=1e-12）break;
		if（i==N?&&?fabs（a[N-1][r]）<=1e-12?）
			continue;

		sum=0;
		for（i=r+1;i<=N-1;i++）
			sum+=a[i][r]*a[i][r];
		d=sqrt（sum）;
		if（a[r+1][r]>1e-12）
			c=-d;
		else
			c=d;
		h=c*c-c*a[r+1][r];

		for（i=0;i<=r;i++）
			u[i]=0;
		u[r+1]=a[r+1][r]-c;
		for（i=r+2;i<=N-1;i++）
			u[i]=a[i][r];

		for（i=0;i<=N-1;i++）
			v[i]=u[i]/h;
		for（i=0;i<=N-1;i++）
		{
			p[i]=0;
			q[i]=0;
			for（j=0;j<=N-1;j++）
			{
				p[i]+=a[j][i]*v[j];
				q[i]+=a[i][j]*v[j];
			}
		}

		t=0;
		for（i=0;i<=N-1;i++）
			t+=p[i]*v[i];
		for（i=0;i<=N-1;i++