北航数值分析作业第一题.doc

数值分析作业第一题

算法设计方案

利用带状Dollittle分解，将A[501][501]转存到数组C[5][501]，以节省存储空间

1、计算λ1和λ501

首先使用幂法求出矩阵的按模最大的特征值λ0：如果λ00，则其必为按模最大值，因此λ501=λ0，然后采用原点平移法，平移量为λ501，使用幂法迭代求出矩阵A-λ501I的按模最大的特征值，其特征值按从小到大排列应为λ1-λ501、λ2-λ501、……、。因此A-λ501I的按模最大的特征值应为λ1-λ501。又因为λ501的值已求得，由此可直接求出λ1。

2、计算λS

λS为矩阵A按模最小的特征值，可以通过反幂法直接求出。

3、计算λik

λik是对矩阵A进行λik平移后，再用反幂法求出按模最小的特征值λmin，λik=λik+λmin。

4、计算矩阵A的条件数计算cond（A）2和行列式det(A)

矩阵A的条件数为，其中λ1和λn分别是矩阵A的模最大和最小特征值，直接利用上面求得的结果直接计算。

矩阵A的行列式可先对矩阵A进行LU分解后，det(A)等于U所有对角线上元素的乘积。

二、源程序：

#includemath.h

#includestdio.h

#includestdlib.h

#includeiostream.h

#defines2

#definer2

intMax(intv1,intv2);

intMin(intv1,intv2);

intmaxt(intv1,intv2,intv3);

voidstorage(doubleC[5][501],doubleb,doublec);

doublemifa(doubleC[5][501]);

voidLU(doubleC[5][501]);

doublefmifa(doubleC[5][501]);

intMax(intv1,intv2)//求两个数的最大值

{ return((v1v2)?v1:v2);

}

intMin(intv1,intv2)//求两个数最小值

{ return((v1v2)?v1:v2);

}

intmaxt(intv1,intv2,intv3)//求三个数最大值

{intt;

if(v1v2)t=v1;

elset=v2;

if(tv3)t=v3;

return(t);

}

/***将矩阵值转存在一个数组里，以节省存储空间***/

voidstorage(doubleC[5][501],doubleb,doublec)

{ inti=0,j=0;

C[i][j]=0,C[i][j+1]=0;

for(j=2;j=500;j++)

C[i][j]=c;

i++;

j=0;

C[i][j]=0;

for(j=1;j=500;j++)

C[i][j]=b;

i++;

for(j=0;j=500;j++)

C[i][j]=(1.64-0.024*(j+1))*sin(0.2*(j+1))-0.64*exp(0.1/(j+1));

i++;

for(j=0;j=499;j++)

C[i][j]=b;

C[i][j]=0;

i++;

for(j=0;j=498;j++)

C[i][j]=c;

C[i][j]=0,C[i][j+1]=0;

}

//用于求解最大的特征值，幂法

doublemifa(doubleC[5][501])

{ intm=0,i,j;

doubleb2,b1=0,sum;

doubleu[501],y[501];

for(i=0;i501;i++)

{u[i]=1.0;

}

do

{ sum=0;

if(m!=0)b1=b2;

m++;

for(i=0;i=500;i++)

sum+=u[i]*u[i];

for(i=0;i=500;i++)

y[i]=u[i]/sqrt(sum);

for(i=0;i=500;i++)

{u[i]=0;

for(j=Max(i-r,0);j=Min(i+s,500);j++)

u[i]=u[i]+C[i-j+s][j]*y[j];

}

b2=0;

for(i=0;i=500;i++)

b2=b2+y[i]*u[i];

}

while(fabs(b2-b1)/fabs(b2)=1.0e-12)