版塊導航: 正在加載中...

調(diào)劑小程序

登錄注冊

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

24小時熱門版塊排行榜

>論壇更新日志 (4593)
>考研 (787)
>導師招生 (658)
>蟲友互識 (614)
>文獻求助 (149)
>休閑灌水 (145)
>考博 (97)
>碩博家園 (81)
>招聘信息布告欄 (63)
>基金申請 (55)
>找工作 (42)
>論文投稿 (38)
>公派出國 (37)
>論文道賀祈福 (30)
>教師之家 (29)
>博后之家 (27)

返回列表

_romantic_鎮(zhèn)

木蟲 (著名寫手)

應(yīng)助: 2 (幼兒園)
金幣: 3339
散金: 387
紅花: 6
帖子: 1221
在線: 677小時
蟲號: 3186948
注冊: 2014-05-07
專業(yè): 化學工程及工業(yè)化學

[求助] 初學者求助，matlab模仿程序已有1人參與

function KineticsEst5
% 動力學ODE方程模型的參數(shù)估計

clear all
clc

k0 = [0  0  0  0];       % 參數(shù)初值
lb = [ -inf  -inf  -inf  -inf ];                % 參數(shù)下限
ub = [+inf  +inf  +inf  +inf ]; % 參數(shù)上限
x0 = [8.5 28.8 27.6];
tspan=[0 0.222 0.333 0.444];
yexp =  [8.5000 28.8000 27.6000
4.8000 23.2000 35.3000
4.2000 21.6000 36.5000
4.0000 21.2000 37.3000];                % yexp: 實驗數(shù)據(jù)[x1 x2 x3]

% 使用函數(shù)fmincon()進行參數(shù)估計
[k,fval,flag] = fmincon(@ObjFunc4Fmincon,k0,[],[],[],[],lb,ub,[],[],x0,yexp);
fprintf('\n使用函數(shù)fmincon()估計得到的參數(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('  The sum of the squares is: %.1e\n\n',fval)
k_fmincon = k;

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

% 以函數(shù)fmincon()估計得到的結(jié)果為初值，使用函數(shù)lsqnonlin()進行參數(shù)估計
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()估計得到的參數(shù)值為:\n')
Output

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

% 以函數(shù)fmincon()估計得到的結(jié)果為初值，使用函數(shù)lsqnonlin()進行參數(shù)估計
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()估計得到的參數(shù)值為:\n')
Output

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

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

回復此樓

» 猜你喜歡

考研交流已經(jīng)有0人回復
【新加坡】納米電子器件項目組有“聯(lián)合培養(yǎng)博士生”名額已經(jīng)有0人回復
物理學I論文潤色/翻譯怎么收費? 已經(jīng)有262人回復
【新加坡】納米電子器件項目組有“聯(lián)合培養(yǎng)博士生”名額已經(jīng)有0人回復
幫我的英語口語老師找學生已經(jīng)有3人回復
計算三階力常數(shù)用thirdorder擴胞時參數(shù)如何選擇？臨近原子數(shù)怎么取值已經(jīng)有5人回復
2026年循環(huán)經(jīng)濟功能材料國際會議（ICFMCE 2026）已經(jīng)有0人回復
2026年第五屆電氣、電子與信息工程國際會議（ISEEIE 2026）已經(jīng)有0人回復

1樓 2021-12-20 20:54:40

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

hzlhm

至尊木蟲 (著名寫手)

應(yīng)助: 387 (碩士)
金幣: 17807
紅花: 53
帖子: 2879
在線: 606.3小時
蟲號: 1517335
注冊: 2011-11-30
性別: GG
專業(yè): 常微分方程與動力系統(tǒng)

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

★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★
_romantic_鎮(zhèn): 金幣+50, ★有幫助 2022-03-15 11:07:58

程序中的問題：
1、tspan = [0.00 : 0.01 : 0.44]與tspan=[0 0.222 0.333 0.444]沒有對應(yīng)
2、f = sum((y(:,1)-yexp(:,1)).^2) + sum((y(:,2)-yexp(:,2)).^2) + sum((y(:,3)-yexp(:,3)).^2) 后面的代碼是多余的，計算的結(jié)果不使用
執(zhí)行結(jié)果：

使用函數(shù)fmincon()估計得到的參數(shù)值為:
k1 = -0.6403
k2 = 1.8125
k3 = 3.6968
k4 = 1.9950
The sum of the squares is: 4.0

贊一下(1人)

回復此樓

QQ：2120156492

2樓2021-12-27 22:26:21

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

相關(guān)版塊跳轉(zhuǎn) 我要訂閱樓主 _romantic_鎮(zhèn) 的主題更新