北航数值分析第一次大作业(幂法反幂法)

北航数值分析第一次大作业(幂法反幂法)
北航数值分析第一次大作业(幂法反幂法)

北航数值分析大作业第一题幂法与反幂法

《数值分析》计算实习题目 第一题: 1. 算法设计方案 (1)1λ,501λ和s λ的值。 1)首先通过幂法求出按模最大的特征值λt1,然后根据λt1进行原点平移求出另一特征值λt2,比较两值大小,数值小的为所求最小特征值λ1,数值大的为是所求最大特征值λ501。 2)使用反幂法求λs ,其中需要解线性方程组。因为A 为带状线性方程组,此处采用LU 分解法解带状方程组。 (2)与140k λλμλ-5011=+k 最接近的特征值λik 。 通过带有原点平移的反幂法求出与数k μ最接近的特征值 λik 。 (3)2cond(A)和det A 。 1)1=n λλ2cond(A),其中1λ和n λ分别是按模最大和最小特征值。 2)利用步骤(1)中分解矩阵A 得出的LU 矩阵,L 为单位下三角阵,U 为上三角阵,其中U 矩阵的主对角线元素之积即为det A 。 由于A 的元素零元素较多,为节省储存量,将A 的元素存为6×501的数组中,程序中采用get_an_element()函数来从小数组中取出A 中的元素。 2.全部源程序 #include #include void init_a();//初始化A double get_an_element(int,int);//取A 中的元素函数 double powermethod(double);//原点平移的幂法 double inversepowermethod(double);//原点平移的反幂法 int presolve(double);//三角LU 分解 int solve(double [],double []);//解方程组 int max(int,int); int min(int,int); double (*u)[502]=new double[502][502];//上三角U 数组 double (*l)[502]=new double[502][502];//单位下三角L 数组 double a[6][502];//矩阵A int main() { int i,k; double lambdat1,lambdat2,lambda1,lambda501,lambdas,mu[40],det;

北航数值分析大作业第二题

数值分析第二次大作业 史立峰 SY1505327

一、 方案 (1)利用循环结构将sin(0.50.2)() 1.5cos( 1.2)() {i j i j ij i j i j a +≠+==(i,j=1,2,……,10)进行赋值,得到需要变换的 矩阵A ; (2)然后,对矩阵A 利用Householder 矩阵进行相似变换,把A 化为上三角矩阵A (n-1)。 对A 拟上三角化,得到拟上三角矩阵A (n-1),具体算法如下: 记A(1)=A ,并记A(r)的第r 列至第n 列的元素为()n r r j n i a r ij ,,1,;,,2,1) ( +==。 对于2,,2,1-=n r 执行 1. 若 ()n r r i a r ir ,,3,2) ( ++=全为零,则令A(r+1) =A(r),转5;否则转2。 2. 计算 () ∑+== n r i r ir r a d 1 2 )( ()( )r r r r r r r r r r d c a d a c ==-=++则取,0sgn ) (,1)(,1若 )(,12r r r r r r a c c h +-= 3. 令 () n T r nr r r r r r r r r R a a c a u ∈-=++) ()(,2)(,1,,,,0,,0 。 4. 计算 r r T r r h u A p /)(= r r r r h u A q /)(= r r T r r h u p t /= r r r r u t q -=ω T r r T r r r r p u u A A --=+ω)()1( 5. 继续。 (3)使用带双步位移的QR 方法计算矩阵A (n-1)的全部特征值,也是A 的全部特征值,具体算法如下: 1. 给定精度水平0>ε和迭代最大次数L 。 2. 记n n ij n a A A ?-==][) 1()1()1(,令n m k ==,1。

实验8 反幂法

《数值分析》实验8 一.实验名称:反幂法 二、实验目的: (1) 掌握求矩阵按模最小特征值及其对应的特征向量的规范化反幂法; (2) 掌握原点移位法。 三、实验要求 (1) 按照题目要求完成实验内容 (2) 写出相应的实验原理与C 语言程序 (3) 给出实验结果,结果分析 (4) 写出相应的实验报告 四、实验题目: 用反幂法求以下矩阵的指定特征值及其特征向量(迭代终止误差为1e-3): 41411014110?? ? ? ??? 接近2的特征值及其特征向量. 程序: #include #include int main(){ int n = 3, i, k, j, max_k = 1000; double y[max_k][n], x[max_k][n], max_x[max_k], z[n], u[n][n], l[n][n], t, err = 1e-5, norm_dx = 1,tt=12;//norm_dx:dx 的一范数,tt:为了求接近tt 的特征值, for (i = 0; i < n; i++) for (j = 0; j < n; j++){ u[i][j] = 0; (j == i) ? (l[i][j] = 1) : (l[i][j] = 0); } double a[][3] = { 4, 1, 4, 1, 10, 1, 4, 1, 10 }; for (i = 0; i < n; i++) { a[i][i]-=tt; x[0][i] = 1; } for (i = 0; i < n; i++){ for (k = i; k < n; k++){ t = 0; for (j = 0; j < i; j++) t += l[i][j] * u[j][k];

北航数值分析报告大作业第八题

北京航空航天大学 数值分析大作业八 学院名称自动化 专业方向控制工程 学号 学生姓名许阳 教师孙玉泉 日期2014 年11月26 日

一.题目 关于x , y , t , u , v , w 的方程组(A.3) ???? ?? ?=-+++=-+++=-+++=-+++79 .0sin 5.074.3cos 5.007.1cos sin 5.067.2cos 5.0y w v u t x w v u t y w v u t x w v u t (A.3) 以及关于z , t , u 的二维数表(见表A-1)确定了一个二元函数z =f (x , y )。 表A-1 二维数表 t z u 0 0.4 0.8 1.2 1.6 2 0 -0.5 -0.34 0.14 0.94 2.06 3.5 0.2 -0.42 -0.5 -0.26 0.3 1.18 2.38 0.4 -0.18 -0.5 -0.5 -0.18 0.46 1.42 0.6 0.22 -0.34 -0.58 -0.5 -0.1 0.62 0.8 0.78 -0.02 -0.5 -0.66 -0.5 -0.02 1.0 1.5 0.46 -0.26 -0.66 -0.74 -0.5 1. 试用数值方法求出f (x , y ) 在区域}5.15.0,8.00|), {≤≤≤≤=y x y x D (上的近似表达式 ∑∑===k i k j s r rs y x c y x p 00 ),( 要求p (x , y )以最小的k 值达到以下的精度 ∑∑==-≤-=10020 7210)],(),([i j i i i i y x p y x f σ 其中j y i x i i 05.05.0,08.0+==。 2. 计算),(),,(* ***j i j i y x p y x f (i =1,2,…,8 ; j =1,2,…,5) 的值,以观察p (x , y ) 逼 近f (x , y )的效果,其中j y i x j i 2.05.0,1.0**+==。

数值分析之幂法及反幂法C语言程序实例

数值分析之幂法及反幂法C 语言程序实例 1、算法设计方案: ①求1λ、501λ和s λ的值: s λ:s λ表示矩阵的按模最小特征值,为求得s λ直接对待求矩阵A 应用反幂法即可。 1λ、501λ:已知矩阵A 的特征值满足关系 1n λλ<< ,要求1λ、及501λ时,可 按如下方法求解: a . 对矩阵A 用幂法,求得按模最大的特征值1m λ。 b . 按平移量1m λ对矩阵A 进行原点平移得矩阵1m B A I λ=+,对矩阵B 用反幂法 求得B 的按模最小特征值2m λ。 c . 321m m m λλλ=- 则:113min(,)m m λλλ=,13max(,)n m m λλλ=即为所求。 ②求和A 的与数5011 140 k k λλμλ-=+最接近的特征值 ik λ(k=0,1,…39): 求矩阵A 的特征值中与k μ最接近的特征值的大小,采用原点平移的方法: 先求矩阵 B=A-k μI 对应的按模最小特征值k β,则k β+k μ即为矩阵A 与k μ最接近的特征值。 重复以上过程39次即可求得ik λ(k=0,1,…39)的值。 ③求A 的(谱范数)条件数2cond()A 和行列式det A : 在(1)中用反幂法求矩阵A 的按模最小特征值时,要用到Doolittle 分解方法,在Doolittle 分解完成后得到的两个矩阵分别为L 和U ,则A 的行列式可由U 阵求出,即:det(A)=det(U)。 求得det(A)不为0,因此A 为非奇异的实对称矩阵,则: max 2()s cond A λλ= ,max λ和s λ分别为模最大特征值与模最小特征值。

北航数值分析计算实习报告一

航空航天大学 《数值分析》计算实习报告 第一大题 学院:自动化科学与电气工程学院 专业:控制科学与工程 学生姓名: 学号: 教师: 电话: 完成日期: 2015年11月6日 航空航天大学 Beijing University of Aeronautics and Astronautics

实习题目: 第一题 设有501501?的实对称矩阵A , ??? ???? ?????????=5011A a b c b c c b c b a 其中,064.0,16.0),501,,2,1(64.0)2.0sin()024.064.1(1 .0-==???=--=c b i e i i a i i 。矩阵A 的特征值为)501,,2,1(???=i i λ,并且有 ||min ||,501 150121i i s λλλλλ≤≤=≤???≤≤ 1.求1λ,501λ和s λ的值。 2.求A 的与数40 1 5011λλλμ-+=k k 最接近的特征值)39,,2,1(???=k k i λ。 3.求A 的(谱数)条件数2)A (cond 和行列式detA 。 说明: 1.在所用的算法中,凡是要给出精度水平ε的,都取12-10=ε。 2.选择算法时,应使矩阵A 的所有零元素都不储存。 3.打印以下容: (1)全部源程序; (2)特征值),,39,...,2,1(,s 5011=k k i λλλλ以及A det ,)A (cond 2的值。 4.采用e 型输出实型数,并且至少显示12位有效数字。

一、算法设计方案 1、求1λ,501λ和s λ的值。 由于||min ||,501 150121i i s λλλλλ≤≤=≤???≤≤,可知绝对值最大特征值必为1λ和501 λ其中之一,故可用幂法求出绝对值最大的特征值λ,如果λ=0,则1λ=λ,否则 501λ=λ。将矩阵A 进行一下平移: I -A A'λ= (1) 对'A 用幂法求出其绝对值最大的特征值'λ,则A 的另一端点特征值1λ或501λ为'λ+λ。 s λ为按模最小特征值,||min ||501 1i i s λλ≤≤=,可对A 使用反幂法求得。 2、求A 的与数40 1 5011λλλμ-+=k k 最接近的特征值)39,...,2,1(=k k i λ。 计算1)1,2,...,50=(i i λ-k μ,其模值最小的值对应的特征值k λ与k μ最接近。因此对A 进行平移变换: )39,,2,1k -A A k k ==(I μ (2) 对k A 用反幂法求得其模最小的特征值'k λ,则k λ='k λ+k μ。 3、求A 的(谱数)条件数2)(A cond 和行列式detA 。 由矩阵A 为非奇异对称矩阵可得: | | )(min max 2λλ=A cond (3) 其中max λ为按模最大特征值,min λ为按模最小特征值,通过第一问我们求得的λ和s λ可以很容易求得A 的条件数。 在进行反幂法求解时,要对A 进行LU 分解得到。因L 为单位下三角阵,行 列式为1,U 为上三角阵,行列式为主对角线乘积,所以A 的行列式等于U 的行列式,为U 的主对角线的乘积。

北航数值分析第一次大作业(高斯gauss lu分解)

一、问题分析及算法描述 编写程序,分别用列主元的Gauss 消去法和LU 分解法求解下面线型代数方程组AX=b 的解,其中A 为N ×N 矩阵,N=50,其中第i(i ≥1)行、第j(i ≥1)列元素 a ij =1 i+j ?1, 右端向量b 的第i(i ≥1)个分量为 b i = 10 i+j ?1N j=1. 列主元素Gauss 消去过程中,要用到两种初等行变换。第一种,交换两行的位置;第二种,用一个数乘某一行加到另一行上。在第k 次消元之前,先对增广矩阵 A (k),b (k) 作第一种行变换,使得a ik (k) 中绝对值最大的元素交换到第k 行的主对角线位置上,然后再使用第二种行变换进行消元。如此往复,最后得到一个上三角系数矩阵,并回代求解解向量。由于每次消元前选取了列主元素,因此与顺序Guass 消元法相比,可提高数值计算的稳定性,且其计算量与顺序Guass 消元法相同。列主元的Gauss 消去法要求系数矩阵A 非奇异。 LU 分解法,即通过一系列初等行变换将系数矩阵A 分解成一个下三角矩阵L 与一个上三角矩阵U 的乘积,进一步通过求解两个三角矩阵得出解向量。若L 为单位下三角矩阵,U 是上三角矩阵,则称为Doolittle 分解;若L 为下三角矩阵,U 是单位上三角矩阵,则称为Crout 分解。若系数矩阵A 的前n-1阶顺序主子式不为零,则Doolittle\Crout 分解具有唯一性。若在每步行变换中选取主元,可提高数值计算稳定性。本算例中采用选主元的Doolittle 分解。 通过分析可知,本算例中待求解线型方程组系数矩阵为非奇异矩阵,且其前n-1阶顺序主子式不为零。方程组的解向量为x = 10,10,?,10 T 。满足列主元高斯消去法以及LU 分解法的基本使用条件。为了验证上述两种方法对本算例的适用性,笔者利用Microsoft Visual C++6.0编写了该算例的列主元高斯消去法以及LU 分解法的程序代码,并进行了运算求解。

北航数值分析大作业第二题精解

目标:使用带双步位移的QR 分解法求矩阵10*10[]ij A a =的全部特征值,并对其中的每一个实特征值求相应的特征向量。已知:sin(0.50.2)() 1.5cos( 1.2)(){i j i j ij i j i j a +≠+== (i,j=1,2, (10) 算法: 以上是程序运作的逻辑,其中具体的函数的算法,大部分都是数值分析课本上的逻辑,在这里特别写出矩阵A 的实特征值对应的一个特征向量的求法: ()[]()() []()[]()111111I 00000 i n n n B A I gause i n Q A I u Bu u λλ-?-?-=-?-?? ?-=????→=??????→= ?? ? 选主元的消元 检查知无重特征值 由于=0i A I λ- ,因此在经过选主元的高斯消元以后,i A I λ- 即B 的最后一行必然为零,左上方变 为n-1阶单位矩阵[]()()11I n n -?-,右上方变为n-1阶向量[]()11n Q ?-,然后令n u 1=-,则 ()1,2,,1j j u Q j n ==???-。

这样即求出所有A所有实特征值对应的一个特征向量。 #include #include #include #define N 10 #define E 1.0e-12 #define MAX 10000 //以下是符号函数 double sgn(double a) { double z; if(a>E) z=1; else z=-1; return z; } //以下是矩阵的拟三角分解 void nishangsanjiaodiv(double A[N][N]) { int i,j,k; int m=0; double d,c,h,t; double u[N],p[N],q[N],w[N]; for(i=0;i

数值分析幂法c语言实现

1.实验目的: 1熟练掌握C 语言程序设计,编程求解问题。 2.运用幂法求解住特征值和特征向量。 2.实验内容: 例题: 用幂法求 A= ??????????0.225.05.025.00.10.15.00.10.1 的特征值和特征向量。 完整代码以及截图如下: #include "stdio.h" #include "math.h" #define M 3 void main() { float fan(),max(),e1,e2,r1,r2; void au(),ex(),print_x(),std(); static float a[M][M]={{1.0,1.0,0.5},{1.0,1.0,0.25},{0.5,0.25,2.0}}; static float u0[M],u1[M],maxn0,maxn1; int i;

printf("*********************************\n"); printf("****** 幂法*********\n"); printf("******求特征值与特征向量*********\n"); printf("*********************************\n\n"); printf("input precision e1,e2:"); scanf("%f,%f",&e1,&e2); printf("\ninput u(%d):",M); for (i=0;ie1 || r2>e2) { printf("%4d",i++); print_x(u0); printf("\n"); ex(u0,u1); } else break; } while (1); } void au(a,u0,u1) float a[][M],u0[],u1[]; { int i,j; for (i=0;i

数值分析幂法与反幂法-matlab程序

数值分析幂法与反幂法 matlab程序 随机产生一对称矩阵,对不同的原点位移和初值(至少取3个)分别使用幂法求计算矩阵的主特征值及主特征向量,用反幂法求计算矩阵的按模最小特征值及特征向量。 要求 1)比较不同的原点位移和初值说明收敛性 2)给出迭代结果,生成DOC文件。 3)程序清单,生成M文件。 解答: >> A=rand(5) %随机产生5*5矩阵求随机矩阵 A = 0.7094 0.1626 0.5853 0.6991 0.1493 0.7547 0.1190 0.2238 0.8909 0.2575 0.2760 0.4984 0.7513 0.9593 0.8407 0.6797 0.9597 0.2551 0.5472 0.2543 0.6551 0.3404 0.5060 0.1386 0.8143 >> B=A+A' %A矩阵和A的转置相加,得到随机对称矩阵B B = 1.4187 0.9173 0.8613 1.3788 0.8044 0.9173 0.2380 0.7222 1.8506 0.5979 0.8613 0.7222 1.5025 1.2144 1.3467 1.3788 1.8506 1.2144 1.0944 0.3929 0.8044 0.5979 1.3467 0.3929 1.6286

B=?? ????? ???? ?? ???6286.13929.03467.15979.08044 .03929.00944 .12144.18506 .13788.13467.12144.15025.17222.08613.05979.08506.17222.02380.09173.08044.03788.18613 .09173 .04187.1 编写幂法、反幂法程序: function [m,u,index,k]=pow(A,u,ep,it_max) % 求矩阵最大特征值的幂法,其中 % A 为矩阵; % ep 为精度要求,缺省为1e-5; % it_max 为最大迭代次数,缺省为100; % m 为绝对值最大的特征值; % u 为对应最大特征值的特征向量; % index ,当index=1时,迭代成功,当index=0时,迭代失败 if nargin<4 it_max=100; end if nargin<3 ep=1e-5; end n=length(A); index=0; k=0; m1=0; m0=0.01; % 修改移位参数,原点移位法加速收敛,为0时,即为幂法 I=eye(n) T=A-m0*I while k<=it_max v=T*u; [vmax,i]=max(abs(v)); m=v(i); u=v/m; if abs(m-m1)

北航数值分析课程第一次大作业讲解

《数值分析A》计算实习题目第一题 一.算法设计方案: 1.矩阵A的存储与检索 将带状线性矩阵A[501][501]转存为一个矩阵MatrixC[5][501] . 由于C语言中数组角标都是从0开始的,所以在数组MatrixC[5][501]中检索A的带内元素a ij的方法是: A的带内元素a ij=C中的元素c i-j+2,j 2.求解λ1,λ501,λs ①首先分别使用幂法和反幂法迭代求出矩阵按摸最大和最小的特征值λmax和λmin。λmin即为λs; 如果λmax>0,则λ501=λmax;如果λmax<0,则λ1=λmax。 ②使用带原点平移的幂法(mifa()函数),令平移量p=λmax,求出对应的按摸最大的特征值λ,max, 如果λmax>0,则λ1=λ,max+p;如果λmax<0,则λ501=λ,max+p。 3.求解A的与数μk=λ1+k(λ501-λ1)/40的最接近的特征值λik (k=1,2,…,39)。 使用带原点平移的反幂法,令平移量p=μk,即可求出与μk最接近的特征值λik。 4.求解A的(谱范数)条件数cond(A)2和行列式d etA。 ①cond(A)2=|λ1/λn|,其中λ1和λn分别是矩阵A的模最大和最小特征值。 ②矩阵A的行列式可先对矩阵A进行LU分解后,detA等于U所有

对角线上元素的乘积。 二.源程序(VS2010环境下,C++语言) #include #include #include #include #include #include #include #define E 1.0e-12 /*定义全局变量相对误差限*/ int max2(int a,int b) /*求两个整型数最大值的子程序*/ { if(a>b) return a; else return b; } int min2(int a,int b) /*求两个整型数最小值的子程序*/ { if(a>b) return b; else return a; } int max3(int a,int b,int c) /*求三整型数最大值的子程序*/ { int t; if(a>b) t=a; else t=b; if(t

BUAA数值分析大作业三

北京航空航天大学2020届研究生 《数值分析》实验作业 第九题 院系:xx学院 学号: 姓名: 2020年11月

Q9:方程组A.4 一、 算法设计方案 (一)总体思路 1.题目要求∑∑=== k i k j s r rs y x c y x p 00 ),(对f(x, y) 进行拟合,可选用乘积型最小二乘拟合。 ),(i i y x 与),(i i y x f 的数表由方程组与表A-1得到。 2.),(* * j i y x f 与1使用相同方法求得,),(* * j i y x p 由计算得出的p(x,y)直接带入),(* * j i y x 求得。

1. ),(i i y x 与),(i i y x f 的数表的获得 对区域D ={ (x,y)|1≤x ≤1.24,1.0≤y ≤1.16}上的f (x , y )值可通过xi=1+0.008i ,yj=1+0.008j ,得到),(i i y x 共31×21组。将每组带入A4方程组,即可获得五个二元函数组,通过简单牛顿迭代法求解这五个二元数组可获得z1~z5有关x,y 的表达式。再将 ),(i i y x 分别带入z1~z5表达式即可获得f(x,y)值。 2.乘积型最小二乘曲面拟合 2.1使用乘积型最小二乘拟合,根据k 值不用,有基函数矩阵如下: ????? ??=k i i k x x x x B 0000 , ????? ??=k j j k y y y y G 0000 数表矩阵如下: ???? ? ? ?=),(),(),(),(0000j i i j y x f y x f y x f y x f U 记C=[rs c ],则系数rs c 的表达式矩阵为: 11-)(-=G G UG B B B C T T T )( 通过求解如下线性方程,即可得到系数矩阵C 。 UG B G G C B B T T T =)()( 2.2计算),(),,(* ***j i j i y x p y x f (i =1,2,…,31 ; j =1,2,…,21) 的值 ),(**j i y x f 的计算与),(j i y x f 相同。将),(**j i y x 代入原方程组,求解响应) ,(* *ij ij u t 进行分片双二次插值求得),(**j i y x f 。),(* *j i y x p 的计算则可以直接将),(**j i y x 代入所求p(x,y)。 二、 源程序 ********* 第三次数值分析大作业Q9************ integer::i, j, K1, L1, n, m dimension X(31), Y(21), T(6), U(6), Z(6, 6), UX(11, 21), TY(11, 21), FXY(11, 21), C(6, 6) dimension z1(31, 21), z2(31, 21), z3(31, 21), z4(31, 21), z5(31, 21) dimension X1(8), Y1(5), FXY1(8, 5), PXY1(8, 5), UX1(8, 5), TY1(8, 5)

北航数值分析第二次大作业--QR分解

《数值分析A》

一、算法设计方案 整个程序主要分为四个函数,主函数,拟上三角化函数,QR分解函数以及使用双步位移求解矩阵特征值、特征向量的函数。因为在最后一个函数中也存在QR分解,所以我没有采用参考书上把矩阵M进行的QR分解与矩阵Ak的迭代合并的方法,而是在该函数中调用了QR分解函数,这样增强了代码的复用性,减少了程序长度;但由于时间关系,对阵中方法的运算速度没有进行深入研究。 1.为了减少QR分解法应用时的迭代次数,首先对给定矩阵进行拟上三角化处理。 2.对经过拟上三角化处理的矩阵进行QR分解。 3.注意到计算特征值与特征向量的过程首先要应用前面两个函数,于是在拟上三角化矩阵的基础上对QR分解函数进行了调用。计算过程中,没有采用goto语句,而是根据流程图采用其他循环方式完成了设计,通过对迭代过程的合并,简化了程序的循环次数,最后在计算特征向量的时候采用了列主元高斯消去法。

二、源程序代码 #include #include #include int i,j,k,l,m; //定义外部变量double d,h,b,c,t,s; double A[10][10],AA[10][10],R[10][10],Q[10][10],RQ[10][10]; double X[10][10],Y[10][10],Qt[10][10],M[10][10]; double U[10],P[10],T[10],W[10],Re[10]={0},Im[10]={0}; double epsilon=1e-12; void main() { void Quasiuppertriangular(double A[][10]); void QRdecomposition(double A[][10]); void DoublestepsQR(double A[][10]); int i,j; for(i=0;i<10;i++) { for(j=0;j<10;j++) { A[i][j]=sin(0.5*(i+1)+0.2*(j+1)); Q[i][j]=0; AA[i][j]=A[i][j]; } A[i][i]=1.5*cos(2.2*(i+1)); AA[i][i]=A[i][i];

数值方法课程设计幂法反幂法计算矩阵特征值和特征向量附Matlab程序

数值方法课程设计幂法反幂法计算矩阵特征值和特征向量附Matlab程序

矩阵的特征值与特征向量的计算 摘要 物理,力学,工程技术中的很多问题在数学上都归结于求矩阵特征值的问题,例如振动问题(桥梁的振动,机械的振动,电磁振动等)、物理学中某些临界值的确定问题以及理论物理中的一些问题。矩阵特征值的计算在矩阵计算中是一个很重要的部分,本文使用幂法和反幂法分别求矩阵的按模最大,按模最小特征向量及对应的特征值。 幂法是一种计算矩阵主特征值的一种迭代法,它最大的优点是方法简单,对于稀疏矩阵比较合适,但有时收敛速度很慢。其基本思想是任取一个非零的初始向量。由所求矩阵构造一向量序列。再经过所构造的向量序列求出特征值和特征向量。 反幂法用来计算矩阵按模最小特征向量及其特征值,及计算对应于一个给定近似特征值的特征向量。本文中主要使用反幂法计算一个矩阵的按模最小特征向量及其对应的特征值。计算矩阵按模最小特征向量的基本思想是将其转化为求逆矩阵的按模最大特征向量。然后经过这个按模最大的特征向量反推出原矩阵的按模最小特征向量。

关键词:矩阵;特征值;特征向量;冥法;反冥法 THE CALCULATIONS OF EIGENVALUE AND EIGENVECTOR OF MATRIX ABSTRACT Physics, mechanics, engineering technology in a lot of problems in mathematics are attributed to matrix eigenvalue problem, such as vibration (vibration of the bridge, mechanical vibration, electromagnetic vibration, etc.) in physics, some critical values determine problems and

北航数值分析大作业第二次

《数值分析》计算实习作业 (第二题)

算法设计方案: 1、对矩阵A 赋值,取计算精度ε=1×10-12; 2、对矩阵A 进行拟上三角化,得到A (n-1),并输出A (n-1); 对矩阵A 的拟上三角化,通过直接调用子函数inftrianglize(A)来实现;拟上三角化得到的矩阵A (n-1)输出至文件solution.txt 中。 3、对A (n-1)进行QR 分解并输出Q 、R 及RQ 矩阵; QR 分解通过直接调用子函数QRdescom(A,Q,R, n)实现。 4、运用QR 方法求所有的特征值,并输出; (1)初始时令m=n ,在m>2的条件下执行; (2)判断如果|A mm-1|<ε,则得到一个特征值,m=m-1,转(4);否则转(3); (3)判断如果|A m-1m-2|<ε,则得到两个特征值,m=m-2,转(4); (4)判断如果m ≤2,转(6);否则转(5); (5)执行相似迭代,转(2); k k T k k k k k k k k k k Q A Q A R Q M I D A D tr A M ==+-=+1)2)det(( (6)求出最后的一个或两个特征值; (7)输出全部的特征值至文件solution.txt 中。 5、输出QR 分解法迭代结束之后的A (n-1)至文件solution.txt 中; 6、通过反幂法求出所有实特征值的特征向量并输出。 首先令B=(A-λi I),其中λi 是实特征值;反幂法通过调用子函数Bpowmethod(B,x1)实现,最终λi 对应的特征向量就是x1;最后将所有的实特征值的特征向量输出。

数值分析试验幂法与反幂法matlab

一、问题的描述及算法设计 (一)问题的描述 我所要做的课题是:对称矩阵的条件数的求解设计 1、求矩阵A 的二条件数 问题 A=?? ?? ? ?????----210121012 2、设计内容: 1)采用幂法求出A 的 错误!未找到引用源。. 2)采用反幂法求出A 的错误!未找到引用源。. 3)计算A 的条件数 ⅡA Ⅱ2* ⅡA -1Ⅱ2=cond2(A )=错误!未找到引用源。/错误!未找到引用源。.(精度要求为10-6) 3、设计要求 1)求出ⅡA Ⅱ2。 2)并进行一定的理论分析。 (二)算法设计 1、幂法算法 (1)取初始向量u )0((例如取u )0(=(1,1,…1)T ),置精度要求ε,置k=1. (2)计算 v )(k =Au )1(-k ,m k =max(v )(k ), u )(k = v )(k / m k (3)若| m k = m 1-k |<ε,则停止计算(m k 作为绝对值最大特征值1λ,u )(k 作为相应的特征向量)否则置k=k+1,转(2) 2、反幂法算法 (1)取初始向量u )0((例如取u )0(=(1,1,…1)T ),置精度要求ε,置k=1. (2)对A 作LU 分解,即A=LU (3)解线性方程组 Ly )(k =u )1(-k ,Uv )(k =y )(k (4)计算 m k =max(v )(k ), u )(k = v )(k / m k (5)若|m k =m 1-k |<ε,则停止计算(1/m k 作为绝对值最小特征值n λ,u )(k 作

为相应的特征向量);否则置k=k+1,转(3).

二、算法的流程图(一)幂法算法的流程图

北航数值分析计算实习报告一

北航数值分析计算实习 报告一 Company number:【0089WT-8898YT-W8CCB-BUUT-202108】

北京航空航天大学 《数值分析》计算实习报告 第一大题 学 院:自动化科学与电气工程学院 专 业: 控制科学与工程 学 生 姓 名: 学 号: 教 师: 电 话: 完 成 日 期: 2015年11月6日 北京航空航天大学 Beijing University of Aeronautics and Astronautics 实习题目: 第一题 设有501501?的实对称矩阵A , 其中,064.0,16.0),501,,2,1(64.0)2.0sin()024.064.1(1.0-==???=--=c b i e i i a i i 。矩阵A 的特征值为)501,,2,1(???=i i λ,并且有 1.求1λ,501λ和s λ的值。 2.求A 的与数40 1 5011λλλμ-+=k k 最接近的特征值)39,,2,1(???=k k i λ。 3.求A 的(谱范数)条件数2)A (cond 和行列式detA 。

说明: 1.在所用的算法中,凡是要给出精度水平ε的,都取12-10=ε。 2.选择算法时,应使矩阵A 的所有零元素都不储存。 3.打印以下内容: (1)全部源程序; (2)特征值),,39,...,2,1(,s 5011=k k i λλλλ以及A det ,)A (cond 2的值。 4.采用e 型输出实型数,并且至少显示12位有效数字。 一、算法设计方案 1、求1λ,501λ和s λ的值。 由于||min ||,501 150121i i s λλλλλ≤≤=≤???≤≤,可知绝对值最大特征值必为1λ和501λ其中之 一,故可用幂法求出绝对值最大的特征值λ,如果λ=0,则1λ=λ,否则501λ=λ。将矩阵A 进行一下平移: I -A A'λ= (1) 对'A 用幂法求出其绝对值最大的特征值'λ,则A 的另一端点特征值1λ或501λ为 'λ+λ。 s λ为按模最小特征值,||min ||501 1i i s λλ≤≤=,可对A 使用反幂法求得。 2、求A 的与数40 1 5011λλλμ-+=k k 最接近的特征值)39,...,2,1(=k k i λ。 计算1)1,2,...,50=(i i λ-k μ,其模值最小的值对应的特征值k λ与k μ最接近。因此对A 进行平移变换: )39,,2,1k -A A k k ==(I μ (2) 对k A 用反幂法求得其模最小的特征值'k λ,则k λ='k λ+k μ。

北航数值分析作业第一题题解

北航数值分析作业第一题: 一、算法设计方案 1.要求计算矩阵的最大最小特征值,通过幂法求得模最大的特征值,进行一定 判断即得所求结果; 2.求解与给定数值接近的特征值,可以该数做漂移量,新数组特征值倒数的绝 对值满足反幂法的要求,故通过反幂法即可求得; 3.反幂法计算时需要方程求解中间过渡向量,需设计Doolite分解求解; 4.|A|=|B||C|,故要求解矩阵的秩,只需将Doolite分解后的U矩阵的对角线相 乘即为矩阵的Det。 算法编译环境:vlsual c++6.0 需要编译函数:幂法,反幂法,Doolite分解及方程的求解 二、源程序如下: #include #include #include #include int Max(int value1,int value2); int Min(int value1,int value2); void Transform(double A[5][501]); double mifa(double A[5][501]); void daizhuangdoolite(double A[5][501],double x[501],double b[501]); double fanmifa(double A[5][501]); double Det(double A[5][501]); /***定义2个判断大小的函数,便于以后调用***/ int Max(int value1,int value2) { return((value1>value2)?value1:value2); } int Min(int value1,int value2) { return ((value1

关于幂法与反幂法的研究

东北大学秦皇岛分校 数值计算课程设计报告 幂法及反幂法 学院数学与统计学院 专业信息与计算科学 学号****** 姓名*** 指导教师*** *** 成绩 教师评语: 指导教师签字: 2014年07月07日

1 绪论 1.1 课题的背景 矩阵特征值的数值算法,在科学和工程技术中很多问题在数学上都归结为矩阵的特征值问题。例如,结构的振动波形和频率可分别由适当矩阵的特征向量和特征值来决定,结构的稳定性由特征值决定;又如机械和机件的振动问题,无线电工及光学系统第电磁振荡问题和物理学中各种临界值都牵涉到特征值计算。所以说研究利用数学软件解决求特征值的问题是非常必要的。 求矩阵特征值的一种方法是从原始矩阵出发,求出其特征多项式及其根,即得到矩阵的特征值。但高次多项式求根问题尚有困难,而且重根的计算精度较低。另外,原始矩阵求特征多项式系数的过程,对舍入误差非常敏感,对最终结果影响很大。所以,从数值计算的观点来看,这种求矩阵特征值的方法不够好。 实际问题中,有时需要的并不是所有的特征根,而是最大最小的实特征根。称模最大的特征根为主特征值。解决特征值计算的算法有很多种,古老的雅可比方法、兰乔斯方法以及较为常用的幂法、QR方法。QR方法是一种变换法,可求全部的特征值;幂法和反幂法是迭代法,只求模最大与模最小的特征值及特征向量。下面主要来研究一下幂法、反幂法,利用MATLAB解决矩阵特征值问题。 幂法是一种计算矩阵主特征值(矩阵按模最大的特征值)及对应特征向量的迭代方法,特别适用于大型稀疏矩阵。反幂法是计算海森伯格阵或三对角阵的对应一个给定近似特征值的特征向量的有效方法之一。 1.2 概念的认识 对于n阶矩阵A,若存在数λ和n维向量x满足:x =,则称λ为矩阵A的特征值, Axλ x为相应的特征向量。 病态矩阵:求解方程组时对数据的小扰动很敏感的矩阵。例如希尔伯特矩阵就是一类著名的病态矩阵。本次课题不对病态矩阵做深入研究。 非亏损矩阵:矩阵存在n个线性无关的特征向量,即有一个完全的特征向量组。 2 MATLAB特征值计算工具简介 查阅MATLAB HELP可以知道,利用eig函数可以快速求解矩阵的特征值和特征向量。可利用该函数对以下所做的幂法及反幂法程序进行检验。

相关文档
最新文档