现在时间是:
当前位置:首 页 >> 矩阵计算>> 教学区>> 文章列表

130213037+左芝+矩阵计算上机作业九

作者:130213037左芝   发布时间:2015-06-23 16:58:49   浏览次数:130

 

3
程序:
function [v,u]=house(x)
n=length(x);
v=zeros(n,1);
a=norm(x,inf);
x=x/a;
b=x(2:n)'*x(2:n);
v(1)=1;
v(2:n)=x(2:n);
if b==0
    u=0;
else
    c=sqrt(x(1)*x(1)+b);
    if x(1)<=0
        v(1)=x(1)-c;
    else
        v(1)=-b/(x(1)+c);
    end
    u=(2*v(1)*v(1))/(b+v(1)*v(1));
    v=v/v(1);
end
验证:
clc
x=[1 0 4 6 3 4]';
x1=[0;3;4];
A=eye(6);
n=size(x1,1);
[v,u]=house(x1)
H=eye(n)-u*v*v'
A([2,5,6],[2,5,6])=H;
y=A*x
ans=y(2)
运行结果:
 
v =
 
    1.0000
   -0.6000
   -0.8000
 
 
u =
 
     1
 
 
H =
 
         0    0.6000    0.8000
    0.6000    0.6400   -0.4800
    0.8000   -0.4800    0.3600
 
 
y =
 
    1.0000
    5.0000
    4.0000
    6.0000
         0
   -0.0000
 
 
ans =
 
     5
 
6、
程序1
%given
function [G]=givens(x)
 
n=length(x);
 
for i=1:n
 
    for j=1:n
 
        if i==j
 
            G(i,j)=1;
 
        else G(i,j)=0;
 
    end
 
    end
 
end
 
        for i=n-1:-1:1
 
            if x(i+1,1)==0
 
                continue;
 
            else if abs(x(i+1,1))>abs(x(i,1))
 
        t=x(i,1)/x(i+1,1);
 
        s=1/sqrt(1+t*t);
 
        c=s*t;
 
          
 
    else
 
        t=x(i+1,1)/x(i,1);
 
        c=1/sqrt(1+t*t);
 
        s=c*t;
 
    end
 
            end    
 
    x(i,1)=c*x(i,1)+s*x(i+1,1);
 
    x(i+1,1)=0;
 
    for j=1:n
 
        a=G(i,j);b=G(i+1,j);
 
        G(i,j)=c*a+s*b;
 
        G(i+1,j)=c*b-s*a;
 
    end  
 
        end
 
end
 
 
程序2
function [cs]=givens2(x,y)
 
n=length(x);
 
for i=1;n
 
    %cs*x=(y,a)'
 
a=sqrt(x(i,1)*x(i,1)+x(i+1,1)*x(i+1,1)-y(i,1)*y(i,1));
 
if a==0
 
    cs(i,1)=1;
 
    cs(i+1,1)=0;
 
else
 
    cs(i,1)=(x(i,1)*y(i,1)+a*x(i+1,1))/(norm(x,2)*norm(x,2));
 
    cs(i+1,1)=(x(i+1,1)*y(i,1)-a*x(i,1))/(norm(x,2)*norm(x,2));
 
    x(i,1)=y(i,1);
 
    x(i,1)=a;
 
end
 
end
 
end
 
 
测试:
%function [G]=givens3(x,y)
 
n=6;
 
G=zeros(n,n);
 
A=rand(n,1);
 
B=rand(n,1);
 
%x,y为单位向量
 
x=A/norm(A);
 
%E=1
 
E=norm(x,2);
 
y=B/norm(B);
 
[G]=givens(x);
 
n=length(x);
 
%验证x规范,即G*x=1,0。。。)’
 
e=G*x;
 
t=0;
 
for i=1:n-1
 
    [cs]=givens2(x,y);
 
    for j=1:n
 
     c=G(i,j);
 
     s=G(i+1,j);
 
     G(i,j)=cs(1,1)*c+cs(2,1)*s;
 
     G(i+1,j)=cs(1,1)*s-cs(2,1)*c;
 
     t=t+y(i,1)*y(i,1);
 
    end
 
end
 
     G
 
     d=G*x;
 
     err=abs(y-d)
 
     if t==1
 
         break;
 
     end
运行结果:
n =
     6

n =
     6

n =
     6

n =
     6

n =
     6

G =
    0.3165    0.0124    0.0326    0.0261    0.0479    0.0410
   -0.0183    0.2303    0.0168    0.0135    0.0247    0.0211
   -0.0041   -0.0538    0.2104   -0.0276   -0.0506   -0.0432
   -0.0009   -0.0119   -0.0492    0.2155   -0.0433   -0.0370
   -0.0002   -0.0026   -0.0109   -0.0504    0.1649   -0.1497
   -0.0002   -0.0024   -0.0100   -0.0465   -0.1359    0.1987

err =
    0.0766
    0.0401
    0.5569
    0.3052
    0.4254
    0.5476
 






上一篇:没有了    下一篇:没有了

Copyright ©2019    计算数学达人 All Right Reserved.

技术支持:自助建站 | 领地网站建设 |短信接口 |燕窝 版权所有 © 2005-2019 lingw.net.粤ICP备16125321号 -5