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

合集下载

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

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

矩阵的特征值与特征向量的计算摘要物理,力学,工程技术中的很多问题在数学上都归结于求矩阵特征值的问题,例如振动问题(桥梁的振动,机械的振动,电磁振动等)、物理学中某些临界值的确定问题以及理论物理中的一些问题。

矩阵特征值的计算在矩阵计算中是一个很重要的部分,本文使用幂法和反幂法分别求矩阵的按模最大,按模最小特征向量及对应的特征值。

幂法是一种计算矩阵主特征值的一种迭代法,它最大的优点是方法简单,对于稀疏矩阵比较合适,但有时收敛速度很慢。

其基本思想是任取一个非零的初始向量。

由所求矩阵构造一向量序列。

再通过所构造的向量序列求出特征值和特征向量。

反幂法用来计算矩阵按模最小特征向量及其特征值,及计算对应于一个给定近似特征值的特征向量。

本文中主要使用反幂法计算一个矩阵的按模最小特征向量及其对应的特征值。

计算矩阵按模最小特征向量的基本思想是将其转化为求逆矩阵的按模最大特征向量。

然后通过这个按模最大的特征向量反推出原矩阵的按模最小特征向量。

关键词:矩阵;特征值;特征向量;冥法;反冥法THE CALCULATIONS OF EIGENVALUE AND EIGENVECTOR OF MATRIXABSTRACTPhysics, 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 theoretical physics in some of the problems. Matrix eigenvalue calculation is a very important part in matrix computation. In this paper, we use the power method and inverse power method to calculate the maximum of the matrix, according to the minimum characteristic vector and the corresponding characteristic value.Power method is an iterative method to calculate the eigenvalues of a matrix. It has the advantage that the method is simple and suitable for sparse matrices, but sometimes the convergence rate is very slow. The basic idea is to take a non - zero initial vector. Construct a vector sequence from the matrix of the matrix. Then the eigenvalues and eigenvectors are obtained by using the constructed vector sequence.The inverse power method is used to calculate the minimum feature vectors and their eigenvalues of the matrix, and to calculate the eigenvalues of the matrix. In this paper, we use the inverse power method to calculate the minimum eigenvalue of a matrix and its corresponding eigenvalues. The basic idea of calculating the minimum characteristic vector of a matrix is to transform it to the maximum characteristic vector of the modulus of the inverse matrix. Then, according to the model, the minimum feature vector of the original matrix is introduced.Key words: Matrix;Eigenvalue;Eigenvector;Iteration methods;目录1 引言 (1)2 相关定理。

数值分析幂法c语言实现

数值分析幂法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 3void 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;i<M;++i){scanf("%f",&u0[i]);}std(u0);maxn0=max(u0);i=0;printf("\n- - - - - - - - - - - - - - - - - -\n");printf(" ............NMD\n");do{au(a,u0,u1);maxn1=max(u1);std(u1);r1=fan(u0,u1);r2=(float)fabs(maxn0-maxn1);maxn0=maxn1;if (r1>e1 || r2>e2){printf("%4d",i++);print_x(u0);printf("\n");ex(u0,u1);}elsebreak;} while (1);}void au(a,u0,u1)float a[][M],u0[],u1[];{int i,j;for (i=0;i<M;i++){u1[i]=0;for (j=0;j<M;j++){u1[i]+=a[i][j]*u0[j];}}}void std(u)float u[];{int i;float t,max();t=max(u);for (i=0;i<M;i++){u[i]=u[i]/t;}}float fan(u0,u1)float u0[],u1[];{float max();int i;float uu[M];for (i=0;i<M;i++){uu[i]=u0[i]-u1[i];}return max(uu);}float max(u)float u[];{int i;float m;m=u[0];for (i=0;i<M;i++){if (u[i]>m){m=u[i];}}return m;}void ex(u0,u1)float u0[],u1[];{int i;for (i=0;i<M;i++){u0[i]=u1[i];}}void print_x(u)float u[];{int i;for (i=0;i<M;i++){printf("%12.6f",u[i]);}}3.运行结果:。

数值分析(10)幂法

数值分析(10)幂法

(2)1 2 , 1 3 , 且矩阵A有n个线性无关的特征向量。
) 由上式可知, y(k 是个摆动序列,当k 充分大时,有 V k k 1) 2 k 1 (1) (2) X x (2 ( x x ) X X 2 k 1 1 1 2 1 2
3 k (3) X X y k [1 x ( 1) x ( ) x V X 2 3 2 3 1
两点说明: 1)如果V y (0) 的选取恰恰使得1 0, 幂法计算仍能 0 进行。因为计算过程中舍入误差的影响,迭代若干
k) (1) X 次后,必然会产生一个向量V y (k , 它在x 1 方向上的分
量不为零,这样,以后的计算就满足所设条件。
k k) k (1) 2)因V y (k 1 x X 1 1 1 1 , 计算过程中可能会出现溢出
k 1 n n i k k i k 1 [1 X 1 ( ) i X i ] [1 X 1 ( ) i X i ] i 2 1 i 2 1 n n i k i k k 1 max [1 X 1 ( ) i X i ] max [1 X 1 ( ) i X i ] i 2 1 i 2 1 n
第四章 代数特征值问题
工程实践中有多种振动问题,如桥梁 或建筑物 的振动,机械机件、飞机机翼的振动,及 一些稳定 性分析和相关分析可转 化为求矩阵特征值与特征向 量的问题。
矩阵A aij

n n
的特征值是A的特征多项式
f ( ) I A 的n个零点.
但高次多项式求根精度低 , 一般不作为求解方 法. 目前的方法是针对矩阵不同的特点给出不同的 有效方法.

《幂法和反幂法》课件

《幂法和反幂法》课件

应用范围比较
总结词
幂法适用于求解特征值和特征向量,而反幂法适用于求解线性方程组和最小二 乘问题。
详细描述
幂法主要用于求解特征值和特征向量,在物理、工程和科学计算等领域有广泛 应用。反幂法适用于求解线性方程组和最小二乘问题,在统计学、机器学习和 数据分析等领域有广泛应用。
优缺点比较
总结词
幂法的优点在于能够求解特征值和特征向量,但缺点是计算复杂度高;反幂法的优点在于计算复杂度低,但缺点 是可能存在数值不稳定性。
幂法的性质
01
02
03
幂法具有高效性
相对于直接计算矩阵的幂 ,幂法可以大大减少计算 量和存储空间。
幂法具有收敛性
在适当的条件下,幂法能 够收敛到正确的矩阵幂的 结果。
幂法具有稳定性
在计算过程中,幂法能够 保持数值的稳定性,避免 误差的累积。
幂法的应用场景
数值分析
用于求解线性方程组、特 征值问题等数值计算问题 。
详细描述
幂法的优点在于能够精确求解特征值和特征向量,适用于需要高精度计算的情况。然而,由于其计算复杂度高, 对于大规模数据集可能效率较低。反幂法的优点在于计算复杂度相对较低,适用于处理大规模数据集。然而,反 幂法可能存在数值不稳定性,对于某些问题可能需要额外的数值稳定化技术。
04
幂法和反幂法的实现
05
幂法和反幂法的应用实 例
幂法在密码学中的应用
加密算法
幂法常被用于构造加密算法,如RSA算法。通过使用幂法,可以 快速地计算大数的幂次,从而实现高效的加密和解密过程。
密钥交换
在Diffie-Hellman密钥交换协议中,幂法被用于生成共享密钥,确 保通信双方安全地交换密钥。
数字签名

幂法,反幂法求解矩阵最大最小特征值及其对应的特征向量

幂法,反幂法求解矩阵最大最小特征值及其对应的特征向量

数值计算解矩阵的按模最大最小特征值及对应的特征向量一.幂法1. 幂法简介:当矩阵A 满足一定条件时,在工程中可用幂法计算其主特征值(按模最大)及其特征向量。

矩阵A 需要满足的条件为:(1) 的特征值为A i n λλλλ,0||...||||21≥≥≥>(2) 存在n 个线性无关的特征向量,设为n x x x ,...,,21 1.1计算过程:i ni i i u x x αα,1)0()0(∑==,有对任意向量不全为0,则有111111*********1111011)()(...u u a u a u λu λαu αA x A Ax x k n n k n k k ni i k i i n i i i k )(k (k))(k αλλλλλα++++=+=+++≈⎥⎦⎤⎢⎣⎡+++======∑∑ 可见,当||12λλ越小时,收敛越快;且当k 充分大时,有1)1111)11111λαλαλ=⇒⎪⎩⎪⎨⎧==+++(k )(k k (k k )(k x x u x u x ,对应的特征向量即是)(k x 1+。

2 算法实现.,, 3,,1 , ).5()5(,,,,||).4();max(,).3()(max(;0,1).2(,).1()()()(停机否则输出失败信息转置若转否则输出若计算最大迭代次数,误差限,初始向量输入矩阵βλβεβλβλε←+←<<-←←=←←k k N k y x Ay x x abs x y k Nx A k k k 3 matlab 程序代码function [t,y]=lpowerA,x0,eps,N) % t 为所求特征值,y 是对应特征向量k=1;z=0; % z 相当于λy=x0./max(abs(x0)); % 规范化初始向量x=A*y; % 迭代格式b=max(x); % b 相当于 βif abs(z-b)<eps % 判断第一次迭代后是否满足要求t=max(x);return;endwhile abs(z-b)>eps && k<Nk=k+1;z=b;y=x./max(abs(x));x=A*y;b=max(x);end[m,index]=max(abs(x)); % 这两步保证取出来的按模最大特征值t=x(index); % 是原值,而非其绝对值。

数值分析 -第7讲_幂法和反幂法

数值分析 -第7讲_幂法和反幂法
数值分析
则存在酉矩阵U使 定理9( Schur定理) 设A ∈ R n×n, r11 r12 L r1n r22 L r2n ∆ = R, U T AU = O rnn 其中rii (i = 1,2,L, n)为A的特征值.
定理10(实Schur分解) 设A ∈ R n×n, 则存在正交矩阵Q使 R11 R12 L R1m R22 L R2m , QT AQ = O Rmm 其中当Rii (i = 1,2,L, m)为一阶时Rii是A的实特征值,当Rii为 二阶时Rii的两个特征值是A的两个共轭复特征值.
xn xn
α1 x1 α1 x1
数值分析
不同范数选取下的特征值的计算
1. 取范数为2-范数时 取范数为2
T T yk −1uk = yk −1 Ayk −1 ⇒
α1 x1T α1 x1 A = λ1 α1 x1 2 α1 x1 2
对应的迭代公式
∀ u0 ∈ R n T η k −1 = uk −1uk −1 yk −1 = uk −1 η k −1 uk = Ayk −1 T β k = yk −1uk ( k = 1, 2,...)
数值分析
实际使用的迭代公式为: 实际使用的迭代公式为:
uk −1 yk −1 = u k −1 u = Ay k −1 k
于是可得
Auk −1 A2uk −2 A k u0 uk = = = L = k −1 uk −1 Auk −2 A u0
uk Ak u0 yk = = k uk A u0
数值分析
定义3 定义3 设A = (aij ) n×n , 令 n ( )i = ∑ | aij | (2) Di = {z | | z − aii |≤ ri , z ∈ C }, (i = 1,L, n) 1 r , j≠i 称Di为复平面上以aii为圆心以ri为半径的Gerschgorin圆盘.

幂法,反幂法求解矩阵最大最小特征值及其对应的特征向量(DOC)

幂法,反幂法求解矩阵最大最小特征值及其对应的特征向量(DOC)

数值计算解矩阵的按模最大最小特征值及对应的特征向量一.幂法1. 幂法简介:当矩阵A 满足一定条件时,在工程中可用幂法计算其主特征值(按模最大)及其特征向量。

矩阵A 需要满足的条件为: (1) 的特征值为A i n λλλλ,0||...||||21≥≥≥>(2) 存在n 个线性无关的特征向量,设为n x x x ,...,,21 1.1计算过程:i ni i i u xx αα,1)0()0(∑==,有对任意向量不全为0,则有1111112211211111111011)()(...u u a u a u λu λαu αA x A Ax x k n n k n k k ni ik i i ni i i k )(k (k))(k αλλλλλα++++=+=+++≈⎥⎦⎤⎢⎣⎡+++======∑∑ 可见,当||12λλ越小时,收敛越快;且当k 充分大时,有1)1111)11111λαλαλ=⇒⎪⎩⎪⎨⎧==+++(k )(k k(k k )(k x x u x u x ,对应的特征向量即是)(k x 1+。

2 算法实现.,, 3,,1 , ).5()5(,,,,||).4();max(,).3()(max(;0,1).2(,).1()()()(停机否则输出失败信息转置若转否则输出若计算最大迭代次数,误差限,初始向量输入矩阵βλβεβλβλε←+←<<-←←=←←k k N k y x Ay x x abs x y k N x A k k k3 matlab 程序代码function [t,y]=lpowerA,x0,eps,N) % t 为所求特征值,y是对应特征向量k=1;z=0; % z 相当于λy=x0./max(abs(x0)); % 规范化初始向量x=A*y; % 迭代格式b=max(x); % b 相当于βif abs(z-b)<eps % 判断第一次迭代后是否满足要求t=max(x);return;endwhile abs(z-b)>eps && k<Nk=k+1;z=b;y=x./max(abs(x));x=A*y;b=max(x);end[m,index]=max(abs(x)); % 这两步保证取出来的按模最大特征值t=x(index); % 是原值,而非其绝对值。

幂法 数值分析

幂法 数值分析
用反幂法计算出矩阵 Bk 的模最小的特征值。 Bk = A -μ k I
(其中,k = 1,2,…,39)。 6. 求出 A 的条件数 cond(A)2 以及 A 的行列式
detA,其中: cond(A)2 = |λ m|/|λ s|; detA = det(LU) = detU。
二、全部源程序:
//1.主程序 #include "stdafx.h" #include "iostream.h" #include "iomanip.h" #include "head.h" #include "math.h"
为:"<<u_min<<endl; cout<<endl;
//求出 A 最大的特征值
//创建初始向量 u0 for ( i = 0; i < n; ++i) u0[i] = 5; u0[0] = 1.0;
//生成矩阵 A 的非零元素 arr_ger( a_U, a_D, a_L, s, r, n, fabs(u_max));
int i, j; //生成上带宽矩阵 a_U 的非 0 元素 for ( i = 0; i < s; ++i ) {
a_U[i] = new double [n-1-i]; } for ( j = 0; j < n-1; ++j ) {
a_U[0][j] = 0.16; } for ( j = 0; j < n-2; ++j ) {
//将带状矩阵 A 压缩为矩阵 C arr_compress( c, a_D, a_U, a_L, s, r, n);

北航数值分析-lec7-幂法和反幂法

北航数值分析-lec7-幂法和反幂法
线性方程组求解
迭代收敛性
反幂法在求解特征值问题中的应用
特征值问题
反幂法主要用于求解矩阵的特征值和特征向量问题。通过迭代过程,反幂法能够找到矩阵的所有特征 值和对应的特征向量。
数值稳定性
反幂法在求解特征值问题时,需要关注数值稳定性问题。由于计算机浮点运算的误差累积,反幂法可 能会产生数值不稳定的解。因此,需要采取适当的策略来提高数值稳定性。
误差分析比较
幂法
由于幂法是通过连续的矩阵乘法来计算矩阵的幂,因此误差会随着计算的次数逐渐 累积。为了减小误差,需要选择合适的计算精度和减少计算次数。
反幂法
反幂法是通过求解线性方程组来计算矩阵的逆和行列式,因此误差主要来自于线性 方程组的求解精度。为了减小误差,需要选择合适的求解方法和提高求解精度。
202X
北航数值分析-lec7-幂法 和反幂法
单击此处添加副标题内容
汇报人姓名 汇报日期
目 录幂法介绍Fra bibliotek反幂法介绍
幂法和反幂法的比较
幂法和反幂法的实现细节
幂法和反幂法的实际应用案例
单击此处输入你的正文,文字是
您思想的提炼,请尽量言简意赅
的阐述观点
contents
单击此处输入你的正文,文字是 您思想的提炼,请尽量言简意赅 的阐述观点
反幂法的实现细节
反幂法是一种迭代算法,用 于求解线性方程组的近似逆。
反幂法的收敛速度取决于矩阵的谱 半径,如果矩阵的谱半径较小,则 反幂法收敛速度较快。
ABCD
反幂法的实现步骤包括:选择初始 矩阵、计算迭代矩阵、更新解矩阵 和判断收敛性。
在实际应用中,反幂法通常用于 求解大规模稀疏线性系统的预处 理和后处理问题。
01

数值幂法及反幂法分析方法

数值幂法及反幂法分析方法

即为A的特征值(近似),vk为相应的特征向量.
考虑主特征值1,由于1
vk1 i
vk
i
其中vk i为vk的
第i个分量.
又由 vk1 i
vk
i
1k
的收敛速度,
比值r 2 确定.r越小,收敛速度越快. 1
这种由非零向量v0及矩阵A的乘幂Ak构造向量
序列vk 以计算A的主特征值i及相应特征向量
的方法叫幂法.
,
y
PX
yk xk
1
P
Pi,
j
1
cos
sin
1
(i行)
sin
cos (j行)
1
X x1, x2 ,, xi ,, x j ,, xn T
Y y1, y2 ,, yi , y j ,, yn T .
P的性质:
1。P为正交阵,PPT I.
2。P为单位阵I只在i,i, j, j, i, j, j,i四个位置元素
第9章矩阵的特征值与特
引言 幂法及反幂法 Jacobi方法 Householder方法 QR算法
§1 引 言
预备知识:
定理1 定理2
如果i i 1,2,, n是矩阵A的特征值,则有
n
n
1。 i aii trA;
i 1
i 1
2。det A 1 2n ;
设A与B为相似矩阵(即存在非奇异矩阵
xi
r
vk1 1k ai xi k
i 1
vk1 1vk ,
vk 1 i vk i
1.
注意:应用幂法进行上机计算时,一般将迭代
向量 规范化:
v0 vk
u0 0
Auk1 k

24讲:91幂法反幂法

24讲:91幂法反幂法

0.998077
11 {0.999012, 0.999977, 0.999989}
0.999034
12 {0.999508, 0.999992, 0.999996}
0.999515
{{1, 2, 3}, {{1, 1, 1}, {1, 0, 0}, {1, -1, 1}}}
模最小特征值 1,特征向量约为{0.9995,0.9999,0.9999}
ai
( i 1
)k 1 (vi
)
j
幂方法
解题步骤即k:充 1i1分1 x大xi( ki( k时 )1lki)m, ,((X(X(iXX(k((kk(1)1k)),)12)))j),jj.j..n11 (1)任给n维初始向量X (0) 0
(2)按X (k) AX (k1) Ak X (0)计算X (k)
Print[k," ",x,"
",y[[1]]/x[[1]]];
y=x/Max[Abs[x]],{k,1,20}]
Eigensystem[A]
反幂法程序运行结果
运行结果为:
1 {0.166667, 0.333333, 0.666667}
0
2 {0.458333, 0.666667, 0.833333}
1

P 1 AP
2

...



n
的充要条件是A具有n个线性无关的特征向量.
反之,如果A Rnn有m个(m n)不同的特征值
1 , 2 ,...m, 则对应的特征向量x1, x2 ,...xm线性无关.
二、幂法运算及程序
设n阶方阵A, 任取初始向量X (0) ,进行迭代计算: X (k1=) AX (k )

c语言幂次方

c语言幂次方

C语言幂次方一、引言幂次方是数学中非常常见的概念,表示一个数的多次乘法运算。

在计算机科学中,C语言提供了一种计算幂次方的方法,即使用指数运算符”**“。

本文将介绍C语言中计算幂次方的方法和相关注意事项。

二、计算幂次方的方法1. 使用指数运算符在C语言中,可以使用指数运算符”**“来计算幂次方。

该运算符与其他基本运算符(如加法、减法、乘法、除法等)相似,但在两个数之间使用它时,左侧的数将作为底数,右侧的数将作为指数,得到的结果即为幂次方的值。

2. 编写自定义函数除了使用指数运算符外,我们还可以编写自定义函数来计算幂次方。

这种方法可以扩展到更复杂的幂次方计算,例如计算浮点数的幂次方、计算负指数幂等。

下面是一个简单的示例代码:#include <stdio.h>double power(double base, int exponent) {double result = 1.0;for (int i = 0; i < exponent; i++) {result *= base;}return result;}int main() {double base;int exponent;printf("请输入底数和指数:");scanf("%lf %d", &base, &exponent);double result = power(base, exponent);printf("%lf 的 %d 次幂为 %lf\n", base, exponent, result);return 0;}三、注意事项1. 底数和指数的类型在计算幂次方时,要注意底数和指数的类型。

如果底数为整数而指数为负数,则可以得到小数的结果,但如果底数为浮点数,则结果可能会不准确。

这是因为浮点数的精度问题导致计算结果的误差。

因此,在计算浮点数的幂次方时,应该进行适当的取整或使用更精确的计算方法。

数值分析幂法和反幂法

数值分析幂法和反幂法
>> [m,u,index,k]=pow(A,u)
m =
2.6813
u =
0.8576
0.6934
0.5623
1.0000
index =
1
k =
49
修改M0=1e-3
m =
2.6814
u =
0.8576
0.6934
0.5623
1.0000
index =
0
k =
1001
修改M0=0 %此时为幂法
m =
4.对结果进行解释说明;
采用方法
及结果
说明
对于幂法和反幂法求解矩阵特征值和特征向量的问题将从问题分析,算法设计和流程图,理论依据,程序及结果进行阐述该问题。
一.问题的分析:
求n阶方阵A的特征值和特征向量,是实际计算中常常碰到的问题,如:机械、结构或电磁振动中的固有值问题等。对于n阶矩阵A,若存在数 和n维向量x满足
end
n=length(A);
index=0;
k=0;
m1=0;
m0=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)<ep;
index=1;
break;
end
m=m+m0;
1
16
0
0.3847
[-0.8996 1.0000 0.2726 -0.2364]
1
16
[1 3 5 7]
0.5
0.3847
[-0.8995 1.0000 0.2726 -0.2364]

C语言程序设计100例之(41):快速幂运算

C语言程序设计100例之(41):快速幂运算

C语⾔程序设计100例之(41):快速幂运算例41 快速幂运算题⽬描述输⼊三个整数 b,p,k(0≤b,p,k<231),求 b^p mod k输⼊格式⼀⾏三个整数 b,p,k输出格式输出 b^p mod k=s (s 为运算结果)输⼊样例2 10 9输出样例2^10 mod 9=7(1)编程思路。

在实际应⽤中,我们经常会⽤到幂运算,例如,a n为a的n次幂。

求a的n次⽅通常采⽤快速幂运算。

下⾯我们来探讨快速幂运算的思路。

由于乘法具有结合律,因此 a4 = a*a * a *a = (a*a) * (a*a) = a2 * a2。

由此可以得到这样的结论:当n为偶数时,a n = a n/2 * a n/2;当n为奇数时,a n = a n/2 * a n/2 * a (其中n/2取整)。

这样,我们可以采⽤⼀种类似于⼆分的思想快速求得a的n次幂。

例如,a9 =a*a*a*a*a*a*a*a*a (⼀个⼀个乘,要乘9次)=a*(a*a)*(a*a)*(a*a)*(a*a)=a*(a2)4= a*((a2)2)2(A平⽅后,再平⽅,再平⽅,再乘上剩下的⼀个a,只乘4次)(2)源程序。

#include <stdio.h>typedef long long LL;LL quickPower(LL a, LL b,LL k){LL p = 1;while (b){if (b&1) p = (p*a)% k;b >>= 1;a = (a%k*a%k)%k;}return p%k;}int main(){LL a,b,k;scanf("%lld%lld%lld",&a,&b,&k);printf("%lld^%lld mod %lld=%lld\n",a,b,k,quickPower(a,b,k));return 0;}习题4141-1 组合数题⽬描述⼩明刚学了组合数。

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

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

一、问题分析及算法描述1. 问题的提出:(1)用幂法、反幂法求矩阵A =[a ij ]20×20的按摸最大和最小特征值,并求出相应的特征向量。

其中 a ij ={sin (0.5i +0.2j ) i ≠j 1.5cos (i +1.2j ) i =j要求:迭代精度达到10−12。

(2)用带双步位移的QR 法求上述的全部特征值,并求出每一个实特征值相应的特征向量。

2. 算法的描述:(1) 幂法幂法主要用于计算矩阵的按摸为最大的特征值和相应的特征向量。

其迭代格式为:{ 任取非零向量u 0=(h 1(0),⋯,h n (0))T|h r (k−1)|=max 1≤j≤n |h r (k−1)| y ⃑ k−1=u ⃑ k−1|h r (k−1)| u ⃑ k =Ay ⃑ k−1=(h 1(k ),⋯,h n (k ))T βk =sgn (h r (k−1))h r (k ) (k =1,2,⋯) 终止迭代的控制选用≤ε。

幂法的使用条件为n ×n 实矩阵A 具有n 个线性无关的特征向量x 1,x 2,⋯,x n ,其相应的特征值λ1,λ2,⋯,λn 满足不等式|λ1|>|λ2|≥|λ3|≥⋯≥|λn |或λ1=λ2=⋯=λm|λ1|>|λm+1|≥|λm+2|≥⋯≥|λn |幂法收敛速度与比值|λ2λ1|或|λm+1λ1|有关,比值越小,收敛速度越快。

(2) 反幂法反幂法用于计算n ×n 实矩阵A 按摸最小的特征值,其迭代格式为:{任取非零向量u 0∈R nηk−1=√u ⃑ k−1T u ⃑ k−1 y ⃑ k−1=u ⃑ k−1ηk−1⁄ Au ⃑ k =y ⃑ k−1 βk =y ⃑ k−1u ⃑ k (k =1,2,⋯) 每迭代一次都要求解一次线性方程组Au ⃑ k =y ⃑ k−1。

当k 足够大时,λn ≈1βk ,y ⃑ k−1可近似的作为矩阵A 的属于λn 的特征向量。

数值分析3.1幂法和反幂法

数值分析3.1幂法和反幂法

第三章 矩阵的特征值与特征向量
3.1 幂法与反幂法 3.2 Jacobi方法
3.3 QR方法
第三章 矩阵的特征值与特征向量
3.1幂法与反幂法
一、乘幂法 二、反幂法
三、带原点位移的反幂法
四、反幂法的特点
第三章 矩阵的特征值与特征向量
3.1幂法与反幂法
一、乘幂法
1、基本思想
2、算法(迭代公式) ◆一般算法
具体算法: (1)使用范数 2
1 X 1 yk , k 1 1 X 1
(2)使用范数
uk A yk 1

k
er u k er y k 1
T
T
k
lim k 1
留为作业自学
具体算法: (1)使用范数 2 1 X 1 yk , k 1 1 X 1
1 2 n
第三章 矩阵的特征值与特征向量
一、乘幂法 1、基本思想 设A有n个线性无关的特征向量 X 1 , X 2 ,, X n ,
AX j j X j , j 1,2,, n
3.1幂法与反幂法
★ 设 1为实数而且是单根: 1 2 n
u0 1 X 1 2 X 2 n X n
具体算法: 按取范数的不同, 迭代公式也不同。 (1)使用范数 2
任取初始向量u0 R n T k 1 u k 1 u k 1 u k 1 yk 1 k 1 (3.4) u k A yk 1 k yk 1T uk k 1,2,
T
精确结果:
X 1 (0,0.5,1) , 1 45
T
max( uk ) 表示 u k 的绝对值最大的分量。 (3)

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

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

数值分析幂法与反幂法matlab程序随机产生一对称矩阵,对不同的原点位移和初值(至少取3个)分别使用幂法求计算矩阵的主特征值及主特征向量,用反幂法求计算矩阵的按模最小特征值及特征向量。

要求1)比较不同的原点位移和初值说明收敛性2)给出迭代结果,生成DOC文件。

3)程序清单,生成M文件。

解答:>> A=rand(5) %随机产生5*5矩阵求随机矩阵A =0.7094 0.1626 0.5853 0.6991 0.14930.7547 0.1190 0.2238 0.8909 0.25750.2760 0.4984 0.7513 0.9593 0.84070.6797 0.9597 0.2551 0.5472 0.25430.6551 0.3404 0.5060 0.1386 0.8143>> B=A+A' %A矩阵和A的转置相加,得到随机对称矩阵BB =1.4187 0.9173 0.8613 1.3788 0.80440.9173 0.2380 0.7222 1.8506 0.59790.8613 0.7222 1.5025 1.2144 1.34671.3788 1.8506 1.2144 1.0944 0.39290.8044 0.5979 1.3467 0.3929 1.6286B=⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡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; endn=length(A);index=0;k=0;m1=0;m0=0.01;% 修改移位参数,原点移位法加速收敛,为0时,即为幂法I=eye(n)T=A-m0*Iwhile k<=it_maxv=T*u;[vmax,i]=max(abs(v));m=v(i);u=v/m;if abs(m-m1)<ep;index=1;break;endm=m+m0;m1=m;k=k+1;endfunction[m,u,index,k]=pow_inv(A,u,ep,it_max)% 求矩阵最大特征值的反幂法,其中% A为矩阵;% ep为精度要求,缺省为1e-5;% it_max为最大迭代次数,缺省为100;% m为绝对值最大的特征值;% u为对应最大特征值的特征向量;% index,当index=1时,迭代成功,当index=0时,迭代失败if nargin<4it_max=100;endif nargin<3ep=1e-5;endn=length(A);index=0;k=0;m1=0;m0=0;% 修改移位参数,原点移位法加速收敛,为0时,即为反幂法I=eye(n);T=A-m0*I;invT=inv(T);while k<=it_maxv=invT*u;[vmax,i]=max(abs(v));m=v(i);u=v/m;if abs(m-m1)<epindex=1;break;endm1=m;k=k+1;endm=1/m;m=m+m0;修改输入的m0的值,所得结果:幂法:反幂法:THANKS !!!致力为企业和个人提供合同协议,策划案计划书,学习课件等等打造全网一站式需求欢迎您的下载,资料仅供参考。

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

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

数值分析之幂法及反幂法C 语言程序实例1、算法设计方案:①求1λ、501λ和s λ的值:s λ:s λ表示矩阵的按模最小特征值,为求得s λ直接对待求矩阵A 应用反幂法即可。

1λ、501λ:已知矩阵A 的特征值满足关系 1n λλ<<L ,要求1λ、及501λ时,可按如下方法求解:a . 对矩阵A 用幂法,求得按模最大的特征值1m λ。

b . 按平移量1m λ对矩阵A 进行原点平移得矩阵1m BA I λ=+,对矩阵B 用反幂法求得B 的按模最小特征值2m λ。

c . 321m m m λλλ=-则:113min(,)m m λλλ=,13max(,)n m m λλλ=即为所求。

②求和A 的与数5011140k 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()scond A λλ=,max λ和s λ分别为模最大特征值与模最小特征值。

2、程序源代码:#include<>#include<>#include<>#define N 501 3e\n",k,k,value_s);}}void main(){float cond;double value_det;printf("Contact me\n");Init_matrix_A(); 3e\n",value_1);printf("λ501=%.13e\n",value_N);value_det=Det_matrix(); 3e\n",value_s);cond=Get_cond_A(); 3e\n",cond);printf("value_det=%.13e\n",value_det); }3、程序运行结果:4、迭代初始向量的选取对计算结果的影响:本次计算实习求矩阵A的具有某些特征的特征值,主要用到的方法是幂法和反幂法,这两种方法从原理上看都是迭代法,因此迭代初始向量的选择对计算结果会产生一定影响,主要表现在收敛速度上。

反幂法程序

反幂法程序

反幂法:#include <stdio.h>#include <math.h>#include <stdlib.h>void Gauss(int h,double a[10][10],doubleb[10],double t[10]);main(){int k,i,j,n,m;doubleq=0,l,y[10],vec[10],c[10][10],u[10],pat=0,a[10][10];printf("\n 请输入矩阵阶数(不要超过10):\n ");scanf("%d",&n);m=n;printf("\n请输入矩阵:\n ");for(i=1;i<m+1;i++){for(j=1;j<m+1;j++){scanf("%lf",&c[i-1][j-1]);}}for(i=1;i<m+1;i++){for(j=1;j<m+1;j++){printf("%f ",c[i-1][j-1]);}printf("\n");}for(i=0;i<n;i++)u[i]=1;for(k=0;k<500;k++) //循环求解特征向量及特征值{for(i=0;i<n;i++){q+=u[i]*u[i];}l=sqrt(q);q=0;for(i=0;i<n;i++) //求特征向量{y[i]=u[i]/l;}Gauss(n,c,y,u); //调用求解线性方程组矩阵,解得向量Upat=0;for(i=0;i<n;i++) //求特征值{pat+=y[i]*u[i];}}for(i=0;i<n;i++)vec[i]=y[i];printf("按模最小的特征值的特征向量为:\n");for(i=0;i<n;i++){printf("%f\n",y[i]);}pat=1/pat;printf("按模最小的特征值为:\n");printf("%f\n",pat);system("pause");//跳出主函数,结束执行}void Gauss(int h,double a[10][10],double b[10],double t[10]) //实现列主元高斯消去法函数{int k,i,j,l,n=0;double g,d[10],c[10][10];n=h;for(i=0;i<n;i++){for(j=0;j<n;j++)c[i][j]=a[i][j];d[i]=b[i];}for(k=0;k<n-1;k++) //列主元消去过程{g=0;for(i=k;i<n;i++) // 找列主对角元以下最大值if(fabs(c[i][k])>g){g=fabs(c[i][k]);l=i; // l为列主对角元以下最大值所在行}for(i=k;i<n;i++) // 交换行,最大值放在主对角位置上{g=c[k][i];c[k][i]=c[l][i];c[l][i]=g;}g=d[k]; //交换d[k]与d[l]d[k]=d[l];d[l]=g;for(i=k+1;i<n;i++) //消去过程{g=c[i][k]/c[k][k];for(j=k;j<n;j++)c[i][j]=c[i][j]-g*c[k][j];d[i]=d[i]-g*d[k];}}t[n-1]=d[n-1]/c[n-1][n-1];for(k=n-2;k>=0;k--) //回带求解{g=0;for(j=k+1;j<n;j++)g=g+c[k][j]*t[j];t[k]=(d[k]-g)/c[k][k];} }。

  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

数值分析之幂法及反幂法C 语言程序实例1、算法设计方案:①求1λ、501λ和s λ的值:s λ:s λ表示矩阵的按模最小特征值,为求得s λ直接对待求矩阵A 应用反幂法即可。

1λ、501λ:已知矩阵A 的特征值满足关系 1n λλ<<,要求1λ、及501λ时,可按如下方法求解: a . 对矩阵A 用幂法,求得按模最大的特征值1m λ。

b . 按平移量1m λ对矩阵A 进行原点平移得矩阵1m BA I λ=+,对矩阵B 用反幂法求得B 的按模最小特征值2m λ。

c . 321m m m λλλ=-则:113min(,)m m λλλ=,13max(,)n m m λλλ=即为所求。

②求和A 的与数5011140k 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()scond A λλ=,max λ和s λ分别为模最大特征值与模最小特征值。

2、程序源代码:#include<stdio.h>#include<stdio.h>#include<math.h>#define N 501 //列#define M 5 //行#define R 2 //下带宽#define S 2 //上带宽#define K 39#define e 1.0e-12 //误差限float A[M][N]; //初始矩阵float u[N]; //初始向量float y[N],yy[N];float maximum,value1,value2,value_1,value_N,value_s,value_abs_max;const float b=0.16f,c=-0.064f;int max_sign,max_position;void Init_matrix_A() //初始化矩阵A{int i;for(i=2;i<N;i++){A[0][i]= c;}for(i=1;i<N;i++){A[1][i]= b;}for(i=0;i<N;i++){A[2][i]= (1.64-0.024*(i+1))*sin(0.2*(i+1))-0.64*exp(0.1/(i+1));}for(i=0;i<N-1;i++){A[3][i]= b;}for(i=0;i<N-2;i++){A[4][i]= c;}}void Init_u() //初始化迭代向量{int i;for(i=0;i<N;i++)u[i]=1.0;}void Get_max() //获得绝对值最大的数值的模{int i;max_position=0;maximum=fabs(u[0]);for(i=1;i<N;i++){if(maximum<fabs(u[i])){max_position=i;maximum=fabs(u[i]);}}if(u[max_position]<0)max_sign=-1;else max_sign=1;}void Get_y() //单位化迭代向量{int i;for(i=0;i<N;i++)y[i]=u[i]/maximum;}void Get_u() //获得新迭代向量{int i;u[0]=A[2][0]*y[0]+A[1][1]*y[1]+A[0][2]*y[2];u[1]=A[3][0]*y[0]+A[2][1]*y[1]+A[1][2]*y[2]+A[0][3]*y[3];u[N-2]=A[4][N-4]*y[N-4]+A[3][N-3]*y[N-3]+A[2][N-2]*y[N-2]+A[1][N-1]*y[N-1];u[N-1]=A[4][N-3]*y[N-3]+A[3][N-2]*y[N-2]+A[2][N-1]*y[N-1];for(i=2;i<N-2;i++)u[i]=A[4][i-2]*y[i-2]+A[3][i-1]*y[i-1]+A[2][i]*y[i]+A[1][i+1]*y[i+1]+A[0][i+2]*y[i+2]; }void Get_value() //获得迭代后特征值{value2=value1;value1=max_sign*u[max_position];}void Check_value() //幂法第二迭代格迭代{Init_u();Get_max();Get_y();Get_u();Get_value();while(1){Get_max();Get_y();Get_u();Get_value();if(fabs((value2-value1)/value1)<e)break;}}void The_value() //获取绝对值最大的特征值λ_501 {Check_value();value_abs_max=value1;}void The_Other_value() //获取特征值λ_1{int i;float value_temp=value1;for(i=0;i<N;i++){A[2][i]-=value_temp;}Check_value();value1+=value_temp;if(value1<value_temp){value_1=value1;value_N=value_temp;}else{value_N=value1;value_1=value_temp;}}int min(int a,int b) //两值中取最小{if(a<b)return a;elsereturn b;}int max(int a,int b) //两值中取最大{if(a<b)return b;elsereturn a;}void Resolve_LU(){int k,i,j,t;float temp;for(k=1;k<=N;k++){for(j=k;j<=min(k+S,N);j++){temp=0;for(t=max(max(1,k-R),j-S);t<=k-1;t++)temp+=A[k-t+S][t-1]*A[t-j+S][j-1];A[k-j+S][j-1]=A[k-j+S][j-1]-temp;}for(i=k+1;i<=min(k+R,N);i++){temp=0;for(t=max(max(1,i-R),k-S);t<=k-1;t++)temp+=A[i-t+S][t-1]*A[t-k+S][k-1];A[i-k+S][k-1]=(A[i-k+S][k-1]-temp)/A[S][k-1];}}}void Back_substitution()//方程组回代过程{int i,t;float temp=0;for(i=2;i<N+1;i++){for(t=max(1,i-R);t<i;t++)y[i-1]-=A[i-t+S][t-1]*y[t-1];}u[N-1]=y[N-1]/A[S][N-1];for(i=N-1;i>0;i--){temp=0;for(t=i+1;t<=min(i+S,N);t++)temp+=A[i-t+S][t-1]*u[t-1];u[i-1]=(y[i-1]-temp)/A[S][i-1];}}double Det_matrix() //求矩阵行列式值{int i;double det=1;Init_matrix_A();Resolve_LU();for(i=0;i<N;i++)det=det*A[2][i];return det;}float Get_norm() //获得迭代向量模{int i;float normal=0;for(i=0;i<N;i++)normal+=u[i]*u[i];normal=sqrt(normal);return normal;}void Get_yy(float normal) //迭代向量单位化{int i;for(i=0;i<N;i++){y[i]=u[i]/normal;yy[i]=y[i];}}void Get_value_s() //获得绝对值最小的特征值{int i;value2=value1;value1=0;for(i=0;i<N;i++)value1+=yy[i]*u[i];value1=1/value1;}void Value_min() //反幂法求绝对值最小的特征值{float norm=0;int count=0;value1=0,value2=0;Init_u();norm=Get_norm();Get_yy(norm);Back_substitution();Get_value_s();while(count<10000){count++;norm=Get_norm();Get_yy(norm);Back_substitution();Get_value_s();if(fabs((value2-value1)/value1)<e)break;}value_s=value1;}float Get_cond_A() //求矩阵条件数{float cond1;cond1=fabs(value_abs_max/value_s);return cond1;}void Value_translation_min() //偏移条件下反幂法求特征值{int i,k;float tr;for(k=1;k<K+1;k++){tr=value_1+k*(value_N-value_1)/40;Init_matrix_A();for(i=0;i<N;i++)A[2][i]-=tr;Resolve_LU();Value_min();value_s+=tr;printf("k=%d =>>>λi%d=%.13e\n",k,k,value_s);}}void main(){float cond;double value_det;printf("Contactme:****************\n");Init_matrix_A(); //初始化矩阵AThe_value(); //获取绝对值最大的特征值λ_501 The_Other_value(); //获取特征值λ_1printf("λ1=%.13e\n",value_1);printf("λ501=%.13e\n",value_N);value_det=Det_matrix(); //求矩阵行列式值Value_min(); //反幂法求绝对值最小的特征值printf("λs=%.13e\n",value_s);cond=Get_cond_A(); //求矩阵条件数Value_translation_min();//偏移条件下反幂法求特征值printf("cond_A=%.13e\n",cond);printf("value_det=%.13e\n",value_det);}3、程序运行结果:4、迭代初始向量的选取对计算结果的影响:本次计算实习求矩阵A的具有某些特征的特征值,主要用到的方法是幂法和反幂法,这两种方法从原理上看都是迭代法,因此迭代初始向量的选择对计算结果会产生一定影响,主要表现在收敛速度上。

相关文档
最新文档