24小时热门版块排行榜    

查看: 2222  |  回复: 14

l1003785517

新虫 (小有名气)

[求助] 用1stopt计算非线性方程组的未知参数

各位虫友能否帮忙看看
        代码如下,其中a1和a2都可以写成k=k0.exp^(-Ea/RT),即反应速率常数,T是变量。我把参数合并是因为会出现过拟合现象,但即使是这样,我算出来还是会出现参数结果不一样,而且我是希望能够求出k0和Ea来的,不知道有没有什么解决的办法,下面的代码里面只放了一个温度下的实验数据。我把另外几个温度的数据也放到代码的后面了,各位虫友看能帮我解决一下这个难题吗,若有可能帮忙分别计算一下a1和a2中的k0和Ea?
        因为后面还要对不同的动力学方程进行计算,所以想问一下这个是不是得用高阶付费的1stopt软件才能做呢?我用的网上注册版1.5的无法得到正确的结果。谢谢!如果不需要的话还望虫友在解答的过程中能发一下原代码,万分感谢,这里我把剩余的金币都奉上!还望相助,谢谢!
//用1stopt计算非线性方程组的未知参数
Title "kinetic values";
Variables x(1:3),y(1:2);
Parameters a(1:5);
ShareModel;
Function y1=a1*x1/(1+a3*x1+a4*x2+a5*x3);
          y2=(a2*x2-a1*x1)/(1+a3*x1+a4*x2+a5*x2);
Data;
//T=275度 x1,x2,x3;y1,y2;
6.31317        0.029302 0.007644 0.0574 0.0345
6.29412        0.044590 0.011466 0.0407 0.0298
6.27507        0.059878 0.015288 0.0241 0.0251
6.26237        0.069433 0.018473 0.0278 0.0203
6.24713        0.080899 0.022295 0.0315 0.0208
6.23316        0.091091 0.026117 0.0268 0.0194
6.22109 0.099372 0.029939 0.0220 0.0160
6.21093        0.107016 0.032487 0.0172 0.0127
//另外几个实验数据
//T=282度
6.14807        0.133133        0.069433        0.2894        0.1089
6.062345        0.186004        0.102557        0.1519        0.0714
6.021705        0.202566        0.126763        0.0144        0.0338
5.991225        0.212758        0.147147        0.0425        -0.0038
5.968365        0.220402        0.162435        0.0705        0.0217
5.937885        0.234416        0.178997        0.0575        0.0297
5.91566        0.244608        0.1911        0.0445        0.0201
5.892165        0.254163        0.205114        0.0316        0.0105
//T=287度
6.0452        0.187278        0.118482        0.4615        0.1918
5.89026        0.277095        0.184093        0.273        0.1397
5.80136        0.326781        0.223587        0.0845        0.0876
5.74675        0.356083        0.249067        0.0829        0.0356
5.69976        0.377741        0.274547        0.0813        0.0243
5.67182        0.386022        0.294294        0.0677        0.0176
5.648325        0.390481        0.313404        0.0541        0.0152
5.614035        0.400673        0.33761        0.0404        0.0129
//T=292度
5.90677        0.220402        0.224224        0.6236        0.1761
5.741035        0.303212        0.307671        0.4216        0.1584
5.559425        0.383474        0.409591        0.2197        0.1408
5.434965        0.437619        0.480298        0.2059        0.1231
5.311775        0.497497        0.543998        0.1921        0.0677
5.26669        0.514696        0.572026        0.1219        0.0332
5.229225        0.527436        0.596869        0.0817        0.0195
5.201285        0.53508        0.617253        0.0585        0.0058@月只蓝@beefly
回复此楼
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
回帖支持 ( 显示支持度最高的前 50 名 )

l1003785517

新虫 (小有名气)

引用回帖:
3楼: Originally posted by 月只蓝 at 2018-07-04 07:36:35
如2楼所说,公式k=k0.exp^(-Ea/RT)与a的关系说明不清楚,而且此公式的写法本身就是错误的,正确的是“k=k0*exp(-Ea/R/T),作为求助者,像已知的常数R=8.314应该给出。

...

是的,作为一个求助者真的不应该出现这种错误,很抱歉。
我在代码中给出的公式,Function y1=a1*x1/(1+a3*x1+a4*x2+a5*x3);
                                                    y2=(a2*x2-a1*x1)/(1+a3*x1+a4*x2+a5*x2);
     a1=k1*exp^(-Ea1/R/T);
     a2=k2*exp^(-Ea2/R/T);
     这是一个阿伦尼乌斯公式,k1,k2都是指前因子,Ea1,Ea2是活化能,这几个参数都是我希望求出来的,但是若果放到公式进去会出现过拟合,所以我才把它们合并了放在代码里。R是常数8.314,T是反应温度,代码中我给出了4个温度下的实验数据。
    希望能够得到您的帮助?万分感谢!
6楼2018-07-04 12:55:03
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

dingd

铁杆木虫 (职业作家)

【答案】应助回帖

★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★
独孤神宇: 金币+1 2018-07-05 14:30:10
l1003785517: 金币+20, ★★★很有帮助 2018-07-05 15:28:26
CODE:
Title "kinetic values";
Constant R=8.314;
VarConstant T=[275,282,287,292];
ConstStr a1=k1*exp(-Ea1/R/T),a2=k2*exp(-Ea2/R/T);
Variables x(1:3),y(1:2);
Parameters a(3:5);
ParameterDomain = [0,];
ShareModel;
Function y1=a1*x1/(1+a3*x1+a4*x2+a5*x3);
         y2=(a2*x2-a1*x1)/(1+a3*x1+a4*x2+a5*x2);
Data; //T=275
6.31317        0.029302 0.007644 0.0574 0.0345
6.29412        0.044590 0.011466 0.0407 0.0298
6.27507        0.059878 0.015288 0.0241 0.0251
6.26237        0.069433 0.018473 0.0278 0.0203
6.24713        0.080899 0.022295 0.0315 0.0208
6.23316        0.091091 0.026117 0.0268 0.0194
6.22109 0.099372 0.029939 0.0220 0.0160
6.21093        0.107016 0.032487 0.0172 0.0127
Data; //T=282
6.14807        0.133133        0.069433        0.2894        0.1089
6.062345        0.186004        0.102557        0.1519        0.0714
6.021705        0.202566        0.126763        0.0144        0.0338
5.991225        0.212758        0.147147        0.0425        -0.0038
5.968365        0.220402        0.162435        0.0705        0.0217
5.937885        0.234416        0.178997        0.0575        0.0297
5.91566        0.244608        0.1911        0.0445        0.0201
5.892165        0.254163        0.205114        0.0316        0.0105
Data; //T=287
6.0452        0.187278        0.118482        0.4615        0.1918
5.89026        0.277095        0.184093        0.273        0.1397
5.80136        0.326781        0.223587        0.0845        0.0876
5.74675        0.356083        0.249067        0.0829        0.0356
5.69976        0.377741        0.274547        0.0813        0.0243
5.67182        0.386022        0.294294        0.0677        0.0176
5.648325        0.390481        0.313404        0.0541        0.0152
5.614035        0.400673        0.33761        0.0404        0.0129
Data; //T=292
5.90677        0.220402        0.224224        0.6236        0.1761
5.741035        0.303212        0.307671        0.4216        0.1584
5.559425        0.383474        0.409591        0.2197        0.1408
5.434965        0.437619        0.480298        0.2059        0.1231
5.311775        0.497497        0.543998        0.1921        0.0677
5.26669        0.514696        0.572026        0.1219        0.0332
5.229225        0.527436        0.596869        0.0817        0.0195
5.201285        0.53508        0.617253        0.0585        0.0058

结果不好,自己好好检查下公式和数据:

Root of Mean Square Error (RMSE): 0.0829320200375007
Sum of Squared Residual: 0.440174076640027
Correlation Coef. (R): 0.0189691245411682
R-Square: 0.000359827685858348
Adjusted R-Square: 0.51886698386634
Determination Coef. (DC): -2.88389955816035
F-Statistic: 0.226197729277559

Parameter                  Best Estimate
--------------------        -------------
k1        4.60347494177944E35
ea1        166179.343552407
k2        9.73280643547253E23
ea2        93748.2902053327
a3        65737.9703181681
a4        4.29079677286443E-8
a5        49386056.7270752
11楼2018-07-04 23:05:52
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
普通回帖

dingd

铁杆木虫 (职业作家)

【答案】应助回帖

★ ★ ★ ★ ★ ★ ★ ★ ★ ★
感谢参与,应助指数 +1
月只蓝: 金币+5, 感谢指导! 2018-07-04 07:39:00
月只蓝: 金币+5, 感谢指导! 2018-07-04 07:39:33
最好把公式写清楚:
Function y1=a1*x1/(1+a3*x1+a4*x2+a5*x3);
          y2=(a2*x2-a1*x1)/(1+a3*x1+a4*x2+a5*x2);

与:
k=k0.exp^(-Ea/RT)
都是什么关系?公式“k=k0.exp^(-Ea/RT)”应该是“k=k0.exp^(-Ea/(R*T))”吧?R又是多少?
2楼2018-07-03 22:40:05
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

月只蓝

主管区长 (职业作家)

如2楼所说,公式k=k0.exp^(-Ea/RT)与a的关系说明不清楚,而且此公式的写法本身就是错误的,正确的是“k=k0*exp(-Ea/R/T),作为求助者,像已知的常数R=8.314应该给出。

发自小木虫Android客户端
MATLAB、MS小问题、普通问题请发帖求助!时间精力有限,恕不接受无偿私信求助。
3楼2018-07-04 07:36:35
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

l1003785517

新虫 (小有名气)

引用回帖:
2楼: Originally posted by dingd at 2018-07-03 22:40:05
最好把公式写清楚:
Function y1=a1*x1/(1+a3*x1+a4*x2+a5*x3);
          y2=(a2*x2-a1*x1)/(1+a3*x1+a4*x2+a5*x2);

与:
k=k0.exp^(-Ea/RT)
都是什么关系?公式“k=k0.exp^(-Ea/RT)”应该是“k=k0.exp^(-E ...

首先感谢您的回复,对我的描述不清楚表示抱歉。
      代码中的公式,Function y1=a1*x1/(1+a3*x1+a4*x2+a5*x3);
                                             y2=(a2*x2-a1*x1)/(1+a3*x1+a4*x2+a5*x2);
     a1=k1.exp^(-Ea1/RT);
     a2=k2.exp^(-Ea2/RT);
     这是一个阿伦尼乌斯公式,k1,k2都是指前因子,Ea1,Ea2是活化能,这几个参数都是我希望求出来的,但是若果放到公式进去会出现过拟合,所以我才把它们合并了放在代码里。R是常数8.314,T是反应温度,代码中我给出了4个温度下的实验数据。
     不知道您是否能帮我解决这个问题呢?万分感谢!
4楼2018-07-04 12:48:46
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

l1003785517

新虫 (小有名气)

引用回帖:
2楼: Originally posted by dingd at 2018-07-03 22:40:05
最好把公式写清楚:
Function y1=a1*x1/(1+a3*x1+a4*x2+a5*x3);
          y2=(a2*x2-a1*x1)/(1+a3*x1+a4*x2+a5*x2);

与:
k=k0.exp^(-Ea/RT)
都是什么关系?公式“k=k0.exp^(-Ea/RT)”应该是“k=k0.exp^(-E ...

不好意思这两个公式应该是这样,我写错了。
a1=k1*exp^(-Ea1/R/T);
a2=k2*exp^(-Ea2/R/T);
5楼2018-07-04 12:52:59
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

dingd

铁杆木虫 (职业作家)

【答案】应助回帖

★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★
l1003785517: 金币+20, 有帮助 2018-07-04 20:21:40
还是说的不明不白:
1:a3,a4和a5与a1,a2采用相同的处理方式?还是视为3个独立的待求参数?
2:对应4个T有四组数据,是想分别求出这四组数据各自对应的一组参数还是这四组数据拥有一组共同的参数?
3:待求参数有范围限制吗?比如k1,k2必须在[0,1]之间。
7楼2018-07-04 17:52:32
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

月只蓝

主管区长 (职业作家)

引用回帖:
6楼: Originally posted by l1003785517 at 2018-07-04 12:55:03
是的,作为一个求助者真的不应该出现这种错误,很抱歉。
我在代码中给出的公式,Function y1=a1*x1/(1+a3*x1+a4*x2+a5*x3);
                                                    y2=(a2*x2-a1*x1)/(1+a3*x1+a4 ...

建议你求助于dingd,以我所知,他是做拟合的大神。

发自小木虫Android客户端
MATLAB、MS小问题、普通问题请发帖求助!时间精力有限,恕不接受无偿私信求助。
8楼2018-07-04 18:31:35
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

l1003785517

新虫 (小有名气)

引用回帖:
7楼: Originally posted by dingd at 2018-07-04 17:52:32
还是说的不明不白:
1:a3,a4和a5与a1,a2采用相同的处理方式?还是视为3个独立的待求参数?
2:对应4个T有四组数据,是想分别求出这四组数据各自对应的一组参数还是这四组数据拥有一组共同的参数?
3:待求参数 ...

额,只有a1和a2才是这样,a3,a4,a5是三个独立的待求参数,我是希望用这四组数据求出同一组参数来,但是我不知道用1stopt该怎么写,范围限制的话就是所有参数都要是非负,麻烦您了!
9楼2018-07-04 20:14:37
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

l1003785517

新虫 (小有名气)

引用回帖:
8楼: Originally posted by 月只蓝 at 2018-07-04 18:31:35
建议你求助于dingd,以我所知,他是做拟合的大神。
...

好的,非常感谢您的帮助!谢谢
10楼2018-07-04 20:15:45
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 l1003785517 的主题更新
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[硕博家园] 湖北工业大学 生命科学与健康学院-课题组招收2026级食品/生物方向硕士 +3 1喜春8 2026-03-17 5/250 2026-03-17 17:18 by ber川cool子
[考研] 085600材料与化工求调剂 +5 绪幸与子 2026-03-17 5/250 2026-03-17 16:40 by laoshidan
[考研] 26考研求调剂 +6 丶宏Sir 2026-03-13 6/300 2026-03-17 16:13 by 醉在风里
[考研] 材料与化工专硕调剂 +5 heming3743 2026-03-16 5/250 2026-03-17 14:03 by 勇敢太监王公公
[考研] 285化工学硕求调剂(081700) +9 柴郡猫_ 2026-03-12 9/450 2026-03-17 10:18 by Sammy2
[考研] [导师推荐]西南科技大学国防/材料导师推荐 +3 尖角小荷 2026-03-16 6/300 2026-03-16 23:21 by 尖角小荷
[考研] 0854控制工程 359求调剂 可跨专业 +3 626776879 2026-03-14 9/450 2026-03-16 17:42 by 626776879
[考研] 一志愿985,本科211,0817化学工程与技术319求调剂 +5 Liwangman 2026-03-15 5/250 2026-03-16 17:10 by 我的船我的海
[考研] 304求调剂 +5 素年祭语 2026-03-15 5/250 2026-03-16 17:00 by 我的船我的海
[考研] 311求调剂 +5 26研0 2026-03-15 5/250 2026-03-16 16:21 by a不易
[考研] 【0703化学调剂】-一志愿华中师范大学-六级475 +5 Becho359 2026-03-11 5/250 2026-03-14 11:35 by 哦哦123
[考研] 266求调剂 +4 学员97LZgn 2026-03-13 4/200 2026-03-14 08:37 by zhukairuo
[考研] 279求调剂 +3 Dizzy123@ 2026-03-10 3/150 2026-03-13 23:02 by JourneyLucky
[考研] 材料与化工求调剂一志愿 985 总分 295 +8 dream…… 2026-03-12 8/400 2026-03-13 22:17 by 星空星月
[考研] 材料工程调剂 +9 咪咪空空 2026-03-12 9/450 2026-03-13 22:05 by 星空星月
[考研] 26调剂/材料科学与工程/总分295/求收留 +9 2026调剂侠 2026-03-12 9/450 2026-03-13 20:46 by 18595523086
[考研] 材料工程调剂 +4 咪咪空空 2026-03-11 4/200 2026-03-13 19:57 by JourneyLucky
[考研] 求调剂 +5 一定有学上- 2026-03-12 5/250 2026-03-13 18:31 by ms629
[考研] 土木第一志愿276求调剂,科研和技能十分丰富,求新兴方向的导师收留 +3 土木小天才 2026-03-12 3/150 2026-03-13 15:01 by JourneyLucky
[考研] 321求调剂(食品/专硕) +3 xc321 2026-03-12 6/300 2026-03-13 08:45 by xc321
信息提示
请填处理意见