24小时热门版块排行榜    

查看: 3361  |  回复: 11

tandz

木虫 (正式写手)

[求助] 如何拟合morse potential的参数?

结合能E与距离r 已知. 请教各位,如何写个script拟合morse potential E(r) = D*(1-exp(a*(r-rm))^2的参数 D, a, rm呢? 本人几乎没有编程水平,请大家帮忙,非常感谢!
回复此楼
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
回帖支持 ( 显示支持度最高的前 50 名 )

月只蓝

主管区长 (职业作家)

【答案】应助回帖

★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ...
感谢参与,应助指数 +1
tandz: 金币+50, ★★★很有帮助, 太感谢了! 2013-04-04 17:03:15
tandz: 金币+20, ★★★★★最佳答案 2013-04-06 17:18:56
csgt0: 金币+2, 谢谢 2013-04-07 15:11:49
MATLAB代码可以吗?
假如已经获得了一组r,E的数据:
r:1 2 3 9
E:1 5 8 18
MATLAB程序如下,结果见附图1。
function feixianxingnihe
clear all;clc
format long

tspan=[1 2 3 9];    %r的数据,在此输入
xexp=[1 5 8 18];    %E的数据,在此输入

k0=[1 1 1];
lb=-[1 1 1]*1e5;
ub=[1 1 1]*1e5;


%-------------------------------------------------------------------------

% 使用函数lsqnonlin()进行参数估计

OPTIONS=optimset('MaxFunEvals',1000);
[k,resnorm,residual,exitflag,output,lambda,jacobian] = ...
    lsqnonlin(@ObjFunc,k0,lb,ub,OPTIONS,tspan,xexp);

ci = nlparci(k,residual,jacobian);
%residual;
fprintf('\n\n使用函数lsqnonlin()估计得到的参数值为:\n')
fprintf('\n\t参数 D = %.16f',k(1))
fprintf('\n\t参数 a = %.16f',k(2))
fprintf('\n\t参数rm = %.16f',k(3))
y=KineticsEqs(tspan,k);
R2=1-sum((xexp-y).^2)./sum((xexp-mean(y)).^2);
fprintf('\n\t相关系数之平方R^2 = %.16f',R2);
figure
plot(tspan,KineticsEqs(tspan,k),'b',tspan,xexp,'or'),legend('计算值','实验值','Location','Best')


%-------------------------------------------------------------------------

function f = ObjFunc(k,tspan,xexp)
f=KineticsEqs(tspan,k)-xexp;

%------------------------------------------------------------------------
function xt = KineticsEqs(t,k)
xt=k(1)*(1-exp(k(2).*(t-k(3)))).^2;

做方程拟合,1stopt代码简单易懂,1stopt代码如下:
Parameters D,a,rm;
Variable r,E;
Function E=D*(1-exp(a*(r-rm)))^2;
Data;
1 1
2 5
3 8
9 18
计算结果见附图2。   
可以看到两种方法计算结果相当接近。

附图1.jpg



附图2.jpg

MATLAB、MS小问题、普通问题请发帖求助!时间精力有限,恕不接受无偿私信求助。
5楼2013-04-04 16:27:04
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
普通回帖

dingd

铁杆木虫 (职业作家)

【答案】应助回帖

★ ★ ★ ★ ★ ★ ★ ★ ★ ★
感谢参与,应助指数 +1
tandz: 金币+10, ★★★很有帮助 2013-04-03 10:25:47
找个1stOpt,两分钟就能搞定。
2楼2013-04-03 10:15:16
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

tandz

木虫 (正式写手)

引用回帖:
2楼: Originally posted by dingd at 2013-04-03 10:15:16
找个1stOpt,两分钟就能搞定。

老板非要我写script, 没有办法啊
3楼2013-04-03 10:25:28
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

tandz

木虫 (正式写手)

自己顶,求助大家
4楼2013-04-03 19:21:49
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

tandz

木虫 (正式写手)

引用回帖:
5楼: Originally posted by 月只蓝 at 2013-04-04 16:27:04
MATLAB代码可以吗?
假如已经获得了一组r,E的数据:
r:1 2 3 9
E:1 5 8 18
MATLAB程序如下,结果见附图1。
function feixianxingnihe
clear all;clc
format long

tspan=;    %r的数据,在此输入
xexp ...

太感谢了!没学过MATLAB,我先试试。看起来效果非常好!
我学过一点点FORTRAN和C语言,不过现在都忘了,呵呵
6楼2013-04-04 17:06:25
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

tandz

木虫 (正式写手)

引用回帖:
5楼: Originally posted by 月只蓝 at 2013-04-04 16:27:04
MATLAB代码可以吗?
假如已经获得了一组r,E的数据:
r:1 2 3 9
E:1 5 8 18
MATLAB程序如下,结果见附图1。
function feixianxingnihe
clear all;clc
format long

tspan=;    %r的数据,在此输入
xexp ...

为什么我不能运行MATLAB程序呢,有几个错误
??? Error using ==> optim\private\lsqncommon
User supplied function failed with the following error:

Error using ==> feval
Undefined command/function 'ObjFunc'.

Error in ==> lsqnonlin at 147
[x,Resnorm,FVAL,EXITFLAG,OUTPUT,LAMBDA,JACOB] = ...
7楼2013-04-06 16:06:21
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

月只蓝

主管区长 (职业作家)

引用回帖:
7楼: Originally posted by tandz at 2013-04-06 16:06:21
为什么我不能运行MATLAB程序呢,有几个错误
??? Error using ==> optim\private\lsqncommon
User supplied function failed with the following error:

Error using ==> feval
Undefined command/func ...

file-new-m文件或者function新建一个m文件,上述代码复制进m文件中,F5运行。我刚又测试过,可以运行。

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

MATLAB、MS小问题、普通问题请发帖求助!时间精力有限,恕不接受无偿私信求助。
8楼2013-04-06 18:08:39
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

tandz

木虫 (正式写手)

送红花一朵
引用回帖:
8楼: Originally posted by 月只蓝 at 2013-04-06 18:08:39
file-new-m文件或者function新建一个m文件,上述代码复制进m文件中,F5运行。我刚又测试过,可以运行。...

运行了!谢谢帮助!
9楼2013-04-07 09:58:22
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

ytudou

银虫 (小有名气)

引用回帖:
5楼: Originally posted by 月只蓝 at 2013-04-04 16:27:04
MATLAB代码可以吗?
假如已经获得了一组r,E的数据:
r:1 2 3 9
E:1 5 8 18
MATLAB程序如下,结果见附图1。
function feixianxingnihe
clear all;clc
format long

tspan=;    %r的数据,在此输入
xexp ...

大神,你看他这个势函数形式没问题吗?morse势不是形如u=D{exp(-2α(r-ro))-2exp(-α(r-ro))}的样子吗?还是说我写出来的这个可以通过变形变成那个?请赐教
要给力啊,MS~LAMMPS~
10楼2013-10-11 14:37:54
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 tandz 的主题更新
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[考研] 085600材料与化工 求调剂 +11 enenenhui 2026-03-13 12/600 2026-03-14 21:48 by laoshidan
[考研] 材料与化工 323 英一+数二+物化,一志愿:哈工大 本人本科双一流 +4 自由的_飞翔 2026-03-13 5/250 2026-03-14 19:39 by hmn_wj
[考研] 290求调剂 +4 @将就将就看 2026-03-10 8/400 2026-03-14 14:23 by 千千运气
[考研] 085600材料与化工 326 求调剂 +5 热爱生活ing 2026-03-09 5/250 2026-03-14 02:39 by JourneyLucky
[考研] 一志愿浙江大学0856材料与化工求调剂 +4 yansheng@211 2026-03-09 5/250 2026-03-14 02:10 by JourneyLucky
[考研] 308求调剂 +4 是Lupa啊 2026-03-09 4/200 2026-03-14 02:06 by tranquil_ya
[考研] 求调剂! +4 朔朔话 2026-03-09 4/200 2026-03-14 01:38 by JourneyLucky
[考研] 306求调剂 +4 唐薏薏 2026-03-09 4/200 2026-03-14 01:19 by JourneyLucky
[考研] 一志愿安徽大学材料工程专硕313分,求调剂的学校 +8 Yu先生 2026-03-10 10/500 2026-03-14 01:04 by JourneyLucky
[考研] 一志愿西南交大,材料专硕317求调剂 +5 lx8568 2026-03-11 5/250 2026-03-13 21:43 by peike
[考研] 301求调剂 +6 Liyouyumairs 2026-03-11 6/300 2026-03-13 20:11 by JourneyLucky
[考研] 【0856】化学工程(085602)313 分,本科学科评估A类院校化学工程与工艺,诚求调剂 +7 小刘快快上岸 2026-03-11 7/350 2026-03-13 16:06 by ruiyingmiao
[考研] 一志愿211化学学硕310分求调剂 +8 努力奋斗112 2026-03-12 9/450 2026-03-13 15:41 by JourneyLucky
[考研] 285求调剂 +4 ytter 2026-03-12 4/200 2026-03-13 14:48 by jxchenghu
[考研] 085600材料与化工 309分请求调剂 +7 dtdxzxx 2026-03-12 8/400 2026-03-13 14:43 by jxchenghu
[考研] 304求调剂(085602一志愿985) +12 化工人999 2026-03-09 12/600 2026-03-13 12:02 by JourneyLucky
[考研] 289求调剂 +3 李政莹 2026-03-12 3/150 2026-03-13 11:02 by 求调剂zz
[考研] 化工学硕306求调剂 +9 42838695 2026-03-12 9/450 2026-03-13 10:16 by houyaoxu
[考研] 321求调剂(食品/专硕) +3 xc321 2026-03-12 6/300 2026-03-13 08:45 by xc321
[考研] 290求调剂 +3 柯淮然 2026-03-10 8/400 2026-03-11 13:48 by 柯淮然
信息提示
请填处理意见