MATLAB求解方程之trapz函数

pulbieup

5 N/ `- k# c: S/ f! X3 J5 h* v( v" Qtrapz 是基于梯形法则的离散点积分函数。 调用形式：
I = trapz(x,y)
) k" _: O6 y4 ?& X2 Q8 h) Z, V其中 x 和 y 分别是自变量和对应函数值，以 sin(x) 在 [0,pi] 积分为例：
' e% j) f9 J- R1 T+ F5 q8 ^  a" Qx = linspace(0,pi,1e3);     %生成 [0,pi] 内的一系列离散点
- E8 m) c; d: sy = sin(x);
I = trapz(x,y)

浮点数误差
由于计算机中都是以二进制形式存储数据，当用十进制数进行计算时，就会存在十进制数二进制数的转换。但是某些十进制数转化为二进制数是一个无限位的小数，这时必须对该无限位小数进行截断计算机才能存储，那么此时就会产生误差，即浮点数误差。
例如十进制的0.9，在转化为二进制时是无限循环小数0.1110011001100110011...。此时须对该无限位小数进行截断才能保存在内存当中，截断后再转换回十进制时，0.9就变成了0.90000000000000002，这就是浮点数误差的产生过程。
由于浮点数误差的存在，当进行数值计算时就会出现一些不可避免的问题，最常见的就是判断两数相等时得到与预期相反的结果。
例：令 a = 0.1+0.2, b = 0.3, 判断 a==b 时，MATLAB 会返回0, 当执行 a-b 时，会发现结果不是精确等于0，而是一个非常小的数5.5511e-17。
或者在矩阵中寻找数的位置(也相当于是判断两数相等)。
  _. n5 V' I* e%例
; C$ \2 a9 l' s7 d, U( |  z>> a = 0.1:0.1:0.5

3 p7 K! N$ @& ~: Y8 @a =

    0.1000    0.2000    0.3000    0.4000    0.5000
# c8 F. s1 C: M% y
>> find(a==0.3)
- w: f6 b) D1 ?" @/ S7 V
ans =

   Empty matrix: 1-by-0
2 i& j9 Y. n  d, P+ b0 F* P) q由于 a 向量中的 0.3 是由 0.1+0.1+0.1 计算得到，在计算过程中就产生了浮点数误差，这也导致在判断 a==0.3 时结果为false，所以 find(a==0.3) 返回一个空矩阵。
% D3 G0 C* W# K! ~在进行数值计算判断两数相等时，最好不要直接判断，而是设立一个容差值，当两个浮点数的差的绝对值小于给定的容差值时，我们就认为这两个浮点数相等。
比如对于上面的例子：
>> a=0.1:0.1:0.5

- l3 T6 @- S4 ca =
% F- D4 p/ ?& D% B& B8 T
  0.1000   0.2000   0.3000   0.4000   0.5000

>> tol=eps(0.3)*10   %设立容差值，一般比这个点的浮点数误差高一到两个数量级即可。eps函数能够求得该点的浮点数误差值。
" p) |" b, s6 l/ A) D: t/ W
. X( O# V1 H2 ntol =
5 D* E9 Y9 L4 `  ^" [
2 V% H% p$ }& n& P+ ^  5.5511e-15
# G3 a* q3 }" K
! u6 r, G$ a- u# D. C: G1 o( H>> find(abs(a-0.3)
ans =
+ O5 ^2 Y3 J0 R" o3 s! y0 A& z' V
  3
" |+ K, _( v1 u+ L$ O

, X8 p9 U( p4 c6 X4 \$ O生成一系列有规律名变量
当循环迭代需要把每次迭代结果进行保存时，如果每次迭代的结果是尺寸不同的矩阵，无法用矩阵进行存储，那么可以利用 eval 和 num2str 这两个函数可以生成一系列例如 a1、a2、a3… 变量对结果进行保存(不推荐这种方法，原因是 eval 这个函数有很多缺点)。
4 T  y+ `. F% g1 @" S# ^eval:：将括号内的字符串视为语句并运行。
num2str： 将数值转换为字符串。
%例
for ii=1:10
    str=['a',num2str(ii),'=1:',num2str(ii)];
% d5 J) d/ D5 v; ^    eval(str)
3 O% Y4 P. w" W# d6 o0 H% fend
这样可以生成一系列变量 a1、a2…a10 对循环结果进行保存。
: r3 v$ F! X9 t$ L' F3 k不推荐使用 eval 函数的原因，帮助文档有详细的解释。
; M2 W/ Z! |* d3 yMATLAB® compiles code the first time you run it to enhance peRFormance for future runs. However, because code in an eval statement can change at run time, it is not compiled.
$ Y6 E! W2 }1 \+ |% |Code within an eval statement can unexpectedly create or assign to a variable already in the current workspace, overwriting existing data.
) _4 I+ m3 W8 t6 zConcatenating strings within an eval statement is often difficult to read. Other language constructs can simplify the syntax in your code.
! N9 o( L( }: k' OMATLAB 对于这类问题有更好的解决办法，利用元胞数组对结果进行存储。元胞数组是 MATLAB 中的特色数据类型，它的元素可以是任意类型的变量，包括不同尺寸或不同维度的矩阵。 对于上面的例子，利用元胞数组：
for ii=1:10
    a{ii}=1:ii;
end
; s8 m+ _" E' ]0 y即生成一系列元胞存储循环结果。这样无论是程序的可读性、运行效率还是后续程序对保存结果调用的方便程度，都远胜于 eval 函数。
除此之外，在处理符号变量时如果需要生成一系列有规律名符号变量，例如生成一个多项式：
1 U, Z. ]; y; a0 ]& H- N7 Hy = x1+2*x2+3*x3+…+100*x100
eval+num2str 能够实现，但更简便的方法还是利用矩阵：
x=sym('x',[1,100]);
5 k* X/ ?" \/ }' F# l2 u" J$ ?8 Ow=(1:100).*x;
y=sum(w)

& U5 @, s! }! x! X' P  I% W

youOK

MATLAB求解方程之trapz函数