24小时热门版块排行榜    

Znn3bq.jpeg
汕头大学海洋科学接受调剂
查看: 4433  |  回复: 12
本帖产生 1 个 计算强帖 ,点击这里进行查看
当前只显示满足指定条件的回帖,点击这里查看本话题的所有回帖

258190169

铜虫 (小有名气)

[求助] matlab-常微分方程参数估计

初始数据浓度和时间
t=[0,10,30,50,70,90,110,130,150,160];
c=[0,0.23211,0.45906,0.68601,0.92328,1.21213,1.32561,1.34624,1.39782,1.398];
微分方程,dc/dt=[4.41/96485-(4.41*k+L*4.41/96485)c]/(1+4.41*k*t)

要求: 1.得到拟合参数:k 和L 以及相对偏差
       2.得到拟合曲线和数据点的图
       3.最好附上院程序
回复此楼
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

feng198800

金虫 (正式写手)

引用回帖:
3楼: Originally posted by dbb627 at 2011-12-16 17:31:51
function PenicilliumEst
clear all;
t=;
y=;
y0=0;

% Nonlinear least square estimate using lsqnonlin()
beta0=;
lb=;ub=;
= ...
    lsqnonlin(@Func,beta0,lb,ub,[],t,y,y0);         
ci = n ...

请问大侠,初值能不能结合最优化(遗传算法或者PSO算法)寻求?如何把这和上面的程序结合起来?
13楼2012-05-25 18:38:20
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
查看全部 13 个回答

258190169

铜虫 (小有名气)

本人编写的程序如下但是无法运行,由于是新手,希望高手帮忙调试一下:
function PenicilliumEst
clear all;
t=[0,10,30,50,70,90,110,130,150,160];
y=[0,0.23211,0.45906,0.68601,0.92328,1.21213,1.32561,1.34624,1.39782,1.398];
y0=0;

% Nonlinear least square estimate using lsqnonlin()
beta0=[0.005 0.001];
lb=[0 0];ub=[inf inf];
[beta,resnorm,residual,exitflag,output,lambda,jacobian] = ...
    lsqnonlin(@Func,beta0,lb,ub,[],t,y);         
ci = nlparci(beta,residual,jacobian);

% =======================================
function f = Func(beta,t,y,y0)        % Define objective function
tspan = [0  max(x)];
[tt yy] = ode45(@ModelEqs,tspan,y0,[],beta);
yc= spline(tt,yy,x);
f1=y-yc
% ==================================
function dydt = ModelEqs(t,y,beta)          % Model equations
dydt = [4.41/96485-(4.41*beta(1)+beta(2)*4.41/96485)*y]/(1+4.41*beta(1)*t)

% result
fprintf('\n Estimated Parameters by Lsqnonlin():\n')
fprintf('\t k1 = %.4f ± %.4f\n',beta(1),ci(1,2)-beta(1))
fprintf('\t k2 = %.4f ± %.4f\n',beta(2),ci(2,2)-beta(2))
fprintf('  The sum of the residual squares is: %.1e\n\n',sum(residual.^2))

% plot of fit results
tspan = [0  max(t)];
[tt yc] = ode45(@modeleqs,tspan,c0,[],beta);
tc=linspace(0,max(t),200);
yc = spline(tt,yc,tc);
plot(t,c,'ro',tc,yca,'r-');
hold on
xlabel('Time');
ylabel('Concentration');
hold off
2楼2011-12-16 16:49:37
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

dbb627

荣誉版主 (著名写手)

【答案】应助回帖

★ ★ ★ ★ ★
感谢参与,应助指数 +1
258190169(金币+20): 多谢大侠 你的QQ是多少可以和你联系一下嘛 2011-12-16 22:12:03
cenwanglai(金币+5, 计算强帖+1): 谢谢给予帮助~ 2011-12-20 09:07:46
引用回帖:
2楼: Originally posted by 258190169 at 2011-12-16 16:49:37:
本人编写的程序如下但是无法运行,由于是新手,希望高手帮忙调试一下:
function PenicilliumEst
clear all;
t=[0,10,30,50,70,90,110,130,150,160];
y=[0,0.23211,0.45906,0.68601,0.92328,1.21213,1.32561, ...

CODE:
function PenicilliumEst
clear all;
t=[0,10,30,50,70,90,110,130,150,160];
y=[0,0.23211,0.45906,0.68601,0.92328,1.21213,1.32561,1.34624,1.39782,1.398];
y0=0;

% Nonlinear least square estimate using lsqnonlin()
beta0=[0.001 0.001];
lb=[0 0];ub=[inf inf];
[beta,resnorm,residual,exitflag,output,lambda,jacobian] = ...
    lsqnonlin(@Func,beta0,lb,ub,[],t,y,y0);         
ci = nlparci(beta,residual,jacobian);
beta
% result
fprintf('\n Estimated Parameters by Lsqnonlin():\n')
fprintf('\t k1 = %.4f ± %.4f\n',beta(1),ci(1,2)-beta(1))
fprintf('\t k2 = %.4f ± %.4f\n',beta(2),ci(2,2)-beta(2))
fprintf('  The sum of the residual squares is: %.1e\n\n',sum(residual.^2))

% plot of fit results
tspan = [0  max(t)];
[tt yc] = ode45(@ModelEqs,tspan,y0,[],beta);
tc=linspace(0,max(t),200);
yca = spline(tt,yc,tc);
plot(t,y,'ro',tc,yca,'r-');
hold on
xlabel('Time');
ylabel('Concentration');
hold off
% =======================================
function f1 = Func(beta,t,y,y0)        % Define objective function
tspan =t;
[tt yy] = ode45(@ModelEqs,tspan,y0,[],beta);
yc= spline(tt,yy,t);
f1=y-yc;
% ==================================
function dydt = ModelEqs(t,y,beta)          % Model equations
dydt = (4.41/96485-(4.41*beta(1)+beta(2)*4.41/96485)*y)/(1+4.41*beta(1)*t);

改的可以运行了,但是初值不合适
The more you learn, the more you know, the more you know, and the more you forget. The more you forget, the less you know. So why bother to learn.
3楼2011-12-16 17:31:51
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

dingd

铁杆木虫 (职业作家)


dbb627(金币+1): 感谢参与 2011-12-17 11:40:34
用1stOpt试试:
CODE:
Variable t,c;
ODEFunction c'=(4.41/96485-(4.41*k+L*4.41/96485)*c)/(1+4.41*k*t);
Data;
0        0
10        0.23211
30        0.45906
50        0.68601
70        0.92328
90        1.21213
110        1.32561
130        1.34624
150        1.39782
160        1.398

均方差(RMSE): 0.272257391318038
残差平方和(SSE): 0.66711678414573
相关系数(R): 0.932733204157737
相关系数之平方(R^2): 0.869991230138359
决定系数(DC): 0.576237758605965

参数                  最佳估算
--------------------        -------------
k        -2.32798002359073
l        514537.340812352


» 本帖已获得的红花(最新10朵)

4楼2011-12-17 09:52:26
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[考研] 一志愿华中农业071010,320求调剂 +5 困困困困坤坤 2026-04-14 6/300 2026-04-14 23:54 by Xurambo2014
[考研] 211本科材料化工求调剂 +19 YHLAH 2026-04-11 23/1150 2026-04-14 22:25 by fenglj492
[考研] 26药学专硕105500求调剂 +6 喽哈加油 2026-04-13 6/300 2026-04-14 16:40 by zhouxiaoyu
[考研] 366求调剂 +11 不知名的小卅 2026-04-11 11/550 2026-04-14 15:50 by zs92450
[考研] 求调剂 +20 MAX怅惘 2026-04-09 22/1100 2026-04-14 14:57 by 独醉梦孤城
[教师之家] 转长聘了 +7 简单化xn 2026-04-13 7/350 2026-04-14 14:50 by xindong
[基金申请] RY:中国产出的科学垃圾论文,绝对数量和比例都世界第一 +6 zju2000 2026-04-14 17/850 2026-04-14 14:34 by jurkat.1640
[考研] 297求调剂 +23 ORCHID1 2026-04-10 26/1300 2026-04-14 13:52 by 陈皮皮
[考研] 085600材料与化工329分求调剂 +24 叶zilin 2026-04-13 25/1250 2026-04-14 09:20 by 试管破裂
[教师之家] 山东双非院校考核超级无底线,领导幸灾乐祸,教师遭殃恐 +3 qut2026 2026-04-11 7/350 2026-04-12 20:24 by qut2026
[考研] 求调剂 +16 张番茄不炒蛋 2026-04-10 17/850 2026-04-12 13:58 by 熬夜成!
[考研] 调剂求助 +6 果然有我 2026-04-11 7/350 2026-04-11 16:22 by 明月此时有
[考研] 283求调剂,工科! +12 苏打水7777 2026-04-08 12/600 2026-04-11 10:28 by 逆水乘风
[考研] 297求调剂 +9 Kwgyz 2026-04-09 9/450 2026-04-11 10:09 by zhq0425
[考研] 284求调剂 +12 archer.. 2026-04-10 13/650 2026-04-11 08:44 by zhq0425
[考研] 考研调剂 +26 硕星赴 2026-04-09 27/1350 2026-04-10 22:24 by 猪会飞
[考研] 070300化学 求调剂 +13 73372112 2026-04-08 13/650 2026-04-09 20:22 by maddjdld
[考研] 283电子信息求调剂 +4 三石WL 2026-04-08 4/200 2026-04-09 10:21 by wp06
[考研] 334求调剂 +16 Riot2025 2026-04-08 17/850 2026-04-09 09:28 by wdyheheeh
[考研] 286求调剂 +19 Faune 2026-04-08 20/1000 2026-04-09 08:36 by 哦哦123
信息提示
请填处理意见