24小时热门版块排行榜    

CyRhmU.jpeg
南方科技大学公共卫生及应急管理学院2026级博士研究生招生报考通知(长期有效)
查看: 93  |  回复: 0
当前主题已经存档。
【有奖交流】积极回复本帖子,参与交流,就有机会分得作者 Shoney 的 30 个金币

Shoney

木虫 (著名写手)

[交流] 【求助】大家帮忙看看这个求最优参数的程序问题在哪儿?

是一个根据实验i数据求最优扩散系数的问题,即有一种物质扩散入一个圆柱体内,实验测得各个时间点的浓度,方程是一个PDE方程。运行结果没提示错误,但结果就是不对,最优的Deff值老是等于初值,不知道问题出在哪儿。(附件是源程序)程序如下:

function ACdisc_Diffusion
clear all;clc;
global Cm;
%动力学数据
tlist=[0 1 2 5 10 20 30 60];    %测量的时间点,单位min
Cm=[0 4.30 5.20 6.17 7.53 8.41 8.80 9.04];    %测得的AC中DMSO平均浓度,单位%w/w

%非线性拟合
C0=0;
Deff0=0.0116;    %单位mm2/min
[Deff,resnorm,resid,exitflag,output,lambda,jacobian]=lsqnonlin(@OptObjFunc,Deff0,[],[],[],tlist,C0,Cm)
ci=nlparci(Deff,resid,jacobian)

%--------------------------------------------------------------------------
%求解PDE
function C=OptObjFunc(Deff,tlist,C0,Cm)   
%PDE问题定义及参数初始化
g='ACdiscg1';    %定义求解域,几何尺寸单位mm
b='ACdiscb1';    %定义边界条件
d=['x/' num2str(Deff)];c='x';a=0;f=0;    %方程的系数

%网格化(三角形网格划分及网格细化)
[p,e,t]=initmesh(g);
[p,e,t]=refinemesh(g,p,e,t);
[p,e,t]=refinemesh(g,p,e,t);
p=jigglemesh(p,e,t);

%绘制PDE三角形网格
pdemesh(p,e,t);

u=parabolic(C0,tlist,b,p,e,t,c,a,f,d);    %计算得到的AC中任意时间任意位置DMSO浓度
%积分求解DMSO平均浓度的计算值
xi=linspace(0,3,100);yi=linspace(0,0.38,100);
Cc2(1,1)=0;    %赋值0,因进行网格转化时,值由线性插值得到;结果Cc2(1,1)≠0(边界值=9.04)
for i=2:length(tlist)
    Cc{1,i}=tri2grid(p,t,u(:,i),xi,yi);    %各个时间点分别进行网格转化
    Cc0=Cc{1,i};
    Cc1(1,i)=dblquad(@Func,0,3,0,0.38,[],@quadl,xi,yi,Cc0);    %积分求各个时间点对应的总浓度
    Cc2(1,i)=Cc1(1,i)/(3*0.38);    %平均浓度
end
C=Cc2-Cm;    %残差
   
%--------------------------------------------------------------------------
function f=Func(x,y,xi,yi,Cc0)
f=interp2(xi,yi,Cc0,x,y);

运行结果如下:
Optimization terminated: first-order optimality less than OPTIONS.TolFun,
and no negative/zero curvature detected in trust region model.
Deff =

    0.0116


resnorm =

    2.7230


resid =

         0   -1.0917   -0.7775    0.4317    0.6481    0.5194    0.2249    0.0000


exitflag =

     1


output =

    firstorderopt: 0
       iterations: 0
        funcCount: 2
     cgiterations: 0
        algorithm: 'large-scale: trust-region reflective Newton'
          message: [1x137 char]


lambda =

    lower: 0
    upper: 0


jacobian =

   All zero sparse: 8-by-1


ci =

  1.0e+007 *

   -9.8973    9.8973

[ Last edited by Shoney on 2009-9-15 at 10:42 ]
回复此楼
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 Shoney 的主题更新
普通表情 高级回复(可上传附件)
信息提示
请填处理意见