| 2 | 1/1 | 返回列表 |
| 查看: 729 | 回復(fù): 1 | |||
| 【懸賞金幣】回答本帖問題,作者szy242424將贈送您 5 個金幣 | |||
szy242424新蟲 (初入文壇)
|
[求助]
matalb軟件 ode23s,非線性最小二乘法模擬動力學(xué)參數(shù) 已有1人參與
|
||
|
function kk1 k0=[236,18,45.6,10,1.2,0.001,0.1,0.001,0.01,0.001,0.1,0.1]; lb=[236,18,45.6,10, 1.2, 0.001, 0.1,0.001,0.01,0.001,0.1,0.1]; ub=[inf,inf,inf,inf, 10.4,1.2, 2.4,1.2,0.8,0.1,4.5,4.5]; data=... [0 0.562205 0 0 34.2775 6 0.618817941 0 0 33.9805 12 0.797936454 0 0 31.941 18 1.554008384 0.9935 0 28.739 24 2.141789344 1.3815 0 26.3835 30 2.543955264 1.5745 1.4795 23.8955 36 3.017273616 1.908 1.7625 21.334 42 3.295696176 2.9885 2.038 19.128 60 3.274041088 4.693 3.262 11.2065 66 3.4070652 5.1295 3.581 8.9005 72 3.753546608 5.6395 4.041 6.4395 84 3.595773824 6.6735 4.71 1.496 ]; x0=data(1,2:end); tspan =[0,6,12,18,24,30,36,42,60,66,72,84]; yexp = [data(2:end,2) data(2:end,3) data(2:end,4) data(2:end,5) ]; ts=data(1:end,1); [k,resnorm,residual,exitflag,output,lambda,jacobian] =lsqnonlin(@ObjFunc,k0,lb,ub,[],tspan,x0,yexp); ci = nlparci(k,residual,jacobian); fprintf('\n\n使用函數(shù)lsqnonlin()估計得到的參數(shù)值為:\n') fprintf('\tk1 = %.9f ± %.9f\n',k(1),ci(1,2)-k(1)) fprintf('\tk2 = %.9f ± %.9f\n',k(2),ci(2,2)-k(2)) fprintf('\tk3 = %.9f ± %.9f\n',k(3),ci(3,2)-k(3)) fprintf('\tk4 = %.9f ± %.9f\n',k(4),ci(4,2)-k(4)) fprintf('\tk5 = %.9f ± %.9f\n',k(5),ci(5,2)-k(5)) fprintf('\tk6 = %.9f ± %.9f\n',k(6),ci(6,2)-k(6)) fprintf('\tk7 = %.9f ± %.9f\n',k(7),ci(7,2)-k(7)) fprintf('\tk8= %.9f ± %.9f\n',k(8),ci(8,2)-k(8)) fprintf('\tk9 = %.9f ± %.9f\n',k(9),ci(9,2)-k(9)) fprintf('\tk10 = %.9f ± %.9f\n',k(10),ci(10,2)-k(10)) fprintf('\tk11= %.9f ± %.9f\n',k(11),ci(11,2)-k(11)) fprintf('\tk12 = %.9f ± %.9f\n',k(12),ci(12,2)-k(12)) fprintf(' The sum of the squares is: %.9e\n\n',resnorm) [ts, ys] = ode23s(@KineticsEqs,ts,x0,[],k); yy = [data(:,2) data(:,3) data(:,4) data(:,5) ]; plot(ts,ys(:,1),'b',tspan,yy(:,1),'bo'); hold on plot(ts,ys(:,2),'r',tspan,yy(:,2),'r*'); plot(ts,ys(:,3),'k',tspan,yy(:,3),'k+') legend('DCW的計算值','DCW的實驗值','SA的計算值','SA的實驗值','AA的計算值','AA的實驗值','X的計算值','X的實驗值') function f = ObjFunc(k,tspan,x0,yexp) % 目標(biāo)函數(shù) [~, Xsim] = ode23s(@KineticsEqs,tspan,x0,[],k); Xsim1=Xsim(:,1);%提取Xsim的第一列 Xsim2=Xsim(:,2);%提取Xsim的第二列 Xsim3=Xsim(:,3);%提取Xsim的第三列 Xsim4=Xsim(:,4);%提取Xsim的第四列 ysim(:,1) = Xsim1(2:end);%微分方程的第一個變量,從第二個解到最后一個解賦值給ysim的第一列 ysim(:,2) = Xsim2(2:end);%微分方程的第二個變量,從第二個解到最后一個解賦值給ysim的第二列 ysim(:,3) = Xsim3(2:end);%微分方程的第三個變量,從第二個解到最后一個解賦值給ysim的第三列 ysim(:,4) = Xsim4(2:end);%微分方程的第四個變量,從第二個解到最后一個解賦值給ysim的第四列 f = [(ysim(:,1)-yexp(:,1)) (ysim(:,2)-yexp(:,2)) (ysim(:,3)-yexp(:,3)) (ysim(:,4)-yexp(:,4))]; function dCdt = KineticsEqs(~,C,k) % ODE模型方程,C1、C2、C3、C4分別為底物濃度、產(chǎn)物丁二酸濃度、產(chǎn)物乙酸濃度、DCW值,t為導(dǎo)數(shù) dC1dt =(0.1933862*C(4)/(C(4)+0.78+C(4)^2/296.9)*(1-C(2)/k(1))^k(2)*(1-C(3)/k(3))^k(4))*C(1); dC2dt =k(5)*dC1dt+k(6)*C(1);%第二個微分方程,等號前面是C(2)對t的微分 dC3dt =k(7)*dC1dt+k(8)*C(1);%第三個微分方程,等號前面是C(3)對t的微分 dC4dt =-((1/k(9))*dC1dt+k(10)*C(1)+(1/k(11))*dC2dt+(1/k(12))*dC3dt); dCdt = [dC1dt; dC2dt;dC3dt;dC4dt];%微分方程組. 各位大神,幫我看一下這個代碼有什么錯誤沒有,目標(biāo)是用ode23s求解微分方程組,然后用非線性最小二乘模擬參數(shù),急求回復(fù),謝謝大家! |
鐵桿木蟲 (職業(yè)作家)
|
參考下面1stOpt計算的結(jié)果: Root of Mean Square Error (RMSE): 0.673324275669162 Sum of Squared Residual: 19.9480855290377 Correlation Coef. (R): 0.983872215020488 R-Square: 0.968004535489321 Parameter Best Estimate -------------------- ------------- k1 2255.98568602798 k2 2006.16856712862 k3 123314855.293705 k4 11.7085151403374 k5 1.20000000001939 k6 0.00942914835572298 k7 0.100000000023396 k8 0.0201718372859985 k9 0.799999999987507 k10 0.0999392603337909 k11 4.49999420593383 k12 0.75747434507548 |
| 2 | 1/1 | 返回列表 |
| 最具人氣熱帖推薦 [查看全部] | 作者 | 回/看 | 最后發(fā)表 | |
|---|---|---|---|---|
|
[考研] 337一志愿華南理工0805材料求調(diào)劑 +4 | mysdl 2026-03-11 | 6/300 |
|
|---|---|---|---|---|
|
[考研] 材料與化工085600調(diào)劑求老師收留 +5 | jiaanl 2026-03-11 | 5/250 |
|
|
[考研] 求調(diào)劑 +5 | yfihxh 2026-03-09 | 5/250 |
|
|
[考研] 求調(diào)劑 +4 | 鶴遨予卿 2026-03-09 | 4/200 |
|
|
[考研] 420求調(diào)劑 +3 | 莫向外求11 2026-03-10 | 3/150 |
|
|
[考研] 327求調(diào)劑 +3 | Ffff03 2026-03-10 | 3/150 |
|
|
[考博] 26申博求助 +3 | 跳躍餅干 2026-03-10 | 4/200 |
|
|
[考研] 求調(diào)劑材料專碩293 +6 | 段_(:з」∠)_ 2026-03-10 | 6/300 |
|
|
[考研] 085602化工求調(diào)劑 +7 | 董boxing 2026-03-10 | 7/350 |
|
|
[考研] 0856材料與化工309分求調(diào)劑 +4 | ZyZy…… 2026-03-10 | 4/200 |
|
|
[考研] 085600材料與化工 326 求調(diào)劑 +4 | 熱愛生活ing 2026-03-09 | 4/200 |
|
|
[考研] 085600材料與化工,一志愿廣州985,求調(diào)劑 +15 | qqyyaill 2026-03-05 | 15/750 |
|
|
[考研] 320求調(diào)劑 +4 | 魏zy 2026-03-08 | 4/200 |
|
|
[考研] 337求調(diào)劑 +3 | 睡醒,。 2026-03-09 | 3/150 |
|
|
[考研] 297求調(diào)劑 +3 | 胡達(dá)靈 2026-03-05 | 5/250 |
|
|
[考研]
|
Sixuan wang 2026-03-06 | 7/350 |
|
|
[考研] 求調(diào)劑 +4 | 呼呼?~+123456 2026-03-06 | 4/200 |
|
|
[考研] 哈爾濱理工大學(xué)2026年研究生調(diào)劑,材料科學(xué)與化學(xué)工程學(xué)院研究生調(diào)劑 +3 | xinliu866 2026-03-06 | 3/150 |
|
|
[考研] 求調(diào)劑,學(xué)校研究所都可以,材料與化工267分 +6 | wmx1 2026-03-05 | 6/300 |
|
|
[考研] 材料調(diào)劑 +4 | L9370 2026-03-05 | 4/200 |
|