数学建模社区-数学中国
标题:
湖泊有机物新陈代谢模型问题
[打印本页]
作者:
度过无风
时间:
2013-1-25 19:34
标题:
湖泊有机物新陈代谢模型问题
%微分方程组求解主程序
( a0 M( R* P; S( z0 u" L. X- H0 Q
clc;clear all;clf;close all;
( e8 h [0 V. n" Q$ ?
%微分方程组求解主程序
! j$ @# }( z/ w# J# F/ p
clc;clear all;clf;close all;
7 ~3 \6 D- u0 y& b
%windows时钟自动计时
1 y2 {# D6 X4 h
T1=clock; %函数clock返回的值是clock=[year month day hour minute seconds]
$ a2 I8 [9 z$ v. }- r" _8 q
disp('计算机正在准备输出湖泊有机物新陈代谢结果，请耐心等待······');
% X( m V, B( I8 Q
[tt,y]=ode45('lbwfun',[0:10:2020],[95.9,0.83,0.003,0.0001,0.0,0.0]);
: v! X) T' I7 G4 q" ]
t=tt(191:end,
' o" f1 I4 h0 b. \
ys=(191:end,1)
: x3 o. b2 {9 ]. g' V
yp=(191:end,2)
8 |: w% q! p) a5 o. c
yh=(191:end,3)
4 v" i) M- h/ ~$ y. i5 a6 G7 }8 F3 W. Z
yr=(191:end,4)
& }9 y1 F! \0 d2 o3 [5 |0 P
yo=(191:end,5)
' u% L* n6 ~ K% R% r
ye=(191:end,6)
; Y3 c0 n4 o; w8 Z6 ~# [! f
T2=clock;
1 t0 R& x6 g: g2 {/ ~0 _
API_elapsed_time=T2-T1;
1 L& c# W+ B6 I- M Q/ L" t
if API_elapsed_time(6)<0
0 Q: T' s- X( j3 n
API_elapsed_time(6)=API_elapsed_time(6)+60;
. F( x' I) U1 `- i0 V2 @
API_elapsed_time(5)=API_elapsed_time(5)-1;
: i0 p3 Z& h( ?$ K% D& V
end
' P3 [$ m* N8 \ ?/ Z0 N
if API_elapsed_time(4)<0
. L* D' c9 ?( @
API_elapsed_time(4)=API_elapsed_time(4)+60;
1 L( L2 I+ }$ j9 t" ^7 G9 O' [
API_elapsed_time(4)=API_elapsed_time(4)-1;
( z7 G1 g. k2 g( W7 r
end
4 E' J% b% s& s: O' w
str=sprintf('湖泊新陈代谢模拟程序共运行 %d 小时 %d 分钟 %.4f 秒',```
: X9 m2 D3 w1 y, w9 S* w) G
API_elapsed_time(4),API_elapsed_time(5),API_elapsed_time(6));
! m8 B+ r$ k* A! F6 x
disp(str);
$ ~& V; l: T5 \" ]4 T3 q
%子函数：微分方程组odefile文件并命名为lbwfun.m
9 L& _" L& ]& X$ A( L" _
function ydot=lbwfun(t,y);
0 r" L. ]2 B2 @( ` y0 H: t
ydot=[121.793*pi*cos(2*pi*t);
- @0 l4 }( ^+ U, {. i' F3 l
y(1)-4.03*y(2);
5 E$ t+ Y1 w! T
0.48*y(2)-17.87*y(3);
, X" Y* W2 @6 ], v0 p
4.85*y(3)-4.65*y(4);
# \0 V$ E0 K: D t/ g
2.55*y(2)+6.12*y(3)+1.95*y(4);
2 Q, w4 j" s3 \9 J0 J& Q
1.10*y(2)+6.90*y(3)+2.70*y(4)];
' h1 D; s/ E$ P5 K, ]3 W* f! \
运行结果（我是按书上写出来的程序）
8 Y- P* a7 [2 O ^0 G5 G! W- h
计算机正在准备输出湖泊有机物新陈代谢结果，请耐心等待······
$ e# H2 E" W& R+ B% M& [: `
??? Error using ==> feval
/ ]1 T0 t: T+ V0 M3 @" T
Undefined function or method 'lbwfun' for input arguments of type 'double'.
* Z5 ^& b F$ G6 \$ @% h. z/ w# o( T8 A
$ {. R% U X7 E' o
Error in ==> odearguments at 110
C2 B* _! J2 k3 P
f0 = feval(ode,t0,y0,args{:}); % ODE15I sets args{1} to yp0.
& Z* C6 o; ^, {* E* }
2 Z+ ^8 M# t' k* k4 ^
Error in ==> ode45 at 173
7 o. m/ @, S- q( c) W
[neq, tspan, ntspan, next, t0, tfinal, tdir, y0, f0, odeArgs, odeFcn, ...
欢迎光临 数学建模社区-数学中国 (http://www.madio.net/)
Powered by Discuz! X2.5
