当前位置:首页 > 10级数值方法计算实习题
?10?7二、已知Wilson矩阵A???8??7解x??1111?。
T787??32??23?565??,且向量b???,则方程组Ax?b有准确?33?6109????5910??31?⑴用Matlab 内部函数求A,A的所有特征值和cond(A)2;
78.17.2??10?7.085.04?65?,解方程组(A??A)(x??x)?b,并求出向量⑵令A??A???85.989.899???99.98??6.994.99?x和?x2,从理论结果和实际计算结果两方面分析方程组Ax?b解的相对误差
?x2x2与A的相对误差?A2A2的关系;
⑶再改变扰动矩阵?A(其元素的绝对值不超过0.005),重复第2问。
解:(1)由Matlab内部函数:
>> A=[10 7 8 7;7 5 6 5;8 6 10 9;7 5 9 10];b=[32;23;33;31]; >> D=det(A) D = 1
>> X=eig(A) X =
0.0102 0.8431 3.8581 30.2887 >> C=cond(A) C =
2.9841e+003
(2) >> B=[10 7 8.1 7.2;7.08 5.04 6 5;8 5.98 9.89 9;6.99 4.99 9 9.98];
>> b=[32;23;33;31]; 求?x 假设X= ?x >> [rank(B),rank([B,b])] >> x=B\\b;x1=[1;1;1;1];X=x-x1 ans = X = 4 4 -6.4761 >> x=B\\b 10.4934 x =-5.4761 -2.4292 11.4934 1.4838 -1.4292
求 ?x2
2.4838
5
>>norm(X) ans =
12.6552 12.655就是 %求
?x2。
?x2x2
>> norm(X)/norm(x1) ans =
6.3276
6.3276 即为
?x2x2
inv(A)
>> d=inv(A);
>> norm(d)*norm(a)*norm(x) ans =
288.4858
>> 1-norm(d*(a)) ans =
-12.3693
>> 288.4858/ -12.3693 ans =
-23.3227
>> norm(a) ans =
0.2244 >> norm(A) ans =
30.2887
>> norm(a)/norm(A) ans =
0.0074
所以
||A?1||||?A||||x||所以=-23.322
1?||A?1(?A)||>>norm(X) ans =
0.0074
||A?1||||?A||||x||所以?x2> ?11?||A(?A)||?A2A2=0.0074
(3)改变?A,取a2=[0 0.002 0.001 0.003;0.001 0.004 0 0.001;0.003 -0.004 -0.001 0;-0.001 -0.002 0 -0.003] B2=a2+A; C2=[B b] b=[32;23;33;31]
>> rank(B2) ans = 4
>> rank(C2) ans = 4
rank(B2)= rank(C2) 所以扰动矩阵有唯一解 >> x2=B2\\b x2 =
1.0649 0.8893 1.0272
0.9859
>> x=B\\b;x1=[1;1;1;1];X=x-x1 (设X= ?x)
6
X2 =
-6.4761 10.4934 -2.4292 1.4838 norm(X2) >> norm(X2) ans =
12.6552 12.6552就是
?x2
>> norm(X)/norm(x1) >> norm(X2)/norm(x1) ans =
6.3276
所以
?x2x2=6.3276
>> norm(a2) ans =
0.0071
>> norm(a)/norm(A) >> norm(a2)/norm(A)
ans =
2.3336e-004
所以
?A2A2=2.3336e-004
norm(d)*norm(a2)*norm(x2) ans =
9.0875
> 1-norm(d*(a2)) ans =
0.6943 norm(X2) ans =
12.6552
9.0875/0.6943 ans =
13.0887
所以||A?1||||?A||||x||1?||A?1(?A)||= 13.0887
?x<||A?1||||?A||||x||21?||A?1(?A)|| 7
三、解三对角线性方程组的追赶法及其应用
⑴编写解三对角线性方程组的追赶法的通用程序,并应用于方程组
?2?1000??x1??1???12?100??x??0?T???2???52111???0?12?10??x3???0?,检验程序的正确性;(解为x?) ,,,,?????????63236?x00?12?10???4?????000?12????0???x5????d2udu?u?ex?3sinx,0?x??1?2?h?⑵求微分方程边值问题?dx的数值解(取步长),dx128?u(0)??2,u(?)?e??3?并与精确解比较(精确解为u?1)。 1?xu?2ui?ui?1duui?1?ui?1d2u,|?,i?1,2,?,n?1 说明:离散化微分方程时,2|x?xi?i?1x?xdxh2dxi2h解:(1)用C++语言编写追赶法求方程解,源代码如下:
#include
using namespace std; int main() {
int k,n=5; double
a[5]={0,-1,-1,-1,-1},b[5]={2.0,2.0,2.0,2.0,2.0},c[5]={-1,-1,-1,-1},r[5]={1,0,0,0,0},u[5],v[5],x[5];
u[0]=r[0]/b[0]; v[0]=c[0]/b[0]; for(k=1;k v[k]=c[k]/(b[k]-v[k-1]*a[k]); for(k=1;k<=n-1;k++) u[k]=(r[k]-u[k-1]*a[k])/(b[k]-v[k-1]*a[k]); x[n-1]=u[n-1]; for(k=n-2;k>=0;k--) x[k]=u[k]-v[k]*x[k+1]; for(int i=0;i<=n-1;i++) cout< 8
共分享92篇相关文档