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

24小时热门版块排行榜    

小木虫论坛-学术科研互动平台 » 专业学科区 » 数学 » 计算数学 » 【求助】用matlab最优化方法进行参数拟合
181/2返回列表12下一页
查看: 3133  |  回复: 17
只看楼主@他人 存档 新回复提醒 (忽略) 收藏    在APP中查看
本帖产生 1 个 数学EPI ，点击这里进行查看

wodaifei

银虫 (小有名气)

[交流] 【求助】用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

» 猜你喜欢

» 本主题相关价值贴推荐，对您同样有帮助:
要想拥有一切，就要让自己变得足够优秀！！！
1楼2011-03-21 08:48:28
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
回帖支持 ( 显示支持度最高的前 50 名 )

leedobb

金虫 (正式写手)
★ ★ ★ ★ ★ ★
小木虫(金币+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

金虫 (正式写手)
★
小木虫(金币+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的回帖

leedobb

金虫 (正式写手)
★
小木虫(金币+0.5):给个红包，谢谢回帖交流
怎么算出来反应常数负号，应该也可以吧
赞一下(1人) 回复此楼
有一天，我打了个瞌睡就到了这里，但我知道我掉入了时光的循环中，虽得以永生，但只有第一个循环有意义。
4楼2011-03-21 15:48:53
 已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

wodaifei

银虫 (小有名气)

哥，你太给力了！！！！！！！！

哥，顶死你！！！你太给力了。。。。。虽然动力学常数存在负值这是不对的，但问题不是你的程序，是我的实验数据。小弟跪谢！！！！！！！！！！！
赞一下 回复此楼
要想拥有一切，就要让自己变得足够优秀！！！
5楼2011-03-21 17:48:04
 已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

lvfuxi

新虫 (初入文坛)

太强大了

很厉害，呵呵
回复此楼
6楼2011-03-25 10:34:29
 已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

dingd

铁杆木虫 (职业作家)
★
小木虫(金币+0.5):给个红包，谢谢回帖
这种动力学微分方程拟合，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)*yC;
            yD'=k3*yA+k6*yB+k8*yC;
            yE'=k4*yA+k7*yB+k9*yC;
Data;
//反应时间        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
赞一下(2人) 回复此楼
7楼2011-12-05 10:57:31
 已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

xiaohsu2006

木虫 (著名写手)
★
小木虫(金币+0.5):给个红包，谢谢回帖
引用回帖:
2楼: Originally posted by leedobb at 2011-03-21 15:36:02:
% 哦，你的题目是这个啊，
% 这只能依靠数值解，matlab肯定没有现成的程序吧。
%
% 你这样做吧。构造一个能量函数（其实就是最小二乘法啦）
%
% F =  (yA(t2)- YYA(t2,k1,k2,k3...,k9) ) ^ 2
% ..嗯， ...

你好，我也是类似的问题，我想问一下你这个算法的原理是什么，思路是怎么样的？构造能量函数我明白，目的是让它最小，那么开始是不先给k赋一组初值，那么接下来是怎么做的？
赞一下(1人) 回复此楼
8楼2011-12-10 20:27:56
 已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

xiaohsu2006

木虫 (著名写手)
引用回帖:
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)* ...

你好，你这是用什么拟合的，matlab吗？拟合的算法是根据什么的？请赐教，不胜感激！
赞一下 回复此楼
9楼2011-12-10 20:29:47
 已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

dingd

铁杆木虫 (职业作家)
★
小木虫(金币+0.5):给个红包，谢谢回帖
在网上搜一下1stOpt。
赞一下(1人) 回复此楼
10楼2011-12-10 21:04:27
 已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 wodaifei 的主题更新
181/2返回列表12下一页
普通表情 龙 兔 虎 猫 高级回复(可上传附件) SCI英文润色翻译-EnPapers美团队-先修改后付款
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[微米和纳米] 银浆叠压烧结粘片的问题有没有好的解决方法 +3 微纳材料 2024-05-30 5/250 2024-06-02 21:17 by manmansboy
[教师之家] 选择 +8 459582015 2024-05-31 11/550 2024-06-02 18:30 by auvauv
[硕博家园] 各位同学能否分享一下实验室的学生劳务发放标准呀？ +19 ma3252788 2024-05-30 19/950 2024-06-02 18:15 by Lcy6666
[考博] 求生物学的博导收留 +5 申博成功功成博� 2024-06-01 5/250 2024-06-02 13:19 by 小虫子咔咔
[考博] 导师不让硕转博，让我去国外读博，能理解吗？ +12 萧山幽谷 2024-05-29 20/1000 2024-06-02 12:01 by yuan0806
[考博] 求25博导，金属增材制造方向 +3 22机械 2024-06-01 3/150 2024-06-02 11:17 by Napoleonsky
[硕博家园] 实验室太吵闹，无法安静学习，怎么办？ +7 utahh 2024-05-31 12/600 2024-06-02 01:23 by dking22
[教师之家] 在大地上我们只过一生---看完我的阿勒泰上头了好几天，完结那天晚上几乎失眠 +11 瞬息宇宙 2024-05-27 13/650 2024-06-01 22:15 by otani
[基金申请] 博后特别资助状态变化 +24 随梦而飞2017 2024-05-30 35/1750 2024-06-01 10:10 by 青岛阳仔
[基金申请] 数理的人才答辩通知发了吗 +7 yzy3327 2024-05-30 7/350 2024-06-01 09:15 by 撒哈拉西部
[教师之家] 中年 (金币+3) +18 459582015 2024-05-28 19/950 2024-06-01 00:41 by 沈婉婷.Girl
[论文投稿] 纠结选哪一个期刊，电化学领域 50+8 Freya163 2024-05-28 10/500 2024-05-31 15:09 by wzykobe
[高分子] MMA预聚体光固化发雾问题求助 +3 惠亚金总 2024-05-29 10/500 2024-05-31 14:59 by 惠亚金总
[博后之家] 2024公派博后申请 +4 326lhpqk 2024-05-27 5/250 2024-05-29 20:03 by @古月胡
[论文投稿] 有没有老师需要发表论文 +4 金老师论文助理- 2024-05-29 4/200 2024-05-29 16:51 by liuyupu132
[论文投稿] 核心初审被拒，理由是“选题的意义不明确，文章写得不像是科技论文”，怎么改 5+3 工藤雷花樱 2024-05-27 8/400 2024-05-29 10:09 by topedit
[论文投稿] EI学报，一审返修后，为啥不再送审，直接终审中? +4 qweasd12345 2024-05-27 6/300 2024-05-29 00:02 by dut_ameng
[基金申请] E10开始送了，希望有好运 +5 sail 2024-05-27 5/250 2024-05-28 18:36 by 芝小芝
[有机交流] 奇怪的物质 100+4 桃桃PXS 2024-05-27 7/350 2024-05-28 10:22 by 091602
[硕博家园] 课题做到一半担心做不出来怎么办 +4 mkpooooo 2024-05-27 7/350 2024-05-27 11:56 by kanmiaolucky
信息提示
关闭
请填处理意见
关闭