找回密码
 注册
关于网站域名变更的通知
查看: 615|回复: 1
打印 上一主题 下一主题

MATLAB求解方程之trapz函数

[复制链接]

该用户从未签到

跳转到指定楼层
1#
发表于 2021-2-19 18:02 | 只看该作者 回帖奖励 |倒序浏览 |阅读模式

EDA365欢迎您登录!

您需要 登录 才可以下载或查看,没有帐号?注册

x
+ H  ?3 Q; Q1 h1 e& {" v
trapz 是基于梯形法则的离散点积分函数。 调用形式:! K( Q- C( p: ]+ S- Y
I = trapz(x,y)
( N  Q3 G; _3 h6 y! l. D  t其中 x 和 y 分别是自变量和对应函数值,以 sin(x) 在 [0,pi] 积分为例:  t7 u5 Q/ Z3 F7 ?8 L" v
x = linspace(0,pi,1e3);     %生成 [0,pi] 内的一系列离散点
& q7 t. i: l/ Q. }* h0 f; F" jy = sin(x);
# U/ J( ~, Q( W& Z, p+ L5 B. UI = trapz(x,y)
# t% d9 A+ D7 S" t7 ?( `' s. x& w2 }' N- c7 ?
浮点数误差8 N# k: [. r; r! J. j7 t
由于计算机中都是以二进制形式存储数据,当用十进制数进行计算时,就会存在十进制数二进制数的转换。但是某些十进制数转化为二进制数是一个无限位的小数,这时必须对该无限位小数进行截断计算机才能存储,那么此时就会产生误差,即浮点数误差。
5 o& A# v! j# k# n9 @/ g5 a例如十进制的0.9,在转化为二进制时是无限循环小数0.1110011001100110011...。此时须对该无限位小数进行截断才能保存在内存当中,截断后再转换回十进制时,0.9就变成了0.90000000000000002,这就是浮点数误差的产生过程。 % G9 |. V* g- m: g$ r# P3 f. f
由于浮点数误差的存在,当进行数值计算时就会出现一些不可避免的问题,最常见的就是判断两数相等时得到与预期相反的结果。 - ?3 F0 E+ H! U, ?# |" X* w
例:令 a = 0.1+0.2, b = 0.3, 判断 a==b 时,MATLAB 会返回0, 当执行 a-b 时,会发现结果不是精确等于0,而是一个非常小的数5.5511e-17。
' G! _/ L5 I! ]或者在矩阵中寻找数的位置(也相当于是判断两数相等)。
2 h  ~. N9 E" v, R- ^7 W2 A%例+ t& @" y2 |4 v+ ?( e, L2 {
>> a = 0.1:0.1:0.5
; U, e% {+ S( U' |  q% W1 g2 {9 j: u' U/ k" q$ }' O3 z' _
a =
1 D( y) y: s+ W  x( R
2 F" A, K/ o6 H$ p8 `    0.1000    0.2000    0.3000    0.4000    0.5000' g/ b! {7 q/ C. H# p
- q+ {/ k- Y; z; h
>> find(a==0.3)
1 A8 T0 R0 x! `& ^" `8 r7 d
7 Z0 Z1 _% g) m3 D) Jans =7 O2 ?, a  J3 C1 T. t5 k' n; Y

! t5 K" K: W7 R; b" {# p8 x7 B   Empty matrix: 1-by-0
3 U5 h0 t6 X4 H, I3 j' {# T" u% M由于 a 向量中的 0.3 是由 0.1+0.1+0.1 计算得到,在计算过程中就产生了浮点数误差,这也导致在判断 a==0.3 时结果为false,所以 find(a==0.3) 返回一个空矩阵。 1 V. W" A+ a; R& h
在进行数值计算判断两数相等时,最好不要直接判断,而是设立一个容差值,当两个浮点数的差的绝对值小于给定的容差值时,我们就认为这两个浮点数相等。
+ C9 ^& z: H* Z比如对于上面的例子:8 C) b: q8 L1 n% A6 o8 N
>> a=0.1:0.1:0.5. {2 M6 ?# ^1 x& ^  ~6 f8 v: Z

$ m% q. [! i! z+ ia =) [; j' @! o& i

& B5 n) a( |: s* p7 P  0.1000   0.2000   0.3000   0.4000   0.50009 v4 Y+ @* @* U/ G" Z) y
0 d. r& @+ ]+ ^. T
>> tol=eps(0.3)*10   %设立容差值,一般比这个点的浮点数误差高一到两个数量级即可。eps函数能够求得该点的浮点数误差值。2 R6 _* h$ i0 ^
' V/ H, N! C! v  P2 ?. |4 ~
tol = 0 f4 I( _, M5 ~# {  ?9 ^; `# B5 M

, h7 _" @* q8 i  5.5511e-15
& I: G' {; P8 o/ I7 K8 j ; w9 e) C# K( H6 _! A, N7 D5 G. P
>> find(abs(a-0.3)
4 b8 K( Y* p1 z; R! d* N4 ^4 Tans = 2 [6 D" U( w5 N2 u3 B9 _! E" G

; \: E3 b% p. ^  A+ |  3
/ {' {- O% {- w+ w0 I$ L* d& H9 @2 |& u' u

: r  a1 I" k$ T+ P  U0 b7 |生成一系列有规律名变量
' f% |2 g7 p6 S- I当循环迭代需要把每次迭代结果进行保存时,如果每次迭代的结果是尺寸不同的矩阵,无法用矩阵进行存储,那么可以利用 eval 和 num2str 这两个函数可以生成一系列例如 a1、a2、a3… 变量对结果进行保存(不推荐这种方法,原因是 eval 这个函数有很多缺点)。
) V) M! Y9 T% b2 T2 zeval::将括号内的字符串视为语句并运行。! T: C+ \  p5 r3 v
num2str: 将数值转换为字符串。! g2 N# ~% E$ E* K- C( L
%例
" Y/ |( c- H+ pfor ii=1:10
$ k- X1 a4 W7 l' [    str=['a',num2str(ii),'=1:',num2str(ii)];# |* p! f. s; h& Q! f6 q
    eval(str)8 h8 x& |; }  F
end( W6 y6 e2 ?7 T
这样可以生成一系列变量 a1、a2…a10 对循环结果进行保存。 0 m  j. X" C; N; m
不推荐使用 eval 函数的原因,帮助文档有详细的解释。& k' @7 p( [0 I! R) G
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.* c2 x" P5 l; g+ u% A
Code within an eval statement can unexpectedly create or assign to a variable already in the current workspace, overwriting existing data.
6 B1 a( {: I+ O" Q7 d4 UConcatenating strings within an eval statement is often difficult to read. Other language constructs can simplify the syntax in your code.
& B* W0 g# ?& l7 o1 uMATLAB 对于这类问题有更好的解决办法,利用元胞数组对结果进行存储。元胞数组是 MATLAB 中的特色数据类型,它的元素可以是任意类型的变量,包括不同尺寸或不同维度的矩阵。 对于上面的例子,利用元胞数组:
" w! i* b5 N8 @% H& Q- `for ii=1:108 f/ Z' c" `" V9 O! {- }
    a{ii}=1:ii;+ F% e& w2 S, r/ K' C4 V
end
4 ~) d1 \  d6 Z) R' {/ X/ Y即生成一系列元胞存储循环结果。这样无论是程序的可读性、运行效率还是后续程序对保存结果调用的方便程度,都远胜于 eval 函数。
) \7 s9 m# P, @+ h$ m7 W除此之外,在处理符号变量时如果需要生成一系列有规律名符号变量,例如生成一个多项式:
- B$ z9 i8 Y1 @& R9 K+ U( ]y = x1+2*x2+3*x3+…+100*x100/ U7 Q  x: Z& b8 m, N8 {, P3 t
eval+num2str 能够实现,但更简便的方法还是利用矩阵:; F: G, t. [, ^9 N* d, c% w) Y
x=sym('x',[1,100]);, p4 x6 x: T; o: J
w=(1:100).*x;
! P7 |4 s2 A% n/ {y=sum(w)# d: }7 p9 q; ^' o) X2 v
+ n+ e8 e1 v3 J# J" \' s# p8 X

: k4 A( y( `) M9 {

该用户从未签到

2#
发表于 2021-2-19 18:42 | 只看该作者
MATLAB求解方程之trapz函数
您需要登录后才可以回帖 登录 | 注册

本版积分规则

关闭

推荐内容上一条 /1 下一条

EDA365公众号

关于我们|手机版|EDA365电子论坛网 ( 粤ICP备18020198号-1 )

GMT+8, 2025-6-20 05:09 , Processed in 0.062500 second(s), 23 queries , Gzip On.

深圳市墨知创新科技有限公司

地址:深圳市南山区科技生态园2栋A座805 电话:19926409050

快速回复 返回顶部 返回列表