 |
 |
 |
 |
 |
 |
||
 |
 |
 |
 |
 |
| 5 | 1/1 | 返回列表 |
| 查看: 1666 | 回復(fù): 9 | ||
| 【懸賞金幣】回答本帖問(wèn)題，作者h(yuǎn)zd250將贈(zèng)送您 10 個(gè)金幣 | ||
| 當(dāng)前只顯示滿足指定條件的回帖，點(diǎn)擊這里查看本話題的所有回帖 | ||
[求助]
matlab擬合反應(yīng)動(dòng)力學(xué) 已有2人參與
|
||
|
剛?cè)腴Tmatlab的小白，最近想做反應(yīng)的動(dòng)力學(xué)，按照B站up主的視頻自己寫了一段代碼，但是運(yùn)行總是出問(wèn)題： 錯(cuò)誤使用 odearguments (第 93 行)；FUNC 必須返回列向量，出錯(cuò) ode45 (第 115 行)，odearguments(FcnHandlesUsed, solver_name, ode, tspan, y0, options, varargin); 出錯(cuò) Kinetics>fun (第 67 行)，[t,x]=ode45(@func,tspan,x0,[],k)； 出錯(cuò) lsqnonlin (第 218 行)，initVals.F = feval(funfcn{3},xCurrent,varargin{:});出錯(cuò) Kinetics (第 18 行)，lsqnonlin(@fun,k0,lb,ub,[],yexp);%非線性最小二乘法。原因:Failure in initial objective function evaluation. LSQNONLIN cannot continue. 下面是我寫的代碼，讀取的Excel表格里有5列*7行的實(shí)驗(yàn)數(shù)據(jù)，勞煩大佬幫我瞅瞅哪里需要改動(dòng)，萬(wàn)分感謝。 function Kinetics %反應(yīng)一：A+B=C+M %r=k*XA*XB-K*XC*XM %反應(yīng)二：A+C=D+M %r=K*XA*XC-K*XD*XM %反應(yīng)三：B+C=E+M %r=K*XB*XC-K*XE*XM %XM=0.175 clc clear all; global a b tspan=[0.5 1 4 6 8 12 16]; yexp=xlsread('reaction.xls'); k0=[0.1 0.01 0.01 0.001 0.001 0.001];%參數(shù)初值 lb=[0 0 0 0 0 0];%下邊界 ub=[+inf +inf +inf +inf +inf +inf];%上邊界 [k,resnorm,residual,exitflag,output,lambda,jacobian]=... lsqnonlin(@fun,k0,lb,ub,[],yexp);%非線性最小二乘法 tspan=[0.5 1 4 6 8 12 16]; a=1; b=a+6; x0=yexp(a,  ;%積分初值[t,x]=ode45(@func,tspan,x0,[],k); t1=linspace(0.5,16,200); ya1=spline(t,x(:,1),t1);%動(dòng)力學(xué)計(jì)算得到的點(diǎn)進(jìn)行樣條插值 ya2=spline(t,x(:,2),t1); ya3=spline(t,x(:,3),t1); ya4=spline(t,x(:,4),t1); ya5=spline(t,x(:,5),t1); for m=1:7 for n=1:5 yy(a+m-1,n)=x(m,n);%每一次的值存入yy矩陣 end end figure(1) plot(tspan,yexp(a:b,1),'k^',t1,ya1,'k-',tspan,yexp(a:b,2),'ro',t1,ya2,'r-',tspan,yexp(a:b,3),'bd',t1,ya3,'b-',... tspan,yexp(a:b,4),'g*',t1,ya4,'g-',tspan,yexp(a:b,5),'yp',t1,ya5,'y-'); legend('','A濃度','','B濃度','','C濃度','','D濃度','','E濃度'); xlabel('t(h)');ylabel('濃度(mol/L)');title('170℃ 0.1wt%催化劑'); t1=linspace(0.5,16,200); z1=spline(t,yy(1:7,1),t1); h1=spline(t,yy(1:7,2),t1); s1=spline(t,yy(1:7,3),t1); b1=spline(t,yy(1:7,4),t1); u1=spline(t,yy(1:7,5),t1); xlswrite('result.xls',[t1' z1' h1' s1' b1' u1'],'sheet1'); xlswrite('result.xls',residual,'sheet2'); Ne = length(yexp(:,2)); %模型適定性判別 Np = length(k); [rho2,F] = rho2_F(k,yexp,resnorm,Ne,Np); ci=nlparci(k,residual,jacobian) fprintf('\t k1,0=%.1f ± %.4f\n',k(1),ci(1,2)-k(1)); fprintf('\t k2,0=%.1f ± %.4f\n',k(2),ci(2,2)-k(2)); fprintf('\t k3,0=%.1f ± %.4f\n',k(3),ci(3,2)-k(3)); fprintf('\t k4,0=%.1f ± %.4f\n',k(4),ci(4,2)-k(4)); fprintf('\t k5,0=%.1f ± %.4f\n',k(5),ci(5,2)-k(5)); fprintf('\t 殘差平方和：%.3f\n',resnorm) fprintf('\t 實(shí)驗(yàn)點(diǎn)數(shù)和自由度分別為 Ne = %d和 Np = %d\n',Ne,Np) fprintf('\t 決定性指標(biāo)ρ^2: %.4f\n',rho2) fprintf('\t F比: %.3f\n\n',F) %================================================================================= function f=fun(k,yexp) f=[]; tspan=[0.5 1 4 6 8 12 16]; a=1; x0=yexp(a, ;[t,x]=ode45(@func,tspan,x0,[],k) d=a+6; yc1=x(:,1); yc2=x(:,2); yc3=x(:,3); yc4=x(:,4); yc5=x(:,5); f11=yexp(a:d,1)-yc1; f12=yexp(a:d,2)-yc2; f13=yexp(a:d,3)-yc3; f14=yexp(a:d,4)-yc4; f15=yexp(a:d,5)-yc5; ff=[f11 f12 f13 f14 f15]; f=[f;ff]; %================================================================================= function dxdt=func(t,x,k) r1=-k(1)*x(1)*x(2)-k(2)*x(1)*x(3)+k(4)*x(3)*0.175+k(5)*x(4)*0.175; r2=-k(1)*x(1)*x(2)-k(3)*x(2)*x(3)+k(4)*x(3)*0.175+k(6)*x(5)*0.175; r3=k(1)*x(1)*x(2)+k(5)*x(4)*0.175+k(6)*x(5)*0.175-k(2)*x(1)*x(3)-k(3)*x(2)*x(3)-k(4)*x(3)*0.175; r4=k(2)*x(1)*x(3)-k(5)*x(4)*0.175; r5=k(3)*x(2)*x(3)-k(6)*x(5)*0.175; dxdt=[r1 r2 r3 r4 r5] %================================================================================= function [rho2,F] = rho2_F(k,yexp,s,Ne,Np) y=yexp.^2; sy = sum(y( );rho2 = 1 - s/sy; %rho2: 決定性指標(biāo) F = (sy - s)*(Ne-Np)/(Np*s); %F：F比 |
回復(fù)此樓» 猜你喜歡
幫我的英語(yǔ)口語(yǔ)老師找學(xué)生
已經(jīng)有4人回復(fù)
-
計(jì)算三階力常數(shù)用thirdorder擴(kuò)胞時(shí)參數(shù)如何選擇？臨近原子數(shù)怎么取值
已經(jīng)有7人回復(fù)
-
物理學(xué)I論文潤(rùn)色/翻譯怎么收費(fèi)?
已經(jīng)有144人回復(fù)
-
2026年循環(huán)經(jīng)濟(jì)功能材料國(guó)際會(huì)議（ICFMCE 2026）
已經(jīng)有0人回復(fù)
-
2026年第五屆電氣、電子與信息工程國(guó)際會(huì)議（ISEEIE 2026）
已經(jīng)有0人回復(fù)
-
哈爾濱理工大學(xué)物理系招收物理學(xué)考研調(diào)劑
已經(jīng)有1人回復(fù)
-
0702一志愿吉大B區(qū)求調(diào)劑
已經(jīng)有2人回復(fù)
-
求調(diào)劑
已經(jīng)有0人回復(fù)
6樓2021-05-18 21:31:33
hzlhm
至尊木蟲 (著名寫手)
- 應(yīng)助: 387 (碩士)
- 金幣: 17808.5
- 紅花: 53
- 帖子: 2879
- 在線: 606.3小時(shí)
- 蟲號(hào): 1517335
- 注冊(cè): 2011-11-30
- 性別: GG
- 專業(yè): 常微分方程與動(dòng)力系統(tǒng)
hzd250
2樓2021-05-17 22:00:00
3樓2021-05-18 10:22:08
4樓2021-05-18 10:52:35
| 最具人氣熱帖推薦 [查看全部] | 作者 | 回/看 | 最后發(fā)表 | |
|---|---|---|---|---|
|
[考研] 考研材料與化工，求調(diào)劑 +4 | 戲精丹丹丹 2026-03-09 | 4/200 |
|
|---|---|---|---|---|
|
|
[考研] 復(fù)試調(diào)劑 +6 | 呼呼？~+123456 2026-03-08 | 8/400 |
|
|
|
[考博] 2026博士申請(qǐng) +6 | 起泡酒 2026-03-08 | 6/300 |
|
|
[考研] 294 英二數(shù)二物化 求調(diào)劑 +6 | 米飯團(tuán)不好吃 2026-03-09 | 6/300 |
|
|
|
[考研] 2026考研求調(diào)劑-材料類-本科211一志愿985-初試301分 +7 | 蟲友233 2026-03-07 | 7/350 |
|
|
|
[考研] 求調(diào)劑，一志愿華中科大0702，數(shù)一英一，293 +4 | 小羅露一二 2026-03-07 | 4/200 |
|
|
|
[考研] 0817化學(xué)工程與技術(shù)312分求調(diào)劑 +7 | T123 tt 2026-03-04 | 7/350 |
|
|
|
[考研] 化學(xué)工程求調(diào)劑 +12 | 化工人999 2026-03-04 | 12/600 |
|
|
|
[考研] 337一志愿華南理工材料求調(diào)劑 +4 | mysdl 2026-03-07 | 4/200 |
|
|
|
[考研] 一志愿鄭大071000分?jǐn)?shù)282求調(diào)劑 +3 | 研研顏 2026-03-05 | 7/350 |
|
|
|
[考研] 一志愿211 085600 280數(shù)二英二求調(diào)劑 +3 | 月山斜 2026-03-06 | 3/150 |
|
|
|
[考研] 一志愿哈爾濱工業(yè)大學(xué)0856材料與化工，前三科206，總分283，求調(diào)劑 +7 | 26考研求調(diào)劑 2026-03-06 | 7/350 |
|
|
|
[考研] 材料專碩290求調(diào)劑 +8 | 杰尼龜aaa 2026-03-04 | 8/400 |
|
|
|
[考博] 2026年博士名額撿漏 +4 | 科研ya 2026-03-04 | 7/350 |
|
|
|
[考研] 332材料求調(diào)劑 +6 | zjy101327 2026-03-05 | 7/350 |
|
|
|
[考研] 334求調(diào)劑 +6 | Trying] 2026-03-05 | 8/400 |
|
|
|
[考研] 材料085600 303求調(diào)劑 +7 | 1bygone 2026-03-04 | 7/350 |
|
|
|
[論文投稿]
 主編終審卡者不動(dòng)
 100+4
|
Stray2021 2026-03-03 | 4/200 |
|
|
|
[考研] 325求調(diào)劑 +5 | 學(xué)家科 2026-03-04 | 5/250 |
|
|
|
[考研] 293求調(diào)劑 +4 | 是樂(lè)渝哇 2026-03-03 | 4/200 |
|