|
|
EDA365欢迎您登录!
您需要 登录 才可以下载或查看,没有帐号?注册
x
本帖最后由 House 于 2020-3-10 17:43 编辑
7 Z& i3 Y9 r( h' z* v8 H& @! Y4 ^, s1 ^; {7 g
MATLAB使用自定义的欧拉法求解常微分方程组
`! _* y# q! Q* K2 F6 z, B4 s3 v0 H" a9 V* s Z
. Y* Y7 ^8 R; |; v2 H+ K, T. D* \
2 g4 G, G6 Q9 X
function F=f(t,Y)
: v2 |2 ^! z3 r! B% 定义待求的微分方程组0 I) b2 p& F2 d# R1 Q* w, ~
x=Y(1);
! g3 }, W; ?" l1 Yy=Y(2);4 Q% p8 F5 `/ `7 x
f1=3*y;3 G7 p' s6 }8 C! J* n' ^5 O
f2=(1-x^2)*y-x;
5 Q* f2 k* M+ }$ D. `F=[f1;f2];
* D# H% h& d, s n& b aend+ q' _5 A% ^0 @
) D `2 W# ]2 l6 q P
o( v; D3 n |7 w# h @%% 定义计算的步长, 设置变量的初始值
- ^( L. [- {6 f9 Vclear;clc;close all
5 v0 E$ _- d5 r; uDelta=0.001; % 定义步长' |6 K/ J# v* V
t=0: Delta :20; % 定义自变量 t
$ _7 a6 d* B7 Y, |, V$ |n=length(t);
" Z" r! a$ o6 M9 aY(:,1)=[2;0]; % 定义 x y 的初始值
( j2 A6 W% b6 D0 F* P+ u4 h0 ~. I4 q& v4 F- [! x: y
%% 自定义欧拉法, 求解微分方程组
; I7 D& ?4 t/ Tfor k=1:n-1
% c2 x+ {: p5 F. _8 O Y(:,k+1)=Y(:,k)+Delta*f(t(k),Y(:,k));
9 [% V5 O0 b% U: T5 ~, H/ g0 U4 t- \end! g$ f5 j# F1 m' f t, ?
x=Y(1,: );( a+ O2 _5 `7 G; ]1 N% E
y=Y(2,: );1 |1 V: T& R2 E' b" r! X9 C
4 Y/ P: Y# ~, e# Z
%% 绘制 x y 的求解结果
* ~1 Y& W; }. w0 _figure
9 d' n! P9 q6 |; ] Y Tset(gcf,'units','normalized','position',[0.15 0.2 0.7 0.6]); % 设置 figure 窗口的位置和尺寸" B0 j' m* Z- B2 M. A: a
subplot(1,2,1)
5 X# e5 I* _- v9 A1 l$ i0 Y3 splot(t,x,'b')
' r8 l0 F. X) wxlabel('t'); c7 W8 H4 y5 g( n
ylabel('x')8 ~$ W* @5 W5 _! N
- c; a* W7 h4 ^& u
subplot(1,2,2)8 U* j" O! B: h, l8 s
plot(t,y,'r')
4 l/ C( `9 k) j! w% kxlabel('t')( |" G$ }6 L1 E/ F) N- m& L+ V
ylabel('y')! s% c; ]/ E4 @6 s' m# [: p' Z
7 I$ z+ N. L2 R' F' b6 q6 D
# i+ U2 u Z4 t3 {7 C/ f9 z( W9 F7 k) D% R
|
|