matlab实现C语言功能
3.6:
一.求实对称矩阵特征值与特征向量的雅可比过关法
1.1定义的雅可比过关法function函数
%本函数采用Jacobi方法计算实对称矩阵的所有特征值和对应特征向量
%返回值D为特征值对角阵,V为对应特征向量构成的正交方阵
%即有V’AV=D,V’*V=I
%采用查找绝对值最大的非对角元素方法
function [D,V]=jacobi(A)
tic;
%检验输入是否合法
b=size(A);
if b(1)~=b(2) %行列不等
error('MATLAB:Jaco:Invalid Matrix,The Matrix input should be a Symmetry Phalanx. See Jaco.');
end
n=max(b);
for i=1:n %非对称
for j=1:n
if abs(A(i,j)-A(j,i))>eps %不能用不等号,因为机器有误差
error('MATLAB:Jaco:Invalid Phalanx,The Phalanx input should be a Symmetric one. See Jaco.');
end
end
end
%实际计算
%初始化D,V为单位矩阵,避免多次分配空间
%并且相乘不会造成影响
D=eye(n);
V=eye(n);
%采用扫描绝对值最大A(p,q)的算法
p=0; %储存最大元素所在行
q=0; %储存最大元素所在列
maxpq=0; %储存绝对值最大元素
for i=1:n-1
for j=i+1:n
if abs(A(i,j))>abs(maxpq)
maxpq=A(i,j);
p=i;
q=j;
end
end
end
while abs(maxpq)>eps
maxpq=0;%务必清零,否则会死循环
phi=atan2(2*A(p,q),A(p,p)-A(q,q))/2;
U=eye(n);
U(p,p)=cos(phi);
U(q,q)=cos(phi);
U(p,q)=-sin(phi);
U(q,p)=sin(phi);
D=U'*A*U;
V=V*U;
A=D;
%改写maxpq
for i=1:n-1
for j=i+1:n
if abs(A(i,j))>abs(maxpq)
maxpq=A(i,j);
p=i;
q=j;
end
end
end
end
toc;
end
1.2运行的题目代码:
A=[10 1 2 3 4;1 9 -1 2 -3;2 -1 7 3 -5;3 2 3 12 -1;4 -3 -5 -1 15];
[D,V]=jacobi(A);
%[V,D]=eig(A,'nobalance')
1.3 实验结果图片:
二.幂乘法求实矩阵最大的特征值与相应特征向量
3.7
2.1定义的幂乘法function函数
function [a,x,n] = pmethod(A,x0,maxit,tol)
if nargin == 3
tol = 1.0e-6;
elseif nargin == 2
maxit = 1000;
tol = 1.0e-6;
end
a0 = 0;
y = A * x0;
a = max(abs(y));
x = y / a;
n = 1;
while (abs(a - a0) > tol) && (n < maxit)
a0 = a;
x0 = x;
y = A * x0;
a = max(abs(y));
x = y / a;
n = n + 1;
end
if (n > maxit)
disp(['迭代步数大于',num2str(maxit),'收敛太慢了!']);
else
disp(['幂法迭代',num2str(n),'步收敛']);
end
function [x,y]=A
[x,y]=eig(A)
eigenvalue=diag(y);%求对角线向量
lamda=max(eigenvalue);%求最大特征值
for i=1:length(A)%求最大特征值对应的序数
if lamda==eigenvalue(i)
break;
end
end
[m,n]=size(x);%得到行数和列数
y_lamda=x(:,i);%求矩阵最大特征值对应的特征向量
%标准化
sum = 0;
for i=1:m
sum = sum + x(i,1);
end
for i=1:m
y_lamda(i,1)= x(i,1)/sum;
end
disp('最大特征值')
lamda
disp('标准化的特征值向量')
y_lamda
2.2运行的题目代码:
maxit = 2000;
A=[0 1 1.5;-5 -0.5 1;-1 2 3.5];
x0=[1;1;1];
tol = 1.0e-7;
[a,x,n] = pmethod(A,x0,maxit,tol);
A=[0 1 1.5;-5 -0.5 1;-1 2 3.5];%待求矩阵
[x,y]=eig(A);%求矩阵的特征值和特征向量,x为特征向量矩阵,y为特征值矩阵。
eigenvalue=diag(y);%求对角线向量
lamda=max(eigenvalue);%求最大特征值
for i=1:length(A)%求最大特征值对应的序数
if lamda==eigenvalue(i)
break;
end
end
[m,n]=size(x);%得到行数和列数
y_lamda=x(:,i);%求矩阵最大特征值对应的特征向量
%标准化
sum = 0;
for i=1:m
sum = sum + x(i,1);
end
for i=1:m
y_lamda(i,1)= x(i,1)/sum;
end
disp('最大特征值')
lamda
disp('标准化的特征值向量')
y_lamda
2.3实验结果图片: