f*******i 发帖数: 8492 | 1 【 以下文字转载自 Physics 讨论区 】
发信人: fantastli (早晨从中午开始), 信区: Physics
标 题: 求一维Lennard-Jones分子模拟的思路,付代码,有包子
发信站: BBS 未名空间站 (Wed Dec 8 21:26:02 2010, 美东)
需要模拟一维空间下Lennard-Jones system的分子运动,现在程序写了一半,但是没有思路,求
达人指点一下,给个思路就可以。
我目前的思路是用离散的时间点,来描述分子的运动。取dt=0.001秒,然后每个阶段时间,来记录
一下分子间的距离(用于后期制作RDF用)。
我是用Lennard-Jones Force在每个时间点的变化来求加速度,然后再求出每个分子运动的轨
迹。
但是,现在的问题是,我不知道用什么条件来约束分子的运动。比如,当两个小球距离很近的时候,
其排斥力应该是很大的,此时应该怎么用程序来描述,之后小球又该如何改变运动,这一点实在是没
有头绪。
我贴出部分代码,matlab:
%Define the properties of the particles
x=zeros(50,1);
y=zeros(50,1);
v=zeros(50,1);
psize=13;
Diameter=1;
Mass=1;
%Define the particle number
N=input('Please input the particle number (N should be no larger than
50): N=');
%Define Sigma
sig=input('Please input the Sigma: sig=');
%Define Epsilon
eps=input('Please input the Epsilon: eps=');
%Define the intial collision time
CCount=0;
%Define the initial velocity and posistion for each particle
for i=1:N
x(i)=2*(i-1);
y(i)=1;
v(i)=2*randn(1)-1;
end
%Define the initial velocity of the first and last particle
v(1)=0;
v(N)=0;
%Define the initial time
t=0;
t_incr=1;
t_prev=0;
dt=0.001; %Use the scatter time line
while CCount<300 %Define the collision times
%Define the acceleration of each particle
for i=2:N-1; %The first and the last particle do not move
force(i)=0;
for j=1:N;
if(j~=i)
force(i)=force(i)+(2*(sig^12/(x(i)-x(j))^13)-
(sig^6/(x(i)-x(j))^7));
end
end
force(i)=force(i)*(-24)*eps; %Define the Lennard-Jones force for
each particle
a(i)=force(i)/Mass;
end
force(1)=0;
force(N)=0;
%For each dt, renew the velocity and position for each particle
for i=2:N-1
v(i)=v(i)+a(i)*dt;
x(i)=x(i)+v(i)*dt;
end
以下部分,我就觉得不对了。
%Detect whether the collision happens or not
for j=1:N-1
if (x(j+1)-x(j))<=Diameter
if j==1
v(2)=-v(2);
else
if j==N-1
v(N-1)=-v(N-1);
else
v_dum=v(j+1);
v(j+1)=v(j);
v(j)=v_dum;
end
end |
j*********m 发帖数: 33 | 2
有思路,求
时间,来记录
运动的轨
离很近的时候,
google
verlet algorithm
【在 f*******i 的大作中提到】 : 【 以下文字转载自 Physics 讨论区 】 : 发信人: fantastli (早晨从中午开始), 信区: Physics : 标 题: 求一维Lennard-Jones分子模拟的思路,付代码,有包子 : 发信站: BBS 未名空间站 (Wed Dec 8 21:26:02 2010, 美东) : 需要模拟一维空间下Lennard-Jones system的分子运动,现在程序写了一半,但是没有思路,求 : 达人指点一下,给个思路就可以。 : 我目前的思路是用离散的时间点,来描述分子的运动。取dt=0.001秒,然后每个阶段时间,来记录 : 一下分子间的距离(用于后期制作RDF用)。 : 我是用Lennard-Jones Force在每个时间点的变化来求加速度,然后再求出每个分子运动的轨 : 迹。
|
f*******i 发帖数: 8492 | 3 那个只能提供更好的速度模拟
我的问题是,用什么方法,或者物理性质来约束小球的运动范围?
假如我一排有5个小球,开始时时均匀分布,假定第一个和最后一个小球不动。中间三
个小球各有一个随
机的初始速度和运动方向。
那么中间三个小球中的任何一个,就会在其他四个小球的lennard-jones的叠加作用力
下进行运动。
那么该如何限定这三个小球的运动范围呢?比如用哪个物理性质来约束他们。
能给个思路吗?
【在 j*********m 的大作中提到】 : : 有思路,求 : 时间,来记录 : 运动的轨 : 离很近的时候, : google : verlet algorithm
|
j*********m 发帖数: 33 | 4
LJ 在近距离的排斥是无穷大的,中间三个球是无法越过最旁边的固定的两个的,只要
你的积分时间足够
小
【在 f*******i 的大作中提到】 : 那个只能提供更好的速度模拟 : 我的问题是,用什么方法,或者物理性质来约束小球的运动范围? : 假如我一排有5个小球,开始时时均匀分布,假定第一个和最后一个小球不动。中间三 : 个小球各有一个随 : 机的初始速度和运动方向。 : 那么中间三个小球中的任何一个,就会在其他四个小球的lennard-jones的叠加作用力 : 下进行运动。 : 那么该如何限定这三个小球的运动范围呢?比如用哪个物理性质来约束他们。 : 能给个思路吗?
|
S*M 发帖数: 10832 | 5 你看看lennard jones potential长啥样不就知道了
【在 f*******i 的大作中提到】 : 那个只能提供更好的速度模拟 : 我的问题是,用什么方法,或者物理性质来约束小球的运动范围? : 假如我一排有5个小球,开始时时均匀分布,假定第一个和最后一个小球不动。中间三 : 个小球各有一个随 : 机的初始速度和运动方向。 : 那么中间三个小球中的任何一个,就会在其他四个小球的lennard-jones的叠加作用力 : 下进行运动。 : 那么该如何限定这三个小球的运动范围呢?比如用哪个物理性质来约束他们。 : 能给个思路吗?
|
f*******i 发帖数: 8492 | 6 我有两点不明白,请指教一下:
1.因为LJ系统里的小球是soft sphere,那么如果在一维空间中,比如我按照顺序定义
为1,2,3,4,5号球,那么3号球,是否会发生在其他四个球的合力作用下,也不会越
过2号与4号球?
2. 根据LJ定义,sigma是两个球势能为0的时候,那么我该如何定义两个小球在距离很
近的时候,排斥力是无穷大? 比如如何定义那个距离? 还有如何定义排斥力?
因为排斥力和吸引力都是两个小球之间距离r的参变量,所以一开始都是连续变化的,
那么根据LJ的势能和LJ force的表达式,是否能够得到当两个小球的距离小于等于
sigma的时候,其LJ force也是逐渐快速变化到无穷大?
【在 j*********m 的大作中提到】 : : LJ 在近距离的排斥是无穷大的,中间三个球是无法越过最旁边的固定的两个的,只要 : 你的积分时间足够 : 小
|
j*********m 发帖数: 33 | 7
1.是
2.一切都在力的表达式里,什么都不用定义...
classical MD,只要学会了Verlet,其他的东西都是很自然的
【在 f*******i 的大作中提到】 : 我有两点不明白,请指教一下: : 1.因为LJ系统里的小球是soft sphere,那么如果在一维空间中,比如我按照顺序定义 : 为1,2,3,4,5号球,那么3号球,是否会发生在其他四个球的合力作用下,也不会越 : 过2号与4号球? : 2. 根据LJ定义,sigma是两个球势能为0的时候,那么我该如何定义两个小球在距离很 : 近的时候,排斥力是无穷大? 比如如何定义那个距离? 还有如何定义排斥力? : 因为排斥力和吸引力都是两个小球之间距离r的参变量,所以一开始都是连续变化的, : 那么根据LJ的势能和LJ force的表达式,是否能够得到当两个小球的距离小于等于 : sigma的时候,其LJ force也是逐渐快速变化到无穷大?
|
f*******i 发帖数: 8492 | 8 请问,我如果应用如下Verlet表达式,是否贴切?
位置更新:r(t+dt) = r(t) + v(t) * dt + 1/2 a(t) * dt^2
速度更新:v(t+dt) = v(t) + 1/2 * [a(t) + a(t+dt)] * dt
【在 j*********m 的大作中提到】 : : 1.是 : 2.一切都在力的表达式里,什么都不用定义... : classical MD,只要学会了Verlet,其他的东西都是很自然的
|
j*********m 发帖数: 33 | 9
你写的是velocity verlet, 一般的verlet只需要积一个方程就够了
【在 f*******i 的大作中提到】 : 请问,我如果应用如下Verlet表达式,是否贴切? : 位置更新:r(t+dt) = r(t) + v(t) * dt + 1/2 a(t) * dt^2 : 速度更新:v(t+dt) = v(t) + 1/2 * [a(t) + a(t+dt)] * dt
|
f*******i 发帖数: 8492 | 10 能否明示一下?
我不是做MD的,所以很多东西都是边查资料,边做的
【在 j*********m 的大作中提到】 : : 你写的是velocity verlet, 一般的verlet只需要积一个方程就够了
|
j*********m 发帖数: 33 | 11
x(t+dt)=2x(t)-x(t-dt)+a(t)dt^2
你可以去wiki看,verlet是最简单的守恒能量的解牛顿方程的方法,而且长时间误差也小
good night, thank you for baozi
【在 f*******i 的大作中提到】 : 能否明示一下? : 我不是做MD的,所以很多东西都是边查资料,边做的
|
f*******i 发帖数: 8492 | 12 多谢指教~
也小
【在 j*********m 的大作中提到】 : : x(t+dt)=2x(t)-x(t-dt)+a(t)dt^2 : 你可以去wiki看,verlet是最简单的守恒能量的解牛顿方程的方法,而且长时间误差也小 : good night, thank you for baozi
|