版塊導(dǎo)航: 正在加載中...

應(yīng)《網(wǎng)絡(luò)安全法》要求，自2017年10月1日起，未進(jìn)行實(shí)名認(rèn)證將不得使用互聯(lián)網(wǎng)跟帖服務(wù)。為保障您的帳號(hào)能夠正常使用，請(qǐng)盡快對(duì)帳號(hào)進(jìn)行手機(jī)號(hào)驗(yàn)證，感謝您的理解與支持！

24小時(shí)熱門版塊排行榜

返回列表

當(dāng)前只顯示滿足指定條件的回帖，點(diǎn)擊這里查看本話題的所有回帖

ronghengyidu

木蟲 (小有名氣)

應(yīng)助: 0 (幼兒園)
金幣: 2854.5
帖子: 218
在線: 62.2小時(shí)
蟲號(hào): 1003427
注冊(cè): 2010-04-22
專業(yè): 生物化工與食品化工

[求助] matlab擬合反應(yīng)動(dòng)力學(xué)參數(shù) 已有1人參與

本人是matlab大白，最近仿照黃華江編著的“實(shí)用化工計(jì)算機(jī)模擬”例題7-5擬合動(dòng)力學(xué)參數(shù)，
         動(dòng)力學(xué)方程為dc1/dt=k1*c1/(k2+k3*c1+k4*c2+k5*c1*c2)，擬合k1,k2,k3,k4,k5,
         matlab程序和動(dòng)力學(xué)數(shù)據(jù)見附件，運(yùn)行出現(xiàn)如下錯(cuò)誤，請(qǐng)各位大拿幫助。


Error: File: KineticsEst.m Line: 47 Column: 60
Unbalanced or unexpected parenthesis or bracket.

Error in run (line 64)
evalin('caller', [script ';']);

Error: File: KineticsEst.m Line: 47 Column: 60
Unbalanced or unexpected parenthesis or bracket.

Error in run (line 64)
evalin('caller', [script ';']);

回復(fù)此樓

» 本帖附件資源列表

歡迎監(jiān)督和反饋：小木蟲僅提供交流平臺(tái)，不對(duì)該內(nèi)容負(fù)責(zé)。
本內(nèi)容由用戶自主發(fā)布，如果其內(nèi)容涉及到知識(shí)產(chǎn)權(quán)問題，其責(zé)任在于用戶本人，如對(duì)版權(quán)有異議，請(qǐng)聯(lián)系郵箱：xiaomuchong@tal.com
附件 1 : KineticsEst.m

2016-10-03 16:39:52, 1.98 K

附件 2 : KineticsData.m

2016-10-03 16:40:16, 365 bytes

1樓 2016-10-03 15:44:36

已閱回復(fù)此樓關(guān)注TA 給TA發(fā)消息送TA紅花 TA的回帖

lu_yu_lan

新蟲 (初入文壇)

應(yīng)助: 3 (幼兒園)
金幣: 101
散金: 5
紅花: 1
帖子: 43
在線: 26.4小時(shí)
蟲號(hào): 4935903
注冊(cè): 2016-08-22

【答案】應(yīng)助回帖

引用回帖:

5樓: Originally posted by ronghengyidu at 2016-10-04 19:30:02
C2初始值是C1初始值的十倍，但是反應(yīng)過程中C1消耗值與C2消耗值一樣...

根據(jù)樓主的意思，微分方程及擬合參數(shù)似乎可以簡(jiǎn)化為：dc1/dt=c1/(k1+k2*c1+k3*c1^2)
相應(yīng)地，數(shù)據(jù)也去掉1500那一列。
OpenLu代碼：

CODE:

!!!using["IMSL","luopt","math"]; //使用命名空間

f(t,c1,dc1:: k1,k2,k3)=

{

   dc1=c1/(k1+k2*c1+k3*c1^2)

};

目標(biāo)函數(shù)(_k1,_k2,_k3 : tyz: tyArray,tA,k1,k2,k3)=

{

   k1=_k1, k2=_k2, k3=_k3, //傳遞優(yōu)化變量，函數(shù)f中要用到k1,k2,k3

   tyz=ode[@f,tA,ra1(150)],

   sum[(tyArray-tyz).^2.0,0]

};

main(::tyArray,tA)=

{

       tyArray=matrix{   //存放實(shí)驗(yàn)數(shù)據(jù)

"

     0     150                

     10    137.3        

     30    107.7        

     60    72.8       

     90    48.8         

     120   22.9        

     150   11.8         

     180   5.4          

     220   1.9       

     240   0.9      

"

    },

    tA=tyArray(all:0), //tA取矩陣的列

    ClearImslErr(),    //清空IMSL錯(cuò)誤輸出

    ERSET(0,0,0),      //關(guān)閉IMSL所有警告

    Opt[@目標(biāo)函數(shù)],    //Opt函數(shù)全局優(yōu)化

    ERSET(0,2,2), ERSET(0,1,0)   //恢復(fù)IMSL警告

};

結(jié)果（k1,k2,k3，目標(biāo)函數(shù)值）：
-28.81117789336092 -0.6456997498853914 1.11232432991073e-003 19.60179363209554

圖形：

CODE:

!!!using("IMSL","win","math");

f(t,c1,dc1:: k1,k2,k3)=

{

   dc1=c1/(k1+k2*c1+k3*c1^2)

};

init(x,t::k1,k2,k3)=

  x = matrix[

"

     0     150                

     10    137.3        

     30    107.7        

     60    72.8       

     90    48.8         

     120   22.9        

     150   11.8         

     180   5.4          

     220   1.9       

     240   0.9 

"

  ],

  t=x(all:0).reshape(),

  k1=-28.81117789336092, k2=-0.6456997498853914, k3=1.11232432991073e-003,

  luShareX2(x, ode[@f,t,ra1(150)]);

ChartWnd[@init];

tuxing.png

贊一下

回復(fù)此樓

6樓2016-10-04 22:50:57

已閱回復(fù)此樓關(guān)注TA 給TA發(fā)消息送TA紅花 TA的回帖

查看全部 8 個(gè)回答

lu_yu_lan

新蟲 (初入文壇)

應(yīng)助: 3 (幼兒園)
金幣: 101
散金: 5
紅花: 1
帖子: 43
在線: 26.4小時(shí)
蟲號(hào): 4935903
注冊(cè): 2016-08-22

【答案】應(yīng)助回帖

★ ★
感謝參與，應(yīng)助指數(shù) +1
月只藍(lán): 金幣+2, 感謝應(yīng)助！ 2016-10-04 21:59:50

樓主的問題，沒有看明白，找到了黃華江編著的“實(shí)用化工計(jì)算機(jī)模擬”例題7-5擬合動(dòng)力學(xué)參數(shù)

7個(gè)微分方程，5個(gè)擬合參數(shù)k1,k2,k3,k4,k5，初值t=0時(shí)，x為0.1883,  0.2507,  0.0467,  0.0899,  0.1804,  0.1394,  0.1046

  q = 8.75 + k5,
  dx1/t=  k5-q*x1- k1*x1*x2-k4*x1*x6*sqrt(0.9),
  dx2/t=  7.0-q*x2 - k1*x1*x2-2*k2*x2*x3,
  dx3/t=  1.75 -q*x3 - k2*x2*x3,
  dx4/t=  -q*x4 + 2*k1*x1*x2-k3*x4*x5,
  dx5/t=  -q*x5 + 3*k2*x2*x3-k3*x4*x5,
  dx6/t=  -q*x6 + 2*k3*x4*x5-k4*x1*x6*sqrt(0.9),
  dx7/t=  -q*x7 + 2*k4*x1*x6*sqrt(0.9)

數(shù)據(jù)：t,x1,x4,x5,x6
0    0.1883 0.0899 0.1804    0.1394
0.0100 0.2047 0.0866 0.1729 0.1297
0.0200 0.2181 0.0856 0.1680 0.1205
0.0300 0.2291 0.0863 0.1647 0.1123
0.0400 0.2382 0.0878 0.1623 0.1053
0.0500 0.2459 0.0899 0.1604 0.0995
0.0600 0.2523 0.0921 0.1588 0.0948
0.0700 0.2576 0.0945 0.1574 0.0911
0.0800 0.2622 0.0968 0.1561 0.0882
0.0900 0.2660 0.0989 0.1548 0.0859
0.1000 0.2692 0.1010 0.1537 0.0842
0.1100 0.2719 0.1028 0.1525 0.0830
0.1200 0.2742 0.1045 0.1515 0.0821
0.1300 0.2761 0.1060 0.1505 0.0814
0.1400 0.2777 0.1074 0.1495 0.0810
0.1500 0.2790 0.1086 0.1487 0.0807
0.1600 0.2801 0.1096 0.1479 0.0805
0.1700 0.2811 0.1106 0.1471 0.0803
0.1800 0.2819 0.1114 0.1465 0.0803
0.1900 0.2825 0.1121 0.1458 0.0803
0.2000 0.2830 0.1127 0.1453 0.0803

用OpenLu求解：

CODE:

!!!using["IMSL","luopt","math"]; //使用命名空間

f(t,x1,x2,x3,x4,x5,x6,x7,dx1,dx2,dx3,dx4,dx5,dx6,dx7 : q : k1,k2,k3,k4,k5)=

{

  q = 8.75 + k5,

  dx1=  k5-q*x1- k1*x1*x2-k4*x1*x6*sqrt(0.9),

  dx2=  7.0-q*x2 - k1*x1*x2-2*k2*x2*x3,

  dx3=  1.75 -q*x3 - k2*x2*x3,

  dx4=  -q*x4 + 2*k1*x1*x2-k3*x4*x5,

  dx5=  -q*x5 + 3*k2*x2*x3-k3*x4*x5,

  dx6=  -q*x6 + 2*k3*x4*x5-k4*x1*x6*sqrt(0.9),

  dx7=  -q*x7 + 2*k4*x1*x6*sqrt(0.9)

};

目標(biāo)函數(shù)(_k1,_k2,_k3,_k4,_k5 : i,s,tyz,tyy : tyArray,tA,k1,k2,k3,k4,k5)=

{

  k1=_k1, k2=_k2, k3=_k3, k4=_k4, k5=_k5,   //傳遞優(yōu)化變量，函數(shù)f中要用到k1,k2,k3,k4,k5

  tyz=ode[@f,tA,ra1(0.1883,  0.2507,  0.0467,  0.0899,  0.1804,  0.1394,  0.1046)],

  tyy=tyArray(zeros), tyy(all:0,1)=subg(tyz,all:0,1), tyy(all:2,4)=subg(tyz,all:4,6),

  sum[(tyy-tyArray).^2.0,0]

};

main(::tyArray,tA)=

{

      tyArray=matrix{   //存放實(shí)驗(yàn)數(shù)據(jù)t,x1,x4,x5,x6

"

    0                0.1883        0.0899    0.1804          0.1394

    0.0100    0.2047    0.0866    0.1729    0.1297

    0.0200    0.2181    0.0856    0.1680    0.1205

    0.0300    0.2291    0.0863    0.1647    0.1123

    0.0400    0.2382    0.0878    0.1623    0.1053

    0.0500    0.2459    0.0899    0.1604    0.0995

    0.0600    0.2523    0.0921    0.1588    0.0948

    0.0700    0.2576    0.0945    0.1574    0.0911

    0.0800    0.2622    0.0968    0.1561    0.0882

    0.0900    0.2660    0.0989    0.1548    0.0859

    0.1000    0.2692    0.1010    0.1537    0.0842

    0.1100    0.2719    0.1028    0.1525    0.0830

    0.1200    0.2742    0.1045    0.1515    0.0821

    0.1300    0.2761    0.1060    0.1505    0.0814

    0.1400    0.2777    0.1074    0.1495    0.0810

    0.1500    0.2790    0.1086    0.1487    0.0807

    0.1600    0.2801    0.1096    0.1479    0.0805

    0.1700    0.2811    0.1106    0.1471    0.0803

    0.1800    0.2819    0.1114    0.1465    0.0803

    0.1900    0.2825    0.1121    0.1458    0.0803

    0.2000    0.2830    0.1127    0.1453    0.0803

"

    },

    tA=tyArray(all:0), //tA取矩陣的列

    ClearImslErr(),    //清空IMSL錯(cuò)誤輸出

    ERSET(0,0,0),      //關(guān)閉IMSL所有警告

    Opt[@目標(biāo)函數(shù)],    //Opt函數(shù)全局優(yōu)化

    ERSET(0,2,2), ERSET(0,1,0)   //恢復(fù)IMSL警告

};

結(jié)果（k1,k2,k3,k4,k5，目標(biāo)函數(shù)值）：
17.60849027475542 73.06467445464085 51.32551162517818 23.02494410245102 6.001239504989786 6.794401934343191e-008

優(yōu)于原書matlab結(jié)果：
17.49, 72.2996, 50.9086, 22.5339, 5.9799， 3.614535245714516e-007

贊一下(2人)

回復(fù)此樓

» 本帖已獲得的紅花（最新10朵）

ronghengyidu

2樓2016-10-03 21:21:39

已閱回復(fù)此樓關(guān)注TA 給TA發(fā)消息送TA紅花 TA的回帖

ronghengyidu

木蟲 (小有名氣)

應(yīng)助: 0 (幼兒園)
金幣: 2854.5
帖子: 218
在線: 62.2小時(shí)
蟲號(hào): 1003427
注冊(cè): 2010-04-22
專業(yè): 生物化工與食品化工

送紅花一朵

引用回帖:

2樓: Originally posted by lu_yu_lan at 2016-10-03 21:21:39
樓主的問題，沒有看明白，找到了黃華江編著的“實(shí)用化工計(jì)算機(jī)模擬”例題7-5擬合動(dòng)力學(xué)參數(shù)

7個(gè)微分方程，5個(gè)擬合參數(shù)k1,k2,k3,k4,k5，初值t=0時(shí)，x為0.1883, 0.2507, 0.0467, 0.0899, 0.1804, 0.1394, 0 ...

非常感謝您的回復(fù)，我只有一個(gè)微分方程，dc1/dt=k1*c1/(k2++k3*c1+k4*c2+k5*c1*c2), 擬合5個(gè)參數(shù)k1,k2,k3,k4,k5。我的t，c1，c2數(shù)據(jù)在附件中，matlab程序也在附件中，這兩個(gè)文件都可以用記事本打開。在這里我把程序和實(shí)驗(yàn)數(shù)據(jù)再貼出來。

function KineticsEst
% 動(dòng)力學(xué)ODE方程模型的參數(shù)估計(jì)

clear all
clc

k0 = [0.5  0.5  0.5  0.5  0.5];       % 參數(shù)初值
lb = [0  0  0  0  0];                % 參數(shù)下限
ub = [+inf  +inf  +inf  +inf  +inf]; % 參數(shù)上限
x0 = [150  1500];
KineticsData;
yexp = ExpData(:,2:3);                % yexp: 實(shí)驗(yàn)數(shù)據(jù)[x1 x2]

% 使用函數(shù)fmincon()進(jìn)行參數(shù)估計(jì)
[k,fval,flag] = fmincon(@ObjFunc4Fmincon,k0,[],[],[],[],lb,ub,[],[],yexp);
fprintf('\n使用函數(shù)fmincon()估計(jì)得到的參數(shù)值為:\n')
fprintf('\tk1 = %.4f\n',k(1))
fprintf('\tk2 = %.4f\n',k(2))
fprintf('\tk3 = %.4f\n',k(3))
fprintf('\tk4 = %.4f\n',k(4))
fprintf('\tk5 = %.4f\n',k(5))
fprintf('  The sum of the squares is: %.1e\n\n',fval)
k_fmincon = k;

% 使用函數(shù)lsqnonlin()進(jìn)行參數(shù)估計(jì)
[k,resnorm,residual,exitflag,output,lambda,jacobian] = ...
lsqnonlin(@ObjFunc4LNL,k0,lb,ub,[],x0,yexp);
ci = nlparci(k,residual,jacobian);
fprintf('\n\n使用函數(shù)lsqnonlin()估計(jì)得到的參數(shù)值為:\n')
Output

% 以函數(shù)fmincon()估計(jì)得到的結(jié)果為初值，使用函數(shù)lsqnonlin()進(jìn)行參數(shù)估計(jì)
k0 = k_fmincon;
[k,resnorm,residual,exitflag,output,lambda,jacobian] = ...
lsqnonlin(@ObjFunc4LNL,k0,lb,ub,[],x0,yexp);
ci = nlparci(k,residual,jacobian);
fprintf('\n\n以fmincon()的結(jié)果為初值，使用函數(shù)lsqnonlin()估計(jì)得到的參數(shù)值為:\n')
Output

% ------------------------------------------------------------------
function f = ObjFunc4Fmincon(k,x0,yexp)
tspan = ExpData(:,1);
[t x] = ode45(@KineticEqs,tspan,x0,[],k);
y(:,1:2) = x(:,1:2);
f = sum((y(:,1)-yexp(:,1)).^2) + sum((y(:,2)-yexp(:,2)).^2));

% ------------------------------------------------------------------
function f = ObjFunc4LNL(k,x0,yexp)
tspan = = ExpData(:,1);
[t x] = ode45(@KineticEqs,tspan,x0,[],k);
y(:,1:2) = x(:,1:2);
f1 = y(:,1) - yexp(:,1);
f2 = y(:,2) - yexp(:,2);
f = [f1; f2];

% ------------------------------------------------------------------
function dxdt = KineticEqs(t,x,k)
dxdt =  ...
k(1)*x(1)./(k(2)+k(3)*x(1)+k(4)*x(2)+k(5)*x(1).*x(2));

實(shí)驗(yàn)數(shù)據(jù)如下：

動(dòng)力學(xué)數(shù)據(jù):
t    c1       c2
0    150       1500
10 137.3 1487.3
30 107.7 1457.7
60 72.8    1422.8
90 48.8    1398.8
120 22.9    1372.9
150 11.8    1361.8
180 5.4    1355.4
220 1.9    1351.9
240 0.9    1350.9

期望樓上能夠指點(diǎn)，如果您能夠給出openLU程序，我也非常感謝，因?yàn)槲业谝淮温犝f這個(gè)程序，怕仿照編程又要出錯(cuò)。非常感謝您在節(jié)假日能夠回答我的問題。

贊一下

回復(fù)此樓

3樓2016-10-04 05:01:01

已閱回復(fù)此樓關(guān)注TA 給TA發(fā)消息送TA紅花 TA的回帖

lu_yu_lan

新蟲 (初入文壇)

應(yīng)助: 3 (幼兒園)
金幣: 101
散金: 5
紅花: 1
帖子: 43
在線: 26.4小時(shí)
蟲號(hào): 4935903
注冊(cè): 2016-08-22

引用回帖:

3樓: Originally posted by ronghengyidu at 2016-10-04 05:01:01
非常感謝您的回復(fù)，我只有一個(gè)微分方程，dc1/dt=k1*c1/(k2++k3*c1+k4*c2+k5*c1*c2), 擬合5個(gè)參數(shù)k1,k2,k3,k4,k5。我的t，c1，c2數(shù)據(jù)在附件中，matlab程序也在附件中，這兩個(gè)文件都可以用記事本打開。在這里我把程序 ...

是不是還有一個(gè)微分方程，dc2/dt=... ...。
否則應(yīng)該怎樣使用c2數(shù)據(jù)呢？

贊一下

回復(fù)此樓

» 本帖已獲得的紅花（最新10朵）

ronghengyidu

4樓2016-10-04 08:37:24

已閱回復(fù)此樓關(guān)注TA 給TA發(fā)消息送TA紅花 TA的回帖

查看全部 8 個(gè)回答

最具人氣熱帖推薦 [查看全部]		作者	回/看	最后發(fā)表

[考研] 26考研一志愿中國石油大學(xué)(華東)305分求調(diào)劑 +6	嘉年新程 2026-03-15	6/300	2026-03-21 17:07 by Dream007008
[考研] 【考研調(diào)劑】化學(xué)專業(yè) 281分，一志愿四川大學(xué)，誠心求調(diào)劑 +10	吃吃吃才有意義 2026-03-19	10/500	2026-03-21 16:43 by macy2011
[考研] 一志愿重慶大學(xué)085700資源與環(huán)境總分308求調(diào)劑 +7	墨墨漠 2026-03-20	7/350	2026-03-21 16:36 by barlinike
[考研] 317求調(diào)劑 +8	申子申申 2026-03-19	14/700	2026-03-21 16:29 by 我的船我的海
[考研] 求助 +4	夢(mèng)里的無言 2026-03-21	5/250	2026-03-21 15:53 by ColorlessPI
[考研] 268求調(diào)劑 +9	簡(jiǎn)單點(diǎn)0 2026-03-17	9/450	2026-03-21 15:37 by lature00
[考研] 070300化學(xué)319求調(diào)劑 +7	錦鯉0909 2026-03-17	7/350	2026-03-21 03:46 by JourneyLucky
[考研] 初始318分求調(diào)劑（有工作經(jīng)驗(yàn)） +3	1911236844 2026-03-17	3/150	2026-03-21 02:33 by JourneyLucky
[考研] 311求調(diào)劑 +5	冬十三 2026-03-18	5/250	2026-03-21 00:16 by JourneyLucky
[考研] 317求調(diào)劑 +5	申子申申 2026-03-19	9/450	2026-03-20 22:26 by JourneyLucky
[考研] 0817 化學(xué)工程 299分求調(diào)劑有科研經(jīng)歷有二區(qū)文章 +22	rare12345 2026-03-18	22/1100	2026-03-20 20:39 by zhukairuo
[考研] 289求調(diào)劑 +6	懷瑾握瑜l 2026-03-20	6/300	2026-03-20 20:30 by 學(xué)員8dgXkO
[考研] 材料學(xué)求調(diào)劑 +4	Stella_Yao 2026-03-20	4/200	2026-03-20 20:28 by ms629
[考研] 0703化學(xué)調(diào)劑，六級(jí)已過，有科研經(jīng)歷 +13	曦熙兮 2026-03-15	13/650	2026-03-20 19:35 by Dream007008
[考研] 廣西大學(xué)家禽遺傳育種課題組2026年碩士招生（接收計(jì)算機(jī)專業(yè)調(diào)劑） +3	123阿標(biāo) 2026-03-17	3/150	2026-03-20 15:58 by 飛行琦
[考研] 一志愿中國海洋大學(xué)，生物學(xué)，301分，求調(diào)劑 +5	1孫悟空 2026-03-17	6/300	2026-03-19 23:46 by zcl123
[考研] 材料考研調(diào)劑 +3	xwt。 2026-03-19	3/150	2026-03-19 11:22 by w沐陽w
[考研] 生物學(xué)071000 329分求調(diào)劑 +3	我愛生物生物愛?/a> 2026-03-17	3/150	2026-03-18 10:12 by macy2011
[考研] 0703化學(xué)336分求調(diào)劑 +6	zbzihdhd 2026-03-15	7/350	2026-03-18 09:53 by zhukairuo
[考研] [導(dǎo)師推薦]西南科技大學(xué)國防/材料導(dǎo)師推薦 +3	尖角小荷 2026-03-16	6/300	2026-03-16 23:21 by 尖角小荷

亭亭五月天在线观看,亭亭五月天在线观看,国产最新av一区二区,国产 高清 中文字幕,99re热久久亚洲综合精品成人,熟妇 一区二区三区,一级做a爰片性色毛片武则天,美女的骚穴视频播放,国产美女午夜免费视频

24小時(shí)熱門版塊排行榜

ronghengyidu

[求助] matlab擬合反應(yīng)動(dòng)力學(xué)參數(shù) 已有1人參與

» 本帖附件資源列表

» 收錄本帖的淘帖專輯推薦

» 猜你喜歡

» 本主題相關(guān)價(jià)值貼推薦，對(duì)您同樣有幫助:

lu_yu_lan

【答案】應(yīng)助回帖

lu_yu_lan

【答案】應(yīng)助回帖

» 本帖已獲得的紅花（最新10朵）

ronghengyidu

lu_yu_lan

» 本帖已獲得的紅花（最新10朵）

亭亭五月天在线观看,亭亭五月天在线观看,国产最新av一区二区,国产高清中文字幕,99re热久久亚洲综合精品成人,熟妇一区二区三区,一级做a爰片性色毛片武则天,美女的骚穴视频播放,国产美女午夜免费视频