应《网络安全法》要求，自2017年10月1日起，未进行实名认证将不得使用互联网跟帖服务。为保障您的帐号能够正常使用，请尽快对帐号进行手机号验证，感谢您的理解与支持！

24小时热门版块排行榜    

小木虫论坛-学术科研互动平台 » 计算模拟区 » 计算模拟 » 求高版本1stopt拟合，
51/1返回列表
查看: 2365  |  回复: 11
只看楼主 @他人 存档 新回复提醒 (忽略) 收藏    在APP中查看
当前只显示满足指定条件的回帖，点击这里查看本话题的所有回帖

yupeifeng

金虫 (小有名气)

[求助] 求高版本1stopt拟合，

给出拟合后的图像及参数、相关系数值，以便确定选择哪种模型。程序已经编好，只需运行后，把拟合后的图像，参数值及相关系数R复制过来就行了，谢谢！！
一、T=40℃
第一个程序：
Constant c0=6.9401,cp=5.8399;
Variable t,c;
ODEFunction c'=-k*(c0-c)^(1/3)*(c-cp)
Data;
t,c
0        6.94
1        6.45
2        6.35
4        6.2
12        5.96
24        5.85
48        5.84

第二个程序
Constant c0=6.9401,cp=5.8399;
Variable t,c;
ODEFunction c'=-k*(c0-c)^(2/3)*(c-cp)^p
Data;
t,c
0        6.94
1        6.45
2        6.35
4        6.2
12        5.96
24        5.85
48        5.84

第三个程序
Constant c0=6.9401,cp=5.8399;
Variable t,c;
ODEFunction c'=-k*(c0-c)^(4/3)*(c-cp)^p
Data;
t,c
0        6.94
1        6.45
2        6.35
4        6.2
12        5.96
24        5.85
48        5.84

第四个程序
Constant c0=6.9401,cp=5.8399;
Variable t,c;
ODEFunction c'=-k*(c-cp)^n
Data;
t,c
0        6.94
1        6.45
2        6.35
4        6.2
12        5.96
24        5.85
48        5.84

二、T=30℃
第一个程序：
Constant c0=6.9401,cp=6.0499;
Variable t,c;
ODEFunction c'=-k*(c0-c)^(1/3)*(c-cp)
Data;
t,c
0        6.94
2        6.39
4        6.3
8        6.2
12        6.1
48        6.05

第二个程序
Constant c0=6.9401,cp=6.0499;
Variable t,c;
ODEFunction c'=-k*(c0-c)^(2/3)*(c-cp)^p
Data;
t,c
0        6.94
2        6.39
4        6.3
8        6.2
12        6.1
48        6.05

第三个程序
Constant c0=6.9401,cp=6.0499;
Variable t,c;
ODEFunction c'=-k*(c0-c)^(4/3)*(c-cp)^p
Data;
t,c
0        6.94
2        6.39
4        6.3
8        6.2
12        6.1
48        6.05

第四个程序
Constant c0=6.9401,cp=6.0499;
Variable t,c;
ODEFunction c'=-k*(c-cp)^n
Data;
t,c
0        6.94
2        6.39
4        6.3
8        6.2
12        6.1
48        6.05

三、T=50℃
第一个程序：
Constant c0=6.7901,cp=5.7899;
Variable t,c;
ODEFunction c'=-k*(c0-c)^(1/3)*(c-cp)
Data;
t,c
0        6.79
1        6.39
2        6.25
4        6.1
12        5.9
22        5.8
48        5.79
第二个程序：
Constant c0=6.7901,cp=5.7899;
Variable t,c;
ODEFunction c'=-k*(c0-c)^(2/3)*(c-cp)^p
Data;
t,c
0        6.79
1        6.39
2        6.25
4        6.1
12        5.9
22        5.8
48        5.79
第三个程序：
Constant c0=6.7901,cp=5.7899;
Variable t,c;
ODEFunction c'=-k*(c0-c)^(4/3)*(c-cp)^p
Data;
t,c
0        6.79
1        6.39
2        6.25
4        6.1
12        5.9
22        5.8
48        5.79
第四个程序：
Constant c0=6.7901,cp=5.7899;
Variable t,c;
ODEFunction c'=-k*(c-cp)^n
Data;
t,c
0        6.79
1        6.39
2        6.25
4        6.1
12        5.9
22        5.8
48        5.79[ Last edited by fegg7502 on 2013-4-11 at 10:34 ]
回复此楼

» 本帖附件资源列表

  • 欢迎监督和反馈：小木虫仅提供交流平台，不对该内容负责。
    本内容由用户自主发布，如果其内容涉及到知识产权问题，其责任在于用户本人，如对版权有异议，请联系邮箱：xiaomuchong@tal.com
  • 附件 1 : 1stopt.doc
    • 2013-04-09 09:21:26, 35 K

» 猜你喜欢

» 本主题相关价值贴推荐，对您同样有帮助:
憧憬未来
1楼 2013-04-09 09:25:30
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

月只蓝

主管区长 (职业作家)
此外，第一个和四个方程，可以先积分为代数方程，再用低版本1stopt或者origin spss excel等拟合也可。
赞一下 回复此楼
MATLAB、MS小问题、普通问题请发帖求助！时间精力有限，恕不接受无偿私信求助。
3楼2013-04-09 11:03:06
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
查看全部 12 个回答

月只蓝

主管区长 (职业作家)

【答案】应助回帖

★ ★
csgt0: 金币+2, 应助指数+1, 谢谢 2013-04-09 11:07:07
手头没有高版本的1stopt，用MATLAB作了一下，因为一般自编的MATLAB程序拟合结果被初值影响很大，第三情况方程的刚性比较大，结果仅供参考。T=40℃的情况，分别见附图1,2,3,4。目前看来，第四个方程比较好。

附图1.jpg



附图2.jpg



附图3.jpg



附图4.jpg
赞一下 回复此楼
MATLAB、MS小问题、普通问题请发帖求助！时间精力有限，恕不接受无偿私信求助。
2楼2013-04-09 10:58:55
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

yupeifeng

金虫 (小有名气)
引用回帖:
3楼: Originally posted by 月只蓝 at 2013-04-09 11:03:06
此外，第一个和四个方程，可以先积分为代数方程，再用低版本1stopt或者origin spss excel等拟合也可。

可不可以让我看看你的程序怎么编写的啊，谢谢了
赞一下 回复此楼
憧憬未来
4楼2013-04-09 21:31:23
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

月只蓝

主管区长 (职业作家)
★ ★
fegg7502: 金币+2, 鼓励交流 2013-04-11 10:35:08
引用回帖:
4楼: Originally posted by yupeifeng at 2013-04-09 21:31:23
可不可以让我看看你的程序怎么编写的啊，谢谢了...

第一个方程
%--------------------------start-----------------------------------------
function ode2
clear all;clc

format long

tspan=[0 1 2 4 12 24 48];
t1=length(tspan),
yexp=[6.45 6.35 6.2 5.96 5.85 5.84]';
t2=length(yexp)



y0=6.94;

k0=[2];
lb=[-1]*1e3;
ub=[1]*1e3;        

yy=[y0 yexp'];

[k,resnorm,residual,exitflag,output,lambda,jacobian] = ...
    lsqnonlin(@ObjFunc,k0,lb,ub,[],tspan,y0,yexp);
ci = nlparci(k,residual,jacobian);
fprintf('\n\n使用函数lsqnonlin()估计得到的参数值为:\n')
fprintf('\t反应速率常数 k = %.4f\n',k(1))
%fprintf('\t反应级数 n = %.4f\n',k(2))
fprintf('  The sum of the squares is: %.1e\n\n',resnorm)
ts=0:1:max(tspan);

[ts ys]=ode45(@KineticsEqs,ts,y0,[],k);
[ttt XXsim] = ode45(@KineticsEqs,tspan,y0,[],k);
y=XXsim(2:end);
xexp=yexp;
R2=1-sum((xexp-y).^2)./sum((xexp-mean(y)).^2);
fprintf('\n\t相关系数之平方R^2 = %.6f',R2);
plot(ts,ys,'b',tspan,yy,'or'),legend('计算值','实验值','Location','best'),


%---------------------------------------------------------
function f = ObjFunc(k,tspan,y0,yexp)           % 目标函数
[t Xsim] = ode45(@KineticsEqs,tspan,y0,[],k);   
ysim = Xsim(2:end);
f=ysim-yexp;
%----------------------------------------------------------

function dydt = KineticsEqs(t,y,k)
c0=6.9401;cp=5.8399;
dydt=-k(1)*(c0-y).^(1/3).*(y-cp);
%--------------------The end-------------------------------------------
赞一下 回复此楼
MATLAB、MS小问题、普通问题请发帖求助！时间精力有限，恕不接受无偿私信求助。
5楼2013-04-09 21:47:34
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
查看全部 12 个回答
普通表情 龙 兔 虎 猫
信息提示
关闭
请填处理意见
关闭