【求助】用matlab最优化方法进行参数拟合 - 数学 - 计算数学

版块导航: 正在加载中...

客户端APP下载

论文辅导

申博辅导

登录注册

应《网络安全法》要求，自2017年10月1日起，未进行实名认证将不得使用互联网跟帖服务。为保障您的帐号能够正常使用，请尽快对帐号进行手机号验证，感谢您的理解与支持！

24小时热门版块排行榜

返回列表

本帖产生 1 个数学EPI ，点击这里进行查看

当前只显示满足指定条件的回帖，点击这里查看本话题的所有回帖

leedobb

金虫 (正式写手)

数学EPI: 1
应助: 77 (初中生)
金币: 152.4
散金: 646
红花: 15
帖子: 856
在线: 270.4小时
虫号: 1199535
注册: 2011-02-06
性别: GG
专业: 高分子科学

★ ★ ★ ★ ★ ★
小木虫(金币+0.5):给个红包，谢谢回帖交流
nono2009(金币+5, 数学EPI+1): 鼓励详细应助 2011-03-22 07:48:24

引用回帖:

Originally posted by wodaifei at 2011-03-21 08:48:28:
各位学哥、学姐，大侠们！！
小弟最近遇到一个难题：具体为 dyA/dt=-(k1+k2+k3+k4)yA
dyB/dt=k1yA-(k5+k6+k7)yB
...

% 哦，你的题目是这个啊，
% 这只能依靠数值解，matlab肯定没有现成的程序吧。
%
% 你这样做吧。构造一个能量函数（其实就是最小二乘法啦）
%
% F =  (yA(t2)- YYA(t2,k1,k2,k3...,k9) ) ^ 2
% ..嗯，算我直接编给你吧。

%
% dyA/dt=-(k1+k2+k3+k4)yA
%                                                       dyB/dt=k1yA-(k5+k6+k7)yB
%                                                       dyC/dt=k2yA+k5yB-(k8+k9)yC
%                                                       dyD/dt=k3yA+k6yB+k8yC
%                                                       dyE/dt=k4yA+k7yB+k9yC
% A       B       C       D       E
%                   0.833    0.0632    0．699    0.187 0.0256    0.0220
%                            1.25       0.055       0.722       0.199 0.0305    0.0281
%                            2.5       0.0502       0.769       0.201 0.0328       0.031
clear all

kk0(1:9)=[1.2,0.9,0.8, 1.3,0.6,1.6,  1.7,1.1,1.0];
t1= 1.25  - 0.833;
t2=2.5 -0.833;
y0 = [0.0632    0.699    0.187 0.0256    0.0220]';
y1 = [ 0.055       0.722       0.199 0.0305    0.0281]';
y2 = [ 0.0502       0.769       0.201 0.0328       0.031]';
dk = 0.0001;
dt = 0.1
kk=kk0
for i=1:100000
  %%% 计算 dF/dk
for j=1:9
%% 计算k+dk
k=kk;k(j)=k(j)+dk;
A = [ -(k(1)+k(2)+k(3)+k(4)) 0 0
k(1)  -(k(5)+k(6)+k(7))  0
   k(2)  k(5)  -(k(8)+k(9)) ];
  [V,D]=eig(A);
  DD(1,1)=D(1,1); DD(2,1) = D(2,2); DD(3,1) = D(3,3);
  DD =DD;
  VD=V^-1;
  Z0 = VD*y0(1:3);
  %%% z = z0*exp( DD(i)*t)
  Z1 = Z0 .*exp(DD*t1);
  Y1 = V *Z1;
  Int_Z1 = Z1./DD-Z0./DD; %int (Z1) = Z1/DD+C
  Y1(4) = k(3)*sum(V(1,1:3)*Int_Z1)+k(6)*sum(V(2,1:3)*Int_Z1) + k(8)*sum(V(3,1:3)*Int_Z1) + y0(4);
  Y1(5) = k(4)*sum(V(1,1:3)*Int_Z1)+k(7)*sum(V(2,1:3)*Int_Z1) + k(9)*sum(V(3,1:3)*Int_Z1) + y0(4);
  Z2 = Z0 .*exp(DD*t2);
  Y2 = V *Z2;
  Int_Z2 = Z2./DD-Z0./DD; %int (Z1) = Z1/DD+C
  Y2(4) = k(3)*sum(V(1,1:3)*Int_Z2)+k(6)*sum(V(2,1:3)*Int_Z2) + k(8)*sum(V(3,1:3)*Int_Z2) + y0(4);
  Y2(5) = k(4)*sum(V(1,1:3)*Int_Z2)+k(7)*sum(V(2,1:3)*Int_Z2) + k(9)*sum(V(3,1:3)*Int_Z2) + y0(4);
  F_up = ((Y2-y2)'*(Y2-y2) + (Y1-y1)'*(Y1-y1)) /2;

%% 计算k=k-dk
  k=kk;k(j)=k(j)-dk;
A = [ -(k(1)+k(2)+k(3)+k(4)) 0 0
k(1)  -(k(5)+k(6)+k(7))  0
   k(2)  k(5)  -(k(8)+k(9)) ];
  [V,D]=eig(A);
  DD(1,1)=D(1,1); DD(2,1) = D(2,2); DD(3,1) = D(3,3);
  DD =DD;
  VD=V^-1;
  Z0 = VD*y0(1:3);
  %%% z = z0*exp( DD(i)*t)
  Z1 = Z0 .*exp(DD*t1);
  Y1 = V *Z1;
  Int_Z1 = Z1./DD-Z0./DD; %int (Z1) = Z1/DD+C
  Y1(4) = k(3)*sum(V(1,1:3)*Int_Z1)+k(6)*sum(V(2,1:3)*Int_Z1) + k(8)*sum(V(3,1:3)*Int_Z1) + y0(4);
  Y1(5) = k(4)*sum(V(1,1:3)*Int_Z1)+k(7)*sum(V(2,1:3)*Int_Z1) + k(9)*sum(V(3,1:3)*Int_Z1) + y0(4);
  Z2 = Z0 .*exp(DD*t2);
  Y2 = V *Z2;
  Int_Z2 = Z2./DD-Z0./DD; %int (Z1) = Z1/DD+C
  Y2(4) = k(3)*sum(V(1,1:3)*Int_Z2)+k(6)*sum(V(2,1:3)*Int_Z2) + k(8)*sum(V(3,1:3)*Int_Z2) + y0(4);
  Y2(5) = k(4)*sum(V(1,1:3)*Int_Z2)+k(7)*sum(V(2,1:3)*Int_Z2) + k(9)*sum(V(3,1:3)*Int_Z2) + y0(4);
  F_down = ((Y2-y2)'*(Y2-y2) + (Y1-y1)'*(Y1-y1)) /2;
  FF= (F_down+F_up)/2;
  DF(j)= (F_up-F_down)/dk/2;
end

for j=1:9
kk(j)= kk(j) - DF(j) *dt;
end

if(mod(i,1000)==1)
'deviation'
FF/10
'k1->k9', kk
end







end

[ Last edited by leedobb on 2011-3-21 at 15:46 ]