| 查看: 2468 | 回复: 4 | |||
牛排小黑娇新虫 (初入文坛)
|
[求助]
势函数解析问题
|
|
我想做Mo-Ag的合金扩散模拟，但是并没有现成的势函数文件，我查阅文献查到了有人拟合过这个合金的势函数。 势函数形式是FS势格式和EAM势差不多的，我看lammps手册这个势函数可以通过一个DYNAMO的程序解析出来，不知道您是不是也做过这个势函数的解析文件，我想让您帮我看一下我的程序是不 是有问题，我之前弄出来的模拟结果只能运行一行。希望尚老师您能抽时间帮我看下，万分感谢！ #include <stdio.h> #include <math.h> #include <stdlib.h> #include <unistd.h> #include <string> #include <stdarg.h> #include <string.h> #include "vector" using namespace std; #define EOK 0 #define ERROR -1 double phi(double r,double cc,double c0,double c1,double c2,double c3,double c4) { double V_r; if(r<=cc) { V_r= (r-cc)*(r-cc)*(c0+c1*r+c2*r*r+c3*pow(r,3)+c4*pow(r,4)); } else V_r=0.0; return (V_r); } double rho(double r,double d,double beta) { double rho_r; if(r<=d) rho_r= (r-d)*(r-d) + beta*beta*pow((r-d),4); else rho_r=0.0; return (rho_r); } double Frho(double rho,double AA) { double F_rho; F_rho = -1.0*AA*pow(rho,0.5); return (F_rho); } int main (void) { int Nr = 10000; double rmax =9.8; 这边取值的标准是？ double dr = rmax/(double)Nr; int Nrho = 10000; double rhomax = 45.5; double drho = rhomax/(double)Nrho; int i; string LAMMPSFilename="AgMo.eam.fs"; FILE *LAMMPSFile = fopen ((char*)LAMMPSFilename.data(), "w"  ;if (!LAMMPSFile) exit (ERROR); // Header for setfl format fprintf (LAMMPSFile, \ "#-> LAMMPS Potential File in DYNAMO 86 setfl Format <-#\n"\ "# Title :Embedded-Atom Method by Arthur F. Woter 1994\n"\ "# Implemented by Nize(2015) nize19890627@gmail.com\n"\ "2 Ag Mo\n"\ "%d %20.20f %d %20.20f %20.20f\n",Nrho,drho,Nr,dr,rmax); fprintf (LAMMPSFile,"47 107.8682 4.0853 FCC\n" ;// Embedding function and density function for Ag for(i=0;i<Nrho;i++) fprintf (LAMMPSFile,"%20.20e\n",Frho((double)i*drho,0.325514)); // Density function for Ag for (i = 0; i < Nr; i++) fprintf (LAMMPSFile, "%20.20e\n",rho((double)i*dr,4.41,-1.293394)); fprintf (LAMMPSFile,"42 95.94 3.1472 BCC\n" ;// Embedding function and density function for Mo for(i=0;i<Nrho;i++) fprintf (LAMMPSFile,"%20.20e\n",Frho((double)i*drho,1.848648)); // Density function for Mo for (i = 0; i < Nr; i++) fprintf (LAMMPSFile, "%20.20e\n",rho((double)i*dr,4.14,0.0)); // Pair potential for Ag-Ag for (i = 0; i < Nr; i++) fprintf (LAMMPSFile, "%20.20e\n", (double)i*dr*phi((double)i*dr,4.76,10.681200,-12.045170,5.203072,-1.013304,0.0742308)); // Pair potential for Ag-Mo for (i = 0; i < Nr; i++) fprintf (LAMMPSFile, "%20.20e\n", (double)i*dr*phi((double)i*dr,4.50,44.406810,-45.490260,15.605110,-1.793704,0.0)); // Pair potential for Mo-Mo for (i = 0; i < Nr; i++) fprintf (LAMMPSFile, "%20.20e\n", (double)i*dr*phi((double)i*dr,3.2572,47.980660,-34.099240,5.832293,0.101749,0.0203934)); fclose (LAMMPSFile); // system("pause" ;return (EOK); 知道的朋友帮我看看。急求 谢谢啦  QQ图片20150626105650.png  QQ图片20150626105951.png  QQ图片20150702114556.jpg |
回复此楼» 收录本帖的淘帖专辑推荐
 分子动力学模拟 |
» 猜你喜欢
综述论文作为代表作会不会影响评审专家的印象分？
已经有12人回复
-
2019年青年基金涵评意见，大家看看几个A，几个B?
已经有4人回复
-
应该是93bebmhtak前后十一个字符比较关键
已经有29人回复
-
奇怪，两个人的filecode固定段从头到尾一模一样
已经有9人回复
-
是这周出结果还是下周出结果？
已经有3人回复
-
小木虫上这么多卖论文的，真有人买论文么？感觉没必要啊
已经有9人回复
-
有时候，自然基金真的不能太认真 （我的申报经验）
已经有4人回复
-
帮忙看看fileCode
已经有13人回复
-
为什么网上很多人说本周 12号出结果
已经有7人回复
-
读博的好处
已经有3人回复
» 本主题相关价值贴推荐，对您同样有帮助:
- 存在这样的函数吗
已经有11人回复
- 有没有人知道塑性功怎么计算
已经有10人回复
- 关于函数共轭的问题，求助
已经有6人回复
- 急求解答，关于力场的电荷参数！！！
已经有4人回复
- 寻找拟合函数（由0增长，趋近于常值，并且0处的导数值为0等诸多条件）
已经有7人回复
- 分子中电荷变化种种
已经有82人回复
- 至少500金币求解答一个问题-磷氧键（主要是P=O双键，还有P-O单键）键能多少？
已经有12人回复
- 请教密度泛函理论中屏蔽库伦势的问题
已经有14人回复
- Gulp 拟合势函数 初值设置的问题
已经有12人回复
- 求概率密度函数的解析式
已经有9人回复
- 【求助】第一性原理与密度泛函理论的关系是什么
已经有15人回复
- 【求助】RDF分析-lammps
已经有12人回复
- 【交流】归一化结构因子ESEL的上限和下限值怎么设置？
已经有3人回复
|
本帖内容被屏蔽 |
2楼2015-09-14 18:51:44
3楼2015-09-24 17:08:17
weihui333
铜虫 (初入文坛)
- 应助: 0 (幼儿园)
- 金币: 134.1
- 帖子: 17
- 在线: 30.8小时
- 虫号: 2731250
- 注册: 2013-10-17
- 性别: GG
- 专业: 化学反应工程
4楼2015-09-26 08:01:55
5楼2016-09-11 14:02:20
