第1章 数值计算中的误差
1.1 实验目的
了解MATLAB基本操作、数值计算中误差的传递，避免误差的方法。
1.2 MATLAB基本操作
矩阵运算
%%行向量
>> x=[-3 1 7] %也可用逗号分隔
x =
-3 1 7
%%列向量
>> y=[0 2 9]' %行向量转置
y =
0
2
9
%%矩阵
>> A=[1 2 3;4 5 6;7 8 9]%分号换行
A =
1 2 3
4 5 6
7 8 9
%%矩阵转置
>> B=A'
B =
1 4 7
2 5 8
3 6 9
%%矩阵加减法运算
>> C=A+B
C =
2 6 10
6 10 14
10 14 18
>> C=A-1
C =
0 1 2
3 4 5
6 7 8
%%矩阵乘法运算
>> C=C*2
C =
0 2 4
6 8 10
12 14 16
>> C=A*C
C =
48 60 72
102 132 162
156 204 252
%%逆
>> A=[1 2;3 5];
>> inv(A)
ans =
-5.0000 2.0000
3.0000 -1.0000
%%矩阵除法运算
>> A=[1 2 3;4 5 6;7 8 0];
>> b=[1 1 1]';
>> x=A\b %左除
x =
-1.0000
1.0000
-0.0000
此处矩阵的除法运算相当于下面方程组的解：
A
x
=
b
Ax=b
Ax=b
即：
x
=
A
−
1
b
x=A^{-1}b
x=A−1b
数组运算
数组的运算符号是点运算，即 .* ./ 或 .\ 还有 .^ 。（加减运算无点运算）
左除和右除的区别在于交换了被除数和除数。
%%数组输入
>> a=0:2:10
a =
0 2 4 6 8 10
>> a=linspace(0,10,6) %构造一维等差数组
a =
0 2 4 6 8 10
%%数组加减法与矩阵加减法运算相同
%%数组乘法
>> a=0:2:10
a =
0 2 4 6 8 10
>> b=0:2:10
b =
0 2 4 6 8 10
>> c=a.*b
c =
0 4 16 36 64 100
%%数组除法
>> a./b
ans =
NaN 1 1 1 1 1
%%数组乘幂
>> c=a.^2
c =
0 4 16 36 64 100
逻辑运算
与：&，或：|，非：~。
>> a=[-1 2 4;5 4 -8];
>> b=[-2 1 -1;3 -1 2];
>> c=(a>0)&(b>0)
c =
2×3 logical array
0 1 0
1 0 0
>> d=(a>0)|(b>0)
d =
2×3 logical array
0 1 1
1 1 1
>> e=~(a>0)
e =
2×3 logical array
1 0 0
0 0 1
基本函数
基本初等函数
sin(x) %正弦
cos(x) %余弦
tan(x) %正切
asin(x) %反正弦
exp(x) %e的指数
log(x) %自然对数
sqrt(x) %平方根
abs(x) %绝对值
与矩阵有关的常用函数
二范数：二范数指矩阵A的2范数，就是A的转置共轭矩阵与矩阵A的积的最大特征根的平方根值，是指空间上两个向量矩阵的直线距离。类似于求棋盘上两点间的直线距离。
条件数：数值分析中，一个问题的条件数是该数量在数值计算中的容易程度的衡量，也就是该问题的适定性。一个低条件数的问题称为良态的，而高条件数的问题称为病态（或者说非良态）的。条件数是线性方程组Ax=b的解对b中的误差或不确定度的敏感性的度量。数学定义为矩阵A的条件数等于A的范数与A的逆的范数的乘积，即cond(A)=‖A‖·‖A的逆‖，对应矩阵的3种范数，相应地可以定义3种条件数。
>> A=[1 -2;-3 4];
>> norm(A) %求二范数
ans =
5.4650
>> cond(A) %求条件数
ans =
14.9330
>> rank(A) %求矩阵的秩
ans =
2
>> zeros(2,3) %生成一个2行3列的零矩阵
ans =
0 0 0
0 0 0
>> ones(2,3) %生成一个2行3列元素为1的矩阵
ans =
1 1 1
1 1 1
>> eye(3) %生成3阶单位阵
ans =
1 0 0
0 1 0
0 0 1
>> length(A) %求向量长度维数（列数）
ans =
2
控制流语句
for循环语句
for 循环变量 = 初值:增量:终值 %增量缺省值为1
语句
end
while循环语句
while 条件 %条件成立时执行下面语句，否则跳过
语句
end
if条件语句（三种形式）
if 条件
语句
end
if 条件
语句1
else
语句2
end
if 条件1
语句1
else if 条件2
语句2
else
语句3
end
M文件
命令文件
假设有一个tset.m文件如下：
n=6;
pp=1;
for i=1:n
pp=pp*i;
end
在命令输入窗口执行以下操作会出现相应结果：
>>test
>>pp
pp =
720
函数文件
function 因变量 = 函数(自变量)
%说明语句
语句
%计算n!的M文件
function tt = prod(n)
tt = 1;
for i = 1:n
tt = tt*i;
end
简单绘图
画函数曲线以及散点图
命令格式：plot(x,y,‘s’)
x代表横坐标，y代表纵坐标，s是可选择参数（默认为蓝色实线），具体意义如下：
例1：画出正余弦曲线
>> title('sin(x)&cos(x)') %图形标题
>> t=0:0.01:2*pi; %变量范围
>> plot(t,sin(t),'K',t,cos(t),'R')
>> xlabel('x');ylabel('y'); %坐标轴题注
结果如下：
例2：画出y=x^2 (0 <= x <= 1)的散点图
>> x=0:0.1:1;
>> y=x.^2;
>> plot(x,y,'.');
>> xlabel('x');ylabel('y');
结果如下：
例3：在区域[-pi,pi]×[-pi,pi]上画以下函数的网格曲面。
z
=
s
i
n
x
c
o
s
y
z=sinxcosy
z=sinxcosy
>> x=-pi:0.3:pi;y=-pi:0.3:pi;
>> [X,Y]=meshgrid(x,y);Z=sin(X).*cos(Y);
>> mesh(X,Y,Z);xlabel('x');ylabel('y');zlabel('z');
结果如下：
1.3 实验1例题：数值计算中的误差分析
例1：防止大数吃小数。设a=1×10308，b=1.1×10308，c=-1.001×10^308，用MATLAB计算a+(b+c)，(a+b)+c。
>> a=1e+308;b=1.1e+308;c=-1.001e+308;
>> r1=(a+b)+c
r1 =
Inf
>> r2=a+(b+c)
r2 =
1.0990e+308
在MATLAB中能表示的最大的正实数是realmax，即1.7977e+308。凡是比这个数大的数，MATLAB就用Inf表示（正无穷）。这个题中就是因为这个原因使得结合律不成立。
例2：分析下面程序，循环语句是否会结束？
a=2;
b=2;
while a+b~=a
b=b/2;
end
按理来说这个程序不会结束，因为b永远不可能为0。但是实际上这个程序很快就会结束，返回b的值为b=2.2204e-016，说明a远大于b的值时，b的值就会被忽略，这时a+b=a，程序结束。
例3：举例说明太小的数不能做除数。
计算以下表达式：
(
2
+
x
)
−
2
x
\frac{(2+x)-2}{x}
x(2+x)−2
>> x=3;y=((2+x)-2)/x
y =
1
>> x=1e-12;y=((2+x)-2)/x
y =
1.0001
>> x=1e-16;y=((2+x)-2)/x
y =
0
该表达式只要x不为零结果就一定等于1。当x越小的时候可以发现结果和1的偏差越大，所以太小的数不宜做除数。
分析下面程序，循环语句是否会结束？
a=2;
b=2;
while a+b~=a
b=b/2;
end
按理来说这个程序不会结束，因为b永远不可能为0。但是实际上这个程序很快就会结束，返回b的值为b=2.2204e-016，说明a远大于b的值时，b的值就会被忽略，这时a+b=a，程序结束。
例3：举例说明太小的数不能做除数。
计算以下表达式：
(
2
+
x
)
−
2
x
\frac{(2+x)-2}{x}
x(2+x)−2
>> x=3;y=((2+x)-2)/x
y =
1
>> x=1e-12;y=((2+x)-2)/x
y =
1.0001
>> x=1e-16;y=((2+x)-2)/x
y =
0
该表达式只要x不为零结果就一定等于1。当x越小的时候可以发现结果和1的偏差越大，所以太小的数不宜做除数。
349
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
