版块导航: 正在加载中...

登录注册

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

24小时热门版块排行榜

返回列表

本帖产生 1 个数学EPI ，点击这里进行查看

当前只显示满足指定条件的回帖，点击这里查看本话题的所有回帖

wodaifei

银虫 (小有名气)

应助: 0 (幼儿园)
金币: 3339.1
散金: 572
红花: 8
帖子: 248
在线: 351.1小时
虫号: 683543
注册: 2008-12-28
性别: GG
专业: 化工系统工程

[交流] 【求助】用matlab最优化方法进行参数拟合已有7人参与

各位学哥、学姐，大侠们！！
小弟最近遇到一个难题：具体为 dyA/dt=-(k1+k2+k3+k4)yA
                                                      dyB/dt=k1yA-(k5+k6+k7)yB
                                                      dyC/dt=k2yA+k5yB-(k8+k9)yC
                                                      dyD/dt=k3yA+k6yB+k8yC
                                                      dyE/dt=k4yA+k7yB+k9yC
其中y为各组分的质量分率; k为反应速率常数
其中实验数据为：反应时间       A       B       C       D       E
                  0.833    0.0632    0．699    0.187 0.0256    0.0220
                           1.25       0.055       0.722       0.199 0.0305    0.0281
                        2.5       0.0502       0.769       0.201 0.0328       0.031
根据实验数据希望能用最优化法回归出其中的九个参数（最好能用matlab），小弟刚接触这方面的知识，遇到了困难，怎么也估计不出来，希望各位学哥、学姐，大侠们能给小弟指点一下！！跪谢！！！不胜感激!!!!!1真的把小弟难住了。。

回复此楼

» 本帖已获得的红花（最新10朵）

ffjtouchlife

» 猜你喜欢

Cas 72-43-5需要30g，定制合成，能接单的留言已经有8人回复
求助:我三月中下旬出站，青基依托单位怎么办？已经有6人回复
北京211副教授，35岁，想重新出发，去国外做博后，怎么样？已经有8人回复
磺酰氟产物，毕不了业了！已经有5人回复
论文终于录用啦！满足毕业条件了已经有25人回复
2026年机械制造与材料应用国际会议（ICMMMA 2026）已经有3人回复
自荐读博已经有3人回复
不自信的我已经有5人回复
投稿Elsevier的杂志（返修），总是在选择OA和subscription界面被踢皮球已经有8人回复

要想拥有一切，就要让自己变得足够优秀！！！

1楼 2011-03-21 08:48:28

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

xiaohsu2006

木虫 (著名写手)

应助: 23 (小学生)
金币: 2520.4
散金: 207
红花: 4
帖子: 1362
在线: 1137.1小时
虫号: 1012350
注册: 2010-05-06
性别: GG
专业: 多相流热物理学

★
小木虫(金币+0.5):给个红包，谢谢回帖

引用回帖:

7楼: Originally posted by dingd at 2011-12-05 10:57:31:
这种动力学微分方程拟合，1stOpt很简单好使啊！
[code]
Variable T, yA, yB, yC, yD, yE;
ODEFunction yA'=-(k1+k2+k3+k4)*yA;
yB'=k1*yA-(k5+k6+k7)*yB;
yC'=k2*yA+k5*yB-(k8+k9)* ...

t/min y1 y2 y3 y4
2.5 3.8 2.5 0.5 0.5
5 5.8 3.3 1.0 0.6
7.5 5.8 2.7 1.2 0.7
10 5.4 2.3 1.2 0.8
15 5.7 2.4 1.2 1.3
30 8.5 2.6 2.3 2.4
模型方程
Y1’=x*k4*y5*y7+k6*y5+k10*y2*y7-0.1942*k10*y1*y4;
Y2’=x*k3*y5*y7+k5*y5-k10*y2*y7+0.1942*k10*y1*y4-k11*y2*y4+2.7624E-6*k11*y3*y7;
Y3’=k7*y5+k11*y2*y4-2.7624E-6*k11*y3*y7;
Y4’=(x-z+y/2)*k3*y5*y7+(2x-z+y/2)*k4*y5*y7+k8*y5+k10*y2*y7-0.1942*k10*y1*y4-3*k11*y2*y4+8.2872E-6*k11*y3*y7;
Y5’=k2*y7-k3*y5*y7-k4*y5*y7-k9*y5;
Y6’=-k2*y6;
Y7’=-(2x-z)*k4*y5*y7-k10*y2*y7+0.1942*k10*y1*y4+k11*y2*y4-2.7624E-6*k11*y3*y7;
其中，自变量为t，因变量为y1,y2,y3,y4,y5,y6,y7,参数为x,y,z,k2,k3,k4,k5,k6,k7,k8,k9,k10,k11
想回归出其中的模型参数，为什么我的在1stopt上运行不出结果啊，能帮试一下吗？非常感谢！下面是我的程序：
Parameter x,y,z,k2,k3,k4,k5,k6,k7,k8,k9,k10,k11;
SharedModel;
Variable t,y(1:7);
Function Y1'=x*k4*y5*y7+k6*y5+k10*y2*y7-0.1942*k10*y1*y4;
Y2'=x*k3*y5*y7+k5*y5-k10*y2*y7+0.1942*k10*y1*y4-k11*y2*y4+2.7624E-6*k11*y3*y7;
Y3'=k7*y5+k11*y2*y4-2.7624E-6*k11*y3*y7;
Y4'=(x-z+y/2)*k3*y5*y7+(2x-z+y/2)*k4*y5*y7+k8*y5+k10*y2*y7-0.1942*k10*y1*y4-3*k11*y2*y4+8.2872E-6*k11*y3*y7;
Y5'=k2*y7-k3*y5*y7-k4*y5*y7-k9*y5;
Y6'=-k2*y6;
Y7'=-(2x-z)*k4*y5*y7-k10*y2*y7+0.1942*k10*y1*y4+k11*y2*y4-2.7624E-6*k11*y3*y7;
Data;
//t, y1, y2, y3, y4
2.5 3.8 2.5 0.5 0.5
5 5.8 3.3 1.0 0.6
7.5 5.8 2.7 1.2 0.7
10 5.4 2.3 1.2 0.8
15 5.7 2.4 1.2 1.3
30 8.5 2.6 2.3 2.4

[ Last edited by xiaohsu2006 on 2011-12-12 at 20:14 ]

赞一下(1人)

回复此楼

11楼2011-12-12 20:04:19

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

查看全部 18 个回答

leedobb

金虫 (正式写手)

数学EPI: 1
应助: 77 (初中生)
金币: 152.4
散金: 646
红花: 15
帖子: 856
在线: 270.4小时
虫号: 1199535
注册: 2011-02-06
性别: GG
专业: 高分子科学

★ ★ ★ ★ ★ ★
小木虫(金币+0.5):给个红包，谢谢回帖交流
nono2009(金币+5, 数学EPI+1): 鼓励详细应助 2011-03-22 07:48:24

引用回帖:

Originally posted by wodaifei at 2011-03-21 08:48:28:
各位学哥、学姐，大侠们！！
小弟最近遇到一个难题：具体为 dyA/dt=-(k1+k2+k3+k4)yA
dyB/dt=k1yA-(k5+k6+k7)yB
...

% 哦，你的题目是这个啊，
% 这只能依靠数值解，matlab肯定没有现成的程序吧。
%
% 你这样做吧。构造一个能量函数（其实就是最小二乘法啦）
%
% F =  (yA(t2)- YYA(t2,k1,k2,k3...,k9) ) ^ 2
% ..嗯，算我直接编给你吧。

%
% dyA/dt=-(k1+k2+k3+k4)yA
%                                                       dyB/dt=k1yA-(k5+k6+k7)yB
%                                                       dyC/dt=k2yA+k5yB-(k8+k9)yC
%                                                       dyD/dt=k3yA+k6yB+k8yC
%                                                       dyE/dt=k4yA+k7yB+k9yC
% A       B       C       D       E
%                   0.833    0.0632    0．699    0.187 0.0256    0.0220
%                            1.25       0.055       0.722       0.199 0.0305    0.0281
%                            2.5       0.0502       0.769       0.201 0.0328       0.031
clear all

kk0(1:9)=[1.2,0.9,0.8, 1.3,0.6,1.6,  1.7,1.1,1.0];
t1= 1.25  - 0.833;
t2=2.5 -0.833;
y0 = [0.0632    0.699    0.187 0.0256    0.0220]';
y1 = [ 0.055       0.722       0.199 0.0305    0.0281]';
y2 = [ 0.0502       0.769       0.201 0.0328       0.031]';
dk = 0.0001;
dt = 0.1
kk=kk0
for i=1:100000
  %%% 计算 dF/dk
for j=1:9
%% 计算k+dk
k=kk;k(j)=k(j)+dk;
A = [ -(k(1)+k(2)+k(3)+k(4)) 0 0
k(1)  -(k(5)+k(6)+k(7))  0
   k(2)  k(5)  -(k(8)+k(9)) ];
  [V,D]=eig(A);
  DD(1,1)=D(1,1); DD(2,1) = D(2,2); DD(3,1) = D(3,3);
  DD =DD;
  VD=V^-1;
  Z0 = VD*y0(1:3);
  %%% z = z0*exp( DD(i)*t)
  Z1 = Z0 .*exp(DD*t1);
  Y1 = V *Z1;
  Int_Z1 = Z1./DD-Z0./DD; %int (Z1) = Z1/DD+C
  Y1(4) = k(3)*sum(V(1,1:3)*Int_Z1)+k(6)*sum(V(2,1:3)*Int_Z1) + k(8)*sum(V(3,1:3)*Int_Z1) + y0(4);
  Y1(5) = k(4)*sum(V(1,1:3)*Int_Z1)+k(7)*sum(V(2,1:3)*Int_Z1) + k(9)*sum(V(3,1:3)*Int_Z1) + y0(4);
  Z2 = Z0 .*exp(DD*t2);
  Y2 = V *Z2;
  Int_Z2 = Z2./DD-Z0./DD; %int (Z1) = Z1/DD+C
  Y2(4) = k(3)*sum(V(1,1:3)*Int_Z2)+k(6)*sum(V(2,1:3)*Int_Z2) + k(8)*sum(V(3,1:3)*Int_Z2) + y0(4);
  Y2(5) = k(4)*sum(V(1,1:3)*Int_Z2)+k(7)*sum(V(2,1:3)*Int_Z2) + k(9)*sum(V(3,1:3)*Int_Z2) + y0(4);
  F_up = ((Y2-y2)'*(Y2-y2) + (Y1-y1)'*(Y1-y1)) /2;

%% 计算k=k-dk
  k=kk;k(j)=k(j)-dk;
A = [ -(k(1)+k(2)+k(3)+k(4)) 0 0
k(1)  -(k(5)+k(6)+k(7))  0
   k(2)  k(5)  -(k(8)+k(9)) ];
  [V,D]=eig(A);
  DD(1,1)=D(1,1); DD(2,1) = D(2,2); DD(3,1) = D(3,3);
  DD =DD;
  VD=V^-1;
  Z0 = VD*y0(1:3);
  %%% z = z0*exp( DD(i)*t)
  Z1 = Z0 .*exp(DD*t1);
  Y1 = V *Z1;
  Int_Z1 = Z1./DD-Z0./DD; %int (Z1) = Z1/DD+C
  Y1(4) = k(3)*sum(V(1,1:3)*Int_Z1)+k(6)*sum(V(2,1:3)*Int_Z1) + k(8)*sum(V(3,1:3)*Int_Z1) + y0(4);
  Y1(5) = k(4)*sum(V(1,1:3)*Int_Z1)+k(7)*sum(V(2,1:3)*Int_Z1) + k(9)*sum(V(3,1:3)*Int_Z1) + y0(4);
  Z2 = Z0 .*exp(DD*t2);
  Y2 = V *Z2;
  Int_Z2 = Z2./DD-Z0./DD; %int (Z1) = Z1/DD+C
  Y2(4) = k(3)*sum(V(1,1:3)*Int_Z2)+k(6)*sum(V(2,1:3)*Int_Z2) + k(8)*sum(V(3,1:3)*Int_Z2) + y0(4);
  Y2(5) = k(4)*sum(V(1,1:3)*Int_Z2)+k(7)*sum(V(2,1:3)*Int_Z2) + k(9)*sum(V(3,1:3)*Int_Z2) + y0(4);
  F_down = ((Y2-y2)'*(Y2-y2) + (Y1-y1)'*(Y1-y1)) /2;
  FF= (F_down+F_up)/2;
  DF(j)= (F_up-F_down)/dk/2;
end

for j=1:9
kk(j)= kk(j) - DF(j) *dt;
end

if(mod(i,1000)==1)
'deviation'
FF/10
'k1->k9', kk
end







end

[ Last edited by leedobb on 2011-3-21 at 15:46 ]

赞一下(6人)

回复此楼

有一天，我打了个瞌睡就到了这里，但我知道我掉入了时光的循环中，虽得以永生，但只有第一个循环有意义。

2楼2011-03-21 15:36:02

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

leedobb

金虫 (正式写手)

数学EPI: 1
应助: 77 (初中生)
金币: 152.4
散金: 646
红花: 15
帖子: 856
在线: 270.4小时
虫号: 1199535
注册: 2011-02-06
性别: GG
专业: 高分子科学

★
小木虫(金币+0.5):给个红包，谢谢回帖交流

引用回帖:

Originally posted by wodaifei at 2011-03-21 08:48:28:
各位学哥、学姐，大侠们！！
小弟最近遇到一个难题：具体为 dyA/dt=-(k1+k2+k3+k4)yA
dyB/dt=k1yA-(k5+k6+k7)yB
...

这是算出来比较好的结果
k = 0.2974 0.3073 -0.2075 0.0272 0.1690 -0.1148 -0.0798 0.4677 0.2665
此时
Y(t2)= 0.0529 0.7138 0.1869 0.0232 0.0235
y(t3)= 0.0311 0.7527 0.1864 0.0163 0.0138
t1用来当初始条件了。

唉，看起来，我太无聊了。

[ Last edited by leedobb on 2011-3-21 at 15:45 ]

赞一下(1人)

回复此楼

有一天，我打了个瞌睡就到了这里，但我知道我掉入了时光的循环中，虽得以永生，但只有第一个循环有意义。

3楼2011-03-21 15:44:50

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