金属铝分子动力学模拟文档格式.docx
- 文档编号:18599801
- 上传时间:2022-12-29
- 格式:DOCX
- 页数:16
- 大小:290.80KB
金属铝分子动力学模拟文档格式.docx
《金属铝分子动力学模拟文档格式.docx》由会员分享,可在线阅读,更多相关《金属铝分子动力学模拟文档格式.docx(16页珍藏版)》请在冰豆网上搜索。
Brenner发现金属单晶晶须拉伸强度与晶须直径呈反比,Fleck在微米级细铜丝的扭转试验中观察到尺寸效应纳米电机系统(NEMS)的出现同迫切要求了解纳米尺度下材料的力学行为,当前从实验上较难获得详细的信息,而分子动力学模拟可以提供相关细节.
分子动力学通过直接模拟原子的运动过程,使我们能够详细了解模拟对象的演化发展历史分子动力学模拟的一个关键在于原子势函数的选取原子势早期一般采用简单的对势,但对势无法正确描述弹性常数,其结果不理想世纪年代提出的镶嵌原子法、有效介质理论更客观地反映了原子间多体作用的本质,可得到较合理的结果.认为体系总能量为
有效地对分子动力学的数据进行后处理也是一个重要的研究方向,因为对数据进行后处理,获取有用的信息也是一个很烦琐的工作。
现有的一些分子动力学软件功能不够强大,只是偏重于某一行业,通用性不高,还有一些软件为自由软件,可维护性不强。
还没有出现集建模、求解、后处理于一体的分子动力学软件。
因此,针对特定的问题进行自行编程显得尤为重要。
通过分子动力学方法,模拟了铝分子拉伸实验中的形变过程.研究了晶体取向裂纹的形变特点和断裂机理,观察到各种形变现象,如位错形核和发射,位错运动,堆垛层错或孪晶的形成,纳米空洞的形成与连接等.计算结果表明,裂纹扩展是塑性过程和弹性过程相结合的过程,其中塑性过程表现为由裂尖发射的位错导致的原子切变行为,而弹性过程的发生则是由无位错区中的原子断键所导致.本研究采用VC++自行编程,VMD动画演示得到原子在拉伸过程中的变化,对所研究的问题进行求解。
2原理
I)其基本原理是使用一个含有有限个分子并且有周期性边界条件的立方盒子,从该体系某一设定的位能模型出发,通过计算机模拟求解微元中全部分子的牛顿运动方程,记录它们在各个不同时刻的位置、速度和受力等,然后统计得到体系的各种热力学、结构和性质,也就是由体系粒子的微观性质求算其宏观性质。
进行分子动力学模拟的首要问题是要得到准确的原子之间的相互作用势函数。
知道原子间正确的相互作用势,从而必须知道相应的电子基态。
计算中根据以下基本假设:
(1)所有粒子的运动都遵循经典牛顿力学规律。
(2)粒子之间的相互作用满足叠加原理。
显然这两条忽略了量子效应和多体作用,与真实物理系统存在一定差别,仍然属于近似计算。
II)VMD作图:
VMD是一个强大的原子作图及动画演示软件,在运用C语言知识对上面求解计算后将得到一个zuobiao.xyz文件,将此文件直接拖进VMD原子作图工具便可以得到原子运动动画。
3方法
第一步:
即模型的设定,也就是势函数的选取。
势函数的研究和物理系统上对物质的描述研究息息相关。
最早是硬球势,即小于临界值时无穷大,大于等于临界值时为零。
常用的是LJ势函数,还有EAM势函数,不同的物质状态描述用不同的势函数。
第二步:
给定初始条件,也就是给定原子的空间位置和速度。
运动方程的求解需要知道粒子的初始位置和速度,不同的算法要求不同的初始条件。
如:
verlet算法需要两组坐标来启动计算,一组零时刻的坐标,一组是前进一个时间步的坐标或者一组零时刻的速度值。
第三步:
利用公式:
计算在第n+1步时所有粒子所处的空间位置。
计算在第n+1步时所有粒子的速度
,动能和速度标度因子:
计算每个时间步系统的总势能
,考察系统达到平衡所需要的时间。
第四步:
计算将速度乘以标度因子的值,并让该值作为下一次计算时,第n+1步粒子的速度。
第五步:
返回第三步,开始第n+2步的模拟计算。
重复几千次,一直到系统达到平衡。
第六步:
对系统进行位移加载,每次加载后使系统平衡1000步,让系统经历一个准静态加载的过程。
每次加载并让系统平衡后,计算宏观量。
4结果
1)程序运行部分结果用VMD作图的结果:
图
(1)
图
(2)
图(3)
2)程序运行结果,势能用matlab作图的结果:
I).利用matlab演示力随应变拉伸的变化:
图(4)
II).利用matlab演示势能随时间的变化:
图(5)
5讨论
图
(1)图
(2)图(3)是铝分子拉伸过程的一组图,分别表示拉伸前,中,后期。
可以发现在拉伸前铝分子的前两排和后两排的原子保持固定不变,分别向上下两方向拉伸,拉伸一次,趋衡一次并且得到每次原子在平衡状态下的坐标和势能。
由图(3)又可知,拉伸次数过多导致铝分子断裂。
由图(4)图(5)可发现力的变化曲线与势能的变化趋于一致,说明势能与力呈一定的相关性。
参考文献
[1]罗熙淳,梁迎春,董申.单晶铝纳米切削过程分子动力学模拟技术研究.中国机械工程,2000;
11(8):
860-8621.
[2]曹莉霞,王崇愚.α2Fe裂纹的分子动力学研究.钢铁研究总院功能材料所,北京 100081
[3]张媛媛,张文飞,王鹏。
不同温度下纳米单晶铜杆拉伸的分子动力学模拟
燕山大学建筑工程与力学学院,河北秦皇岛066004
附录:
采用VC++编程:
#include<
stdio.h>
math.h>
#include<
stdlib.h>
doublefen;
doublerij;
doubleA=0.1221,ebsl=1.316,r0=2.86378,wendu=10;
intN=300,m,s,t,n,i,j;
doublesm=0.0000861738569;
//k
doubleq=2.516,heli,xheli,eksj;
doubleli[300][2];
doublesheli1,sheli2,p=8.612;
doublea[300][2],v[300][2],lii[300][2],shineng[1000],shineng3;
doubledt=0.5,mal=27,dn,bdyz,kk;
voidlih();
//力的函数的声明
voidlih()//计算力
{
for(i=0;
i<
N;
i++)
{
shineng[i]=0;
li[i][0]=0.0;
li[i][1]=0.0;
}
{
fen=0.0;
for(j=0;
j<
j++)
{
if(j!
=i)
{
rij=pow((a[i][0]-a[j][0]),2)+pow((a[i][1]-a[j][1]),2);
rij=sqrt(rij);
fen+=exp(-2*q*(rij/r0-1));
}
xheli=sqrt(fen);
if(j!
shineng[i]=-ebsl*xheli;
shineng[i]+=A*exp(-p*(rij/r0-1));
for(j=0;
{
rij=pow((a[i][0]-a[j][0]),2)+pow((a[i][1]-a[j][1]),2);
rij=sqrt(rij);
sheli1=A*exp((-p)*(rij/r0-1))*(p/r0)/rij;
sheli2=ebsl*exp((-2)*q*(rij/r0-1))*(2*q)/r0/rij/xheli;
li[i][0]+=(sheli1-sheli2)*(a[i][0]-a[j][0]);
li[j][0]+=-(sheli1-sheli2)*(a[i][0]-a[j][0]);
li[i][1]+=(sheli1-sheli2)*(a[i][1]-a[j][1]);
li[j][1]+=-(sheli1-sheli2)*(a[i][1]-a[j][1]);
}
}
voidmain()//主函数
{
doublea1=sqrt
(2)*4.050/2,a2=sqrt(3/2)*4.050/2;
inti,j,k=0;
FILE*fp;
for(i=0;
20;
for(j=0;
15;
if(i%2==0)
a[k][0]=j*a1;
a[k][1]=a2*i;
else
a[k][0]=j*a1+a1/2;
k++;
}
fp=fopen("
zuobiao.dat"
"
w"
);
fprintf(fp,"
\n\n"
for(i=0;
k;
fprintf(fp,"
Al%.18f%.18f%.12f\n"
a[i][0],a[i][1],0.0);
fclose(fp);
//初始速度
v[i][0]=0.0;
v[i][1]=0.0;
lih();
//调用函数
//输出势能
fp=fopen("
shineng.dat"
AL势能=%.16f\n"
shineng[N]);
fclose(fp);
//趋衡过程
for(s=0;
s<
100;
s++)//趋衡100次
a[i][0]=a[i][0]+dt*v[i][0]+0.5*dt*dt*li[i][0]/mal;
a[i][1]=a[i][1]+dt*v[i][1]+0.5*dt*dt*li[i][1]/mal;
lii[i][0]=li[i][0];
lii[i][1]=li[i][1];
lih();
//重新计算力
//算速度
v[i][0]+=0.5*dt*(li[i][0]+lii[i][0])/mal;
v[i][1]+=0.5*dt*(li[i][1]+lii[i][1])/mal;
dn=0.0;
for(n=0;
n<
n++)
dn+=0.5*mal*(v[n][0]*v[n][0]+v[n][1]*v[n][1]);
//标度因子
bdyz=sqrt(0.5*(2*n-3)*wendu*sm/dn);
//对每个原子的速度进行标度
v[n][0]=v[n][0]*bdyz;
v[n][1]=v[n][1]*bdyz;
a"
for(i=0;
AL%.16f%.16f%.16f\n"
//拉伸阶段
for(t=0;
t<
500;
t++)//拉伸500次
for(i=N-30;
a[i][1]=a[i][1]+0.0005*15*a2;
30;
a[i][1]=a[i][1]-0.0005*15*a2;
50;
s++)//趋衡50次
for(i=30;
N-30;
li[i][0]=0.0;
shineng3=0.0;
xheli=sqrt(fen);
shineng3=-ebsl*xheli;
shineng3+=A*exp(-p*(rij/r0-1));
li[j][0]+=-(sheli1-sheli2)*(a[i][0]-a[j][0]);
li[i][1]+=(sheli1-sheli2)*(a[i][1]-a[j][1]);
li[j][1]+=-(sheli1-sheli2)*(a[i][1]-a[j][1]);
shineng[t]=shineng3;
//标度因子
//对每个原子的速度进行标度:
//输出坐标
printf("
第%d次"
t+1);
//输出势能
shineng[t]);
}//拉伸的括号
}//main函数的
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 金属 分子 动力学 模拟