24小时热门版块排行榜    

CyRhmU.jpeg
查看: 1720  |  回复: 5

Mr_杰克

新虫 (初入文坛)

[求助] [求助] 1stopt 微分方程代码求助已有2人参与

这两天在用1stopt研究一个动力学微分方程的拟合求动力学参数,看起来应该是比较简单的方程,但是研究了很长时间,自己也根据网上能搜到的例子编了两个程序跑了下,但是结果相差比较大,我也不知道程序编的对不对,请懂的大神指点一二,跪谢了!我使用的是试用版的1stopt 5.0

方程是这样的:
dx/dt=-p1*x;
dy/dt=p1*x-p2*y;
初始值:t=0,此时x=3.5,y=0
测量的数据是t 和 y。具体数据请见程序中。p1和p2应为正。

程序1:(dx/dt=-p1*x   这个方程,我先积分成了x=x0*exp(-p1*t),也不知道对不对?)
//y'=dy/dt=3.5*p1*exp(-p1*t)-p2*y
Constant V =0; //Initial Values
Constant dt=1;
Parameter p(2);
Variable t, y;
StartProgram [Pascal];
Procedure MainModel;
var i,t: integer;
    y0 : double;
Begin
    y0 := V;
    for i := 0 to DataLength - 1 do begin
        y0 := (3.5*p1*exp(-p1*t)-p2*y0)*dt+y0;
        y := y0;

    end;
End;
EndProgram;
Data;
//0     0
10     0.010515
14     0.023065
17     0.037836
20     0.062683
22     0.063288
24     0.064385
26     0.063359
28     0.062942
30     0.061384
32     0.060358
35     0.065102
45     0.059101
60     0.043515

这个程序的输出报告和图:[求助] 1stopt 微分方程代码求助
                                          [求助] 1stopt 微分方程代码求助-1

关于程序1,我的主要疑惑是:1.不知道这个程序编的对不对?
                                                2.关于 dt 这个数值,在程序中的作用,我比较担心的是,我的t是不连续取值的,dt会不会使得t在程序中变成连续取值了?从图上看我的这一担心就更明显,因为横坐标只取到了12



程序2:
InitialODEValue t=0,x=3.5,y=0;
Parameters p1,p2;
Variable t[Output],y[Output];
ODEFunction x'=-p1*x;
            y'=p1*x-p2*y;

Data;
//0              0
10                0.010514645
14                0.023065468
17                0.037836003
20                0.062682935
22                0.06328842
24                0.064385481
26                0.063358784
28                0.062942392
30                0.061383892
32                0.060358109
35                0.065102196
45                0.059100765
60                0.043515354

这个程序的输出报告和图:     [求助] 1stopt 微分方程代码求助-2
                                               [求助] 1stopt 微分方程代码求助-3

关于程序2,我的主要疑惑是:1.不知道这个程序编的对不对?
                                                2.关于ODEFunction 语句的使用是否正确?ODEFunction语句对微分函数的定义默认是怎么样的?在函数中我都没定义t,他是怎么计算的?也涉及到t的取值问题,ODEFunction是不是自动把t按照1,2,3,4,5……取值了?


总结一下问题:1.请高手指导我这个动力学微分方程应该如何编程,能给代码就最好了;(授人以鱼,但是是最主要的。。。)
                        2.请指出我自己编写程序的问题,以及关于两个程序中t的取值,是否连续?还是间歇,如何确定(授人以渔)     
                        3.关于ODEFunction的定义,关于程序1中dt的定义。


谢谢!!!
回复此楼

» 猜你喜欢

» 本主题相关价值贴推荐,对您同样有帮助:

已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

月只蓝

主管区长 (职业作家)

【答案】应助回帖

感谢参与,应助指数 +1
去看看这位大神关于 微分方程拟合的回帖吧:http://muchong.com/bbs/space.php?uid=291104&view=post
MATLAB、MS小问题、普通问题请发帖求助!时间精力有限,恕不接受无偿私信求助。
2楼2015-04-08 10:40:32
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

Mr_杰克

新虫 (初入文坛)

引用回帖:
2楼: Originally posted by 月只蓝 at 2015-04-08 10:40:32
去看看这位大神关于 微分方程拟合的回帖吧:http://muchong.com/bbs/space.php?uid=291104&view=post

谢谢版大关注!您是使用Matlab的是吧?

嗯,他的帖子也看了不少,但还是有不少疑惑,我已经站短他了,希望大神本人能现身说法吧。谢谢!
3楼2015-04-08 10:53:12
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

dingd

铁杆木虫 (职业作家)

【答案】应助回帖

感谢参与,应助指数 +1
CODE:
InitialODEValue t=0,x=3.5,y=0;
Parameters p1,p2;
Variable t,y;
ODEFunction x'=-p1*x;
            y'=p1*x-p2*y;
Data;
//0              0
10                0.010514645
14                0.023065468
17                0.037836003
20                0.062682935
22                0.06328842
24                0.064385481
26                0.063358784
28                0.062942392
30                0.061383892
32                0.060358109
35                0.065102196
45                0.059100765
60                0.043515354

结果稳定唯一,但效果不是很好:

Root of Mean Square Error (RMSE):0.0119983163868176
Sum of Square Error:0.00187147474953628
Correlation Coef. (R): 0.77730265287911
R-Square: 0.604199414172902
Determination Coef. (DC): 0.514455569529589
F-Statistic: 22.6549813151864

Parameter                  Best Estimate
--------------------        -------------
p1        0.00131507520277041
p2        0.0675323910740055
[求助] 1stopt 微分方程代码求助-4
c2.jpg

4楼2015-04-08 11:32:00
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

Mr_杰克

新虫 (初入文坛)

引用回帖:
4楼: Originally posted by dingd at 2015-04-08 11:32:00
InitialODEValue t=0,x=3.5,y=0;
Parameters p1,p2;
Variable t,y;
ODEFunction x'=-p1*x;
            y'=p1*x-p2*y;
Data;
//0              0
10                0.010514645
14                0.023 ...

谢谢!跟我的程序2一样,结果也一样看来我编的是对的?那能麻烦大神解答下程序1为什么跑出来的结果不一样吗?
5楼2015-04-08 12:05:38
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

Mr_杰克

新虫 (初入文坛)

引用回帖:
4楼: Originally posted by dingd at 2015-04-08 11:32:00
InitialODEValue t=0,x=3.5,y=0;
Parameters p1,p2;
Variable t,y;
ODEFunction x'=-p1*x;
            y'=p1*x-p2*y;
Data;
//0              0
10                0.010514645
14                0.023 ...

另外再请教一个问题,1stopt能否像Origin非线性拟合时那样给出拟合的参数的标准误?
6楼2015-04-08 12:07:39
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 Mr_杰克 的主题更新
信息提示
请填处理意见