数值分析matlab上机作业答案.doc-资料库

lemongrass1987-3100856-cc0ddd5c413b9140543efdf5d763c2d3.doc.pdf-第1页.png

第1页 / 共28页

lemongrass1987-3100856-cc0ddd5c413b9140543efdf5d763c2d3.doc.pdf-第2页.png

第2页 / 共28页

lemongrass1987-3100856-cc0ddd5c413b9140543efdf5d763c2d3.doc.pdf-第3页.png

第3页 / 共28页

lemongrass1987-3100856-cc0ddd5c413b9140543efdf5d763c2d3.doc.pdf-第4页.png

第4页 / 共28页

lemongrass1987-3100856-cc0ddd5c413b9140543efdf5d763c2d3.doc.pdf-第5页.png

第5页 / 共28页

lemongrass1987-3100856-cc0ddd5c413b9140543efdf5d763c2d3.doc.pdf-第6页.png

第6页 / 共28页

lemongrass1987-3100856-cc0ddd5c413b9140543efdf5d763c2d3.doc.pdf-第7页.png

第7页 / 共28页

lemongrass1987-3100856-cc0ddd5c413b9140543efdf5d763c2d3.doc.pdf-第8页.png

第8页 / 共28页

一、给定向量 x≠0，计算初等反射阵 Hk。 1．程序功能：给定向量 x≠0，计算初等反射阵 Hk。 2．基本原理：若      的分量不全为零，则由  ,  x ,  , x n ) 2 2 ) x ( x  1     u x      1   2   H I   ( sign x 1 e  1 2 2 uu u   u 2   (  x 1 )   I  uu 1 T T 2 2 确定的镜面反射阵 H 使得    ；当 (1   时，由 n k )  ( sign x )( k n  i k  x k ( 2 x i 1/ 2 ))   k , x k 1  ,  ,  ,0, u ( k ) ( k ) T ) u    k ( k  x k ( u    k        ) k (0,  1 (  k 2 H I    k ) x n ) = T  n R  k ( u ( k ) T )  k  u 1 ( k ) ( k ) ( u T ) 有 H x k  ( , x x 1 2 ,  , x  k k  , 1 ,0,  ,0) T  R n 算法： (1)输入 x，若 x 为零向量，则报错 (2)将 x 规范化，    如果 M＝0，则报错同时转出停机否则 x i  , iMx i  ,2,1 , n (3)计算 2x ，如果 1 ，则   (4)   ( 1 ) (5)计算 u x u ,  (1)  x  1  (6) H I   uu 1 T (7) y   ( M ,0,  ,0) (8)按要求输出,结束

3．变量说明： x －输入的 n 维向量； n － n 维向量 x 的维数； M － M 是向量 x 的无穷范数，即 x 中绝对值最大的一项的绝对值； p － Householder 初等变换阵的系数ρ； u － Householder 初等变换阵的向量 U s －向量 x 的二范数； x －输入的 n 维向量； n － n 维向量 x 的维数； p － Householder 初等变换阵的系数ρ； u － Householder 初等变换阵的向量 U k －数 k，H*x=y，使得 y 的第 k+1 项到最后项全为零； 4．程序代码： (1) function [p,u]=holder2(x) %HOLDER2 给定向量 x≠0,计算 Householder 初等变换阵的 p,u %程序功能：函数 holder2 给定向量 x≠0,计算 Householder 初等变换阵的 p,u； %输入：n 维向量 x； %输出：[p,u]。p 是 Householder 初等变换阵的系数ρ， u 是 Householder 初等变换阵的向量 U。 % n=length(x); % 得到 n 维向量 x 的维数； % 初始化 p,u； p=1;u=0; % 得到向量 x 的无穷范数，即 x 中绝对值最大的一项的绝对值； M=max(abs(x)); if M==0 % 如果 x=0，提示出错,程序终止; disp('Error: M=0'); return; else x=x/M; end; s=norm(x); if x(1)<0 s=-s; end u=x; u(1)=s+x(1); p=s*u(1); if n==1 u=0; end % 规范化 % 求 x 的二范数 % 首项为负，s 值要变号 % 除首项外，其余各项 x,u 相同 % 计算 u 的首项 % 计算 p % 若 x 是 1×1 维向量，则 u=0 (2) function H=holderk(x,k) %HOLDERK 给定向量 x≠0,数 k，计算初等反射阵 Hk，使 HkX=Y，其中 Y 的第 k+1 项到最后项全为零；

%程序功能：函数 holderk 给定向量 x≠0，数 k，计算初等反射阵 Hk，使 HkX=Y，%程序功能：函数 holder2 给定向量 x≠0,计算 Householder 初等变换阵的 p,u； %输入：n 维向量 x,数 k； %输出：H。H 是 Householder 初等变换阵，H*x=y，使得 y 的第 k+1 项到最后项全为零； %引用函数：holder2； n=length(x); if k>n % 得到 n 维向量 x 的维数； %如果 k 值溢出，报错； disp('Error: k>n'); % 初始化 H，并使 H(1:k,1:k)=I； % 得到计算 Householde 初等变换阵的系数ρ、向量 U； end H=eye(n); [p,u]=holder2(x(k:n)); H(k:n,k:n)=eye(n-k+1)-p\u*u'; % 计算 H(k:n,k:n)=I-p\u*u'； 5．使用示例: 情形 1: X 为零向量 >> x=[0,0,0,0]'; >> H=holderk(x,1) Error: M=0 H = 1 0 0 0 0 1 0 0 0 0 1 0 0 0 0 1 情形 2： K 值溢出： >> x=[1,2,3,4]'; >> H=holderk(x,5) Error: k>n 情形 3： K 值为 1： >> x=[2,3,4,5]'; >> H=holderk(x,1) H = -0.2722 -0.4082 -0.5443 -0.6804 -0.4082 0.8690 -0.1747 -0.2184 -0.5443 -0.1747 0.7671 -0.2911 -0.6804 -0.2184 -0.2911 0.6361 检验: >> det(H)

ans = -1.0000 >> H*x ans = -7.3485 0.0000 0.0000 0.0000 情形 4： (1)K 值为 3： >> x=[4,3,2,1]'; >> H=holderk(x,3) H = 1.0000 0 0 0 0 1.0000 0 0 0 0 -0.8944 -0.4472 0 0 -0.4472 0.8944 检验: >> det(H) ans = -1 >> H*x ans = 4.0000 3.0000 -2.2361 0 (2)K 值为 2: >> x=[4,3,2,1]'; >> H=holderk(x,2)

H = 1.0000 0 0 0 0 -0.8018 -0.5345 -0.2673 0 0 -0.5345 0.8414 -0.0793 -0.2673 -0.0793 0.9604 >> det(H) ans = -1.0000 >> H*x ans = 4.0000 -3.7417 0.0000 0 二、设 A 为 n 阶矩阵，编写用 Householder 变换法对矩阵 A 作正交分解的程序。 1．程序功能：给定 n 阶矩阵 A，通过本程序用 Householder 变换法对矩阵 A 作正交分解，得出 A＝QR 2．基本原理：任一实列满秩的 m×n 矩阵 A，可以分解成两个矩阵的乘积，即 A＝QR，其中 Q 是具有法正交列向量的 m×n 矩阵，R 是非奇异的 n 阶上三角阵。算法： (1)输入 n 阶矩阵 A (2)对  (3)计算上三角阵 R，仍然存储在 A AH A k n k ，求 Househoulder 初等反射阵的 )1  k k u , A  : , k k ( ) ( ) ( k 。 k  ,2,1  , n  1 j  , kk  ,1 , n t j  1  k )1      a n  ki  ) ( k ij ( k ij  a ut j (4)计算正交阵 Q Q  HQ I )1( k ) ( k k   au i ( k ij )    k  i )1 k  Q (

i ,2,1  , n t i  1  k , kk  ( )1 k   ij  kj  ,1  ( ) k q ij (5)按要求输出,结束 j q n ) ( k uq ij   k j   , n ut i  k j 3．变量说明： A －输入的 n 阶矩阵，同时用于存储上三角阵 R； n －矩（方）阵 A 的阶数； Q － Q 是具有法正交列向量的 n 阶矩阵； p,u －向量 A（k:n,k），对应初等反射阵的ρ,u k，jj，ii －循环变量； t1 －计算上三角阵 R 的系数 tj； t2 －计算正交矩阵 Q 的系数 ti； 4．程序代码： function [Q,A]=qrhh(A) %QRHH 用 Householder 变换法对 n 阶矩阵 A 作正交分解 A=QR； %程序功能：函数 qrhh 用 Householder 变换法对矩阵 A 作正交分解 A=QR； %输入：n 阶矩阵 A； %输出：[Q,A]。Q 是具有法正交列向量的 n 阶矩阵， % %引用函数： % [n,n]=size(A); Q=eye(n); for k=1:n-1 %求矩（方）阵 A 的阶数； %构造正交矩阵 Q(1)=I； holder2；示例 [p,u]=holder2(x); A（即 R）是非奇异的 n 阶上三角阵，仍用输入的矩阵 A 存储。 [p,u]=holder2(A(k:n,k));%向量 A（k:n,k），对应初等反射阵的ρ,u for jj=k:n %计算上三角阵 R（仍存贮于 A） t1=dot(u,A(k:n,jj))/p; %利用向量内积求和 A(k:n,jj)=A(k:n,jj)-t1*u; end for ii=1:n %计算正交矩阵 Q t2=dot(u,Q(ii,k:n))/p; %利用向量内积求和 Q(ii,k:n)=Q(ii,k:n)-t2*u'; end end 5．使用示例: (1)A 为 3 阶矩阵： >> A=[1 2 3; 2 3 0; 3 4 5]

A = 1 2 3 2 3 4 3 0 5 >> [q,r]=qrhh(A) -0.2673 -0.5345 -0.8018 0.8729 0.2182 -0.4364 0.4082 -0.8165 0.4082 q = r = -3.7417 0 -0.0000 -5.3452 0.6547 0.0000 -4.8107 0.4364 3.2660 检验： >> q*r ans = 1.0000 2.0000 3.0000 2.0000 3.0000 4.0000 3.0000 0.0000 5.0000 (2)A 为 4 阶矩阵： >> A=[1 2 3 4; 2 3 0 1; 3 4 5 6;1 6 8 0] A = 1 2 3 1 2 3 4 6 3 0 5 8 4 1 6 0 >> [q,r]=qrhh(A) q =

-0.2582 -0.5164 -0.7746 -0.2582 0.0597 -0.1045 -0.2688 0.9556 -0.2660 0.8434 -0.4662 -0.0222 -0.9268 -0.1049 0.3323 0.1399 r = -3.8730 0 0 0 -6.7132 4.4647 -0.0000 0.0000 -6.7132 -6.1968 6.4805 -3.3070 0 -1.4783 -3.0178 -1.8187 检验： >> q*r ans = 1.0000 2.0000 3.0000 1.0000 2.0000 3.0000 4.0000 6.0000 3.0000 -0.0000 5.0000 8.0000 4.0000 1.0000 6.0000 0 数值求解正方形域上的 Poisson 方程边值问题    2 u     2 x     (0, ) u y    2 u   2 y   (1, ) u y  ( , ) f x y  xy ,0  , x y  1  ( ,0) u x  ( ,1) 1 u x 

资料库

数值分析matlab上机作业答案.doc

相关推荐

课程资源

热门标签

最新资料