EDA365电子论坛网

标题: MATLAB求解方程之trapz函数 [打印本页]

作者: pulbieup    时间: 2021-2-19 18:02
标题: MATLAB求解方程之trapz函数

! P# I+ I9 q8 z' ]trapz 是基于梯形法则的离散点积分函数。 调用形式:
8 a7 e/ }5 E/ |, D( o( [I = trapz(x,y)7 C( u, ?, O" `. l
其中 x 和 y 分别是自变量和对应函数值,以 sin(x) 在 [0,pi] 积分为例:
, u5 _" K4 T/ j. Z: ?8 l8 j- ?x = linspace(0,pi,1e3);     %生成 [0,pi] 内的一系列离散点6 [. y+ `" |3 r& V0 q; F
y = sin(x);
  h* O, c& m# n1 [1 P# |. \I = trapz(x,y)$ P' s: Q/ E; ^- @6 g

* b. T! V3 _  b9 I浮点数误差
+ e, e+ L0 K% m  |由于计算机中都是以二进制形式存储数据,当用十进制数进行计算时,就会存在十进制数二进制数的转换。但是某些十进制数转化为二进制数是一个无限位的小数,这时必须对该无限位小数进行截断计算机才能存储,那么此时就会产生误差,即浮点数误差。
# E1 I0 ]' |2 `; V( k6 U9 ^! t例如十进制的0.9,在转化为二进制时是无限循环小数0.1110011001100110011...。此时须对该无限位小数进行截断才能保存在内存当中,截断后再转换回十进制时,0.9就变成了0.90000000000000002,这就是浮点数误差的产生过程。 ; e1 t' ^+ U: ]& n- e
由于浮点数误差的存在,当进行数值计算时就会出现一些不可避免的问题,最常见的就是判断两数相等时得到与预期相反的结果。 8 _9 Z/ s# x  X$ {3 w  d  I
例:令 a = 0.1+0.2, b = 0.3, 判断 a==b 时,MATLAB 会返回0, 当执行 a-b 时,会发现结果不是精确等于0,而是一个非常小的数5.5511e-17。
& L8 k( z5 @# Z' ^1 H  P8 I; p2 ~或者在矩阵中寻找数的位置(也相当于是判断两数相等)。
5 q3 J( d0 b& v1 u$ F( @%例; g8 B2 T6 g3 v9 T/ [8 Q
>> a = 0.1:0.1:0.52 i% D" z1 F  k9 S/ n$ g7 d
7 v" b, h7 v7 }* B6 b
a =. m/ {& k$ S; y' o3 j4 k/ p

( N. p  M! h) A- N  `) n! h( C$ B    0.1000    0.2000    0.3000    0.4000    0.5000& ?% ]; t; Q* g2 I/ [& g
. U0 [) {8 o: H% f
>> find(a==0.3)1 r* ?. w7 `, g' c- ~; q5 \* x
4 W5 f6 M: a$ g0 X( m* E2 \
ans =
+ ]; \: H5 V9 g$ Y  C5 L3 Q. P" N+ _$ i3 s/ @
   Empty matrix: 1-by-0
/ f! n' R! F  y! C- n3 _由于 a 向量中的 0.3 是由 0.1+0.1+0.1 计算得到,在计算过程中就产生了浮点数误差,这也导致在判断 a==0.3 时结果为false,所以 find(a==0.3) 返回一个空矩阵。 7 X4 U9 j' d9 @4 A/ T" z
在进行数值计算判断两数相等时,最好不要直接判断,而是设立一个容差值,当两个浮点数的差的绝对值小于给定的容差值时,我们就认为这两个浮点数相等。 : t2 n8 v, V+ n& R& q' E
比如对于上面的例子:6 w# _  L; u" U
>> a=0.1:0.1:0.5  @* Y1 [+ Y) |& [& _1 _8 ]

4 z" ^4 E$ B* La =
  X. b/ v  {8 F. e# m5 S
' @& N: L6 g" R7 N8 Q7 V  0.1000   0.2000   0.3000   0.4000   0.5000& Z( |* n/ ]& @9 Y6 n1 {" l% x

/ G2 ^: j0 n) x1 u>> tol=eps(0.3)*10   %设立容差值,一般比这个点的浮点数误差高一到两个数量级即可。eps函数能够求得该点的浮点数误差值。& U, C( h% T9 v7 l' z
# z# Q9 Z0 [& d+ E
tol = & o7 N1 r# o' R, q  |
% \0 Y5 y! e1 h: L* p$ V
  5.5511e-15
- R: a& {( N8 y( s8 u; t
2 X% S6 w- I. ^0 Y9 E>> find(abs(a-0.3)
/ h/ |0 I) O, @) I' Z9 N) kans = 6 z/ b6 u. x) l, E' T
; h+ W$ I: n+ _& A$ s$ m
  3" o; D  s% u7 @/ V: r2 z

' i  o- C- z3 x8 h! b$ ~- O8 G8 v: \" P: b! [8 A% K. M9 x
生成一系列有规律名变量
3 u7 w" G+ a! _" R* z/ A5 e当循环迭代需要把每次迭代结果进行保存时,如果每次迭代的结果是尺寸不同的矩阵,无法用矩阵进行存储,那么可以利用 eval 和 num2str 这两个函数可以生成一系列例如 a1、a2、a3… 变量对结果进行保存(不推荐这种方法,原因是 eval 这个函数有很多缺点)。
0 u0 I4 f$ F. Deval::将括号内的字符串视为语句并运行。
' \9 w$ X2 J* Z# M* cnum2str: 将数值转换为字符串。
% K' J; `; J8 W! X9 E% B%例- p. y! _4 k1 b4 y
for ii=1:10
; T! w. s/ `5 v) [, F( g    str=['a',num2str(ii),'=1:',num2str(ii)];
# ]# l4 Y6 U* d  {! B    eval(str)" L5 L" K2 m' R: B2 w, I
end
" K% Y& h& {4 U( o8 t/ o这样可以生成一系列变量 a1、a2…a10 对循环结果进行保存。
, c- h  ~" ^; |# ^7 E% B不推荐使用 eval 函数的原因,帮助文档有详细的解释。5 a+ S9 L9 ~  u
MATLAB® 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.: @+ _% p- [1 n% ]  o$ J
Code within an eval statement can unexpectedly create or assign to a variable already in the current workspace, overwriting existing data.9 d5 ?' e' i4 p" B5 R  @
Concatenating strings within an eval statement is often difficult to read. Other language constructs can simplify the syntax in your code.1 [: Q2 `1 n7 r1 {
MATLAB 对于这类问题有更好的解决办法,利用元胞数组对结果进行存储。元胞数组是 MATLAB 中的特色数据类型,它的元素可以是任意类型的变量,包括不同尺寸或不同维度的矩阵。 对于上面的例子,利用元胞数组:
6 r- v6 N7 q! ifor ii=1:10" a! i9 p  R+ t+ |
    a{ii}=1:ii;
* Q# \( d# y' W9 n0 X1 X% _end
$ k& `8 [) L0 G% h' F: h即生成一系列元胞存储循环结果。这样无论是程序的可读性、运行效率还是后续程序对保存结果调用的方便程度,都远胜于 eval 函数。
: Z8 }; q& ~" O除此之外,在处理符号变量时如果需要生成一系列有规律名符号变量,例如生成一个多项式:" P. Q/ v2 ~9 j4 L9 u
y = x1+2*x2+3*x3+…+100*x100, C, ~6 j; c) l* n" b
eval+num2str 能够实现,但更简便的方法还是利用矩阵:
: P2 v( B2 |! ux=sym('x',[1,100]);! ?3 \3 }" y% [! Q$ H
w=(1:100).*x;
9 u; y  I8 I& o" l! q* h" G) f# \y=sum(w)
6 p9 ^# a2 i5 J. h1 F" h4 T# h; I) y9 ]/ E

) U' V! b) ?; Z; a4 e
作者: youOK    时间: 2021-2-19 18:42
MATLAB求解方程之trapz函数




欢迎光临 EDA365电子论坛网 (https://bbs.eda365.com/) Powered by Discuz! X3.2