From b3aedc986adbc1af7a867f0c721adfd7855bd1e4 Mon Sep 17 00:00:00 2001 From: sjplimp <sjplimp@f3b2605a-c512-4ea7-a41b-209d697bcdaa> Date: Fri, 12 Oct 2007 16:20:03 +0000 Subject: [PATCH] git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@1028 f3b2605a-c512-4ea7-a41b-209d697bcdaa --- tools/eam_generate/Al_Zhou.c | 111 ++++++++++++++++++++++++++++ tools/eam_generate/Cu_Mishin1.c | 123 ++++++++++++++++++++++++++++++++ tools/eam_generate/Cu_Zhou.c | 111 ++++++++++++++++++++++++++++ tools/eam_generate/README | 8 +++ tools/eam_generate/W_Zhou.c | 111 ++++++++++++++++++++++++++++ 5 files changed, 464 insertions(+) create mode 100644 tools/eam_generate/Al_Zhou.c create mode 100644 tools/eam_generate/Cu_Mishin1.c create mode 100644 tools/eam_generate/Cu_Zhou.c create mode 100644 tools/eam_generate/README create mode 100644 tools/eam_generate/W_Zhou.c diff --git a/tools/eam_generate/Al_Zhou.c b/tools/eam_generate/Al_Zhou.c new file mode 100644 index 0000000000..01b5f4515c --- /dev/null +++ b/tools/eam_generate/Al_Zhou.c @@ -0,0 +1,111 @@ +#include <stdio.h> +#include <math.h> +#include <stdlib.h> +#include <unistd.h> +#include <stdarg.h> +#include <string.h> + +#define EOK 0 +#define ERROR -1 + +#define re 2.886166 +#define fe 1.392302 +#define rhoe 20.226537 +#define alpha 6.942419 +#define beta 3.702623 +#define A 0.251519 +#define B 0.313394 +#define kappa 0.395132 +#define lambda 0.790264 +#define Fn0 -2.806783 +#define Fn1 -0.276173 +#define Fn2 0.893409 +#define Fn3 -1.637201 +#define F0 -2.83 +#define F1 0. +#define F2 0.929508 +#define F3 -0.68232 +#define eta 0.779208 +#define Fe -2.829437 + +double V (double r) { + return ( + ( A*exp (-alpha * (r/re-1.) ) ) / (1. + pow (r/re-kappa, 20.)) + - + ( B*exp (-beta * (r/re-1.) ) ) / (1. + pow (r/re-lambda, 20.)) + ); +} + +double rho (double r) { + return ( + (fe * exp (-beta * (r/re-1.))) / (1. + pow (r/re-lambda, 20.)) + ); +} + +double F (double rho_) { + double rhon = .85*rhoe, + rho0 = 1.15*rhoe; + if (rho_ < rhon) + return ( + Fn0 * pow (rho_/rhon-1., 0.) + + Fn1 * pow (rho_/rhon-1., 1.) + + Fn2 * pow (rho_/rhon-1., 2.) + + Fn3 * pow (rho_/rhon-1., 3.) + ); + else if (rhon <= rho_ && rho_ < rho0) + return ( + F0 * pow (rho_/rhoe-1., 0.) + + F1 * pow (rho_/rhoe-1., 1.) + + F2 * pow (rho_/rhoe-1., 2.) + + F3 * pow (rho_/rhoe-1., 3.) + ); + else if (rho0 <= rho_) + return ( + Fe*(1. - log( pow (rho_/rhoe, eta) ) ) * pow (rho_/rhoe, eta) + ); +} + +int main (void) { + int Nr = 10001; + double rmax = 4.041*2.5; //3.9860*4.; + double dr = rmax/(double)Nr; + int Nrho = 10001; + double rhomax = rho (0.); + double drho = rhomax/(double)Nrho; + + int atomic_number = 1; + double mass = 26.982; + double lattice_constant = 4.041; + char lattice_type[] = "FCC"; + + int i; + + char LAMMPSFilename[] = "Al_Zhou.eam"; + FILE *LAMMPSFile = fopen (LAMMPSFilename, "w"); + if (!LAMMPSFile) exit (ERROR); + + // Header for setfl format + fprintf (LAMMPSFile, \ + "#-> LAMMPS Potential File in DYNAMO 86 setfl Format <-#\n"\ + "# Zhou Al Acta mater(2001)49:4005\n"\ + "# Implemented by G. Ziegenhain (2007) gerolf@ziegenhain.com\n"\ + "%d Al\n"\ + "%d %20.20f %d %20.20f %20.20f\n"\ + "%d %20.20f %20.20f %s\n", + atomic_number, + Nrho, drho, Nr, dr, rmax, + atomic_number, mass, lattice_constant, lattice_type); + + // Embedding function + for (i = 0; i < Nrho; i++) + fprintf (LAMMPSFile, "%20.20f\n", F ((double)i*drho)); + // Density function + for (i = 0; i < Nr; i++) + fprintf (LAMMPSFile, "%20.20f\n", rho ((double)i*dr)); + // Pair potential + for (i = 0; i < Nr; i++) + fprintf (LAMMPSFile, "%20.20f\n", V ((double)i*dr) * (double)i*dr); + + fclose (LAMMPSFile); + return (EOK); +} diff --git a/tools/eam_generate/Cu_Mishin1.c b/tools/eam_generate/Cu_Mishin1.c new file mode 100644 index 0000000000..896593cbd1 --- /dev/null +++ b/tools/eam_generate/Cu_Mishin1.c @@ -0,0 +1,123 @@ +#include <stdio.h> +#include <math.h> +#include <stdlib.h> +#include <unistd.h> +#include <stdarg.h> +#include <string.h> + +#define rc 5.50679// A +#define h 0.50037// A +#define E1 2.01458*1e2//eV +#define E2 6.59288*1e-3//eV +#define r01 0.83591//A +#define r02 4.46867//A +#define a1 2.97758//1/A +#define a2 1.54927//1/A +#define d 0.86225*1e-2//A +#define rs1 2.24//A +#define rs2 1.8//A +#define rs3 1.2//A +#define S1 4.//eV/A^4 +#define S2 40.//eV/A^4 +#define S3 1.15*1e3//eV/A^4 +#define a 3.80362 +#define r03 -2.19885//A +#define r04 -2.61984*1e2//A +#define b1 0.17394//A^-2 +#define b2 5.35661*1e2//1/A +#define F0 -2.28235//eV +#define F2 1.35535//eV +#define q1 -1.27775//eV +#define q2 -0.86074//eV +#define q3 1.78804//eV +#define q4 2.97571//eV +#define Q1 0.4 +#define Q2 0.3 + +#define EOK 0 +#define ERROR -1 + + +double M (double r, double r0, double alpha) { + return exp (-2.*alpha*(r-r0)) - 2.*exp (-alpha*(r-r0)); +} + +double psi (double x) { + if (x >= 0.) return 0.; + else return pow (x, 4.) / (1.+pow (x, 4.) ); +} + +double H (double x) { + if (x >= 0) return 1.; + else return 0.; +} + +double V (double r) { + return ( E1*M(r,r01,a1) + E2*M(r,r02,a2) + d ) * psi ( (r-rc)/h ) - + H (rs1-r)*S1*pow (rs1-r, 4.) - + H (rs2-r)*S2*pow (rs2-r, 4.) - + H (rs3-r)*S3*pow (rs3-r, 4.); +} + +double rho (double r) { + return ( a*exp (-b1*pow (r-r03, 2.)) + exp (-b2*(r-r04)) ) * psi ( (r-rc)/h ); +} + +double F (double rho_) { + if (rho_ < 1.) + return F0 + .5*F2*pow (rho_-1., 2.) + + q1*pow (rho_-1.,3.) + + q2*pow (rho_-1.,4.) + + q3*pow (rho_-1.,5.) + + q4*pow (rho_-1.,6.); + else if (rho_ > 1.) + return ( F0 + .5*F2*pow (rho_-1., 2.) + q1*pow (rho_-1., 3.) + Q1*pow (rho_-1., 4.)) / + ( 1. + Q2*pow (rho_-1., 3.) ); + else exit (ERROR); + +} + +int main (void) { + int Nr = 10001; + double rmax = 9.; + double dr = rmax/(double)Nr; + int Nrho = 10001; + double rhomax = rho (0.); + double drho = rhomax/(double)Nrho; + + int atomic_number = 1; + double mass = 63.55; + double lattice_constant = 3.615; + char lattice_type[] = "FCC"; + + int i; + + char LAMMPSFilename[] = "Cu_Mishin1.eam"; + FILE *LAMMPSFile = fopen (LAMMPSFilename, "w"); + if (!LAMMPSFile) exit (ERROR); + + // Header for setfl format + fprintf (LAMMPSFile, \ + "#-> LAMMPS Potential File in DYNAMO 86 setfl Format <-#\n"\ + "# Mishin Cu EAM1 PRB(2001)63:224106\n"\ + "# Implemented by G. Ziegenhain (2007) gerolf@ziegenhain.com \n"\ + "%d Cu\n"\ + "%d %20.20f %d %20.20f %20.20f\n"\ + "%d %20.20f %20.20f %s\n", + atomic_number, + Nrho, drho, Nr, dr, rc, + atomic_number, mass, lattice_constant, lattice_type); + + // Embedding function + for (i = 0; i < Nrho; i++) + fprintf (LAMMPSFile, "%20.20f\n", F ((double)i*drho)); + // Density function + for (i = 0; i < Nr; i++) + fprintf (LAMMPSFile, "%20.20f\n", rho ((double)i*dr)); + // Pair potential + for (i = 0; i < Nr; i++) + fprintf (LAMMPSFile, "%20.20f\n", V ((double)i*dr) * (double)i*dr); + + fclose (LAMMPSFile); + return (EOK); +} diff --git a/tools/eam_generate/Cu_Zhou.c b/tools/eam_generate/Cu_Zhou.c new file mode 100644 index 0000000000..aa687bc5b7 --- /dev/null +++ b/tools/eam_generate/Cu_Zhou.c @@ -0,0 +1,111 @@ +#include <stdio.h> +#include <math.h> +#include <stdlib.h> +#include <unistd.h> +#include <stdarg.h> +#include <string.h> + +#define EOK 0 +#define ERROR -1 + +#define re 2.556162 +#define fe 1.554485 +#define rhoe 22.150141 +#define alpha 7.669911 +#define beta 4.090619 +#define A 0.327584 +#define B 0.468735 +#define kappa 0.431307 +#define lambda 0.86214 +#define Fn0 -2.17649 +#define Fn1 -0.140035 +#define Fn2 0.285621 +#define Fn3 -1.750834 +#define F0 -2.19 +#define F1 0. +#define F2 0.702991 +#define F3 0.683705 +#define eta 0.921150 +#define Fe -2.191675 + +double V (double r) { + return ( + ( A*exp (-alpha * (r/re-1.) ) ) / (1. + pow (r/re-kappa, 20.)) + - + ( B*exp (-beta * (r/re-1.) ) ) / (1. + pow (r/re-lambda, 20.)) + ); +} + +double rho (double r) { + return ( + (fe * exp (-beta * (r/re-1.))) / (1. + pow (r/re-lambda, 20.)) + ); +} + +double F (double rho_) { + double rhon = .85*rhoe, + rho0 = 1.15*rhoe; + if (rho_ < rhon) + return ( + Fn0 * pow (rho_/rhon-1., 0.) + + Fn1 * pow (rho_/rhon-1., 1.) + + Fn2 * pow (rho_/rhon-1., 2.) + + Fn3 * pow (rho_/rhon-1., 3.) + ); + else if (rhon <= rho_ && rho_ < rho0) + return ( + F0 * pow (rho_/rhoe-1., 0.) + + F1 * pow (rho_/rhoe-1., 1.) + + F2 * pow (rho_/rhoe-1., 2.) + + F3 * pow (rho_/rhoe-1., 3.) + ); + else if (rho0 <= rho_) + return ( + Fe*(1. - log ( pow (rho_/rhoe, eta) )) * pow (rho_/rhoe, eta) + ); +} + +int main (void) { + int Nr = 10001; + double rmax = 3.615*2.5; + double dr = rmax/(double)Nr; + int Nrho = 10001; + double rhomax = rho (0.); + double drho = rhomax/(double)Nrho; + + int atomic_number = 1; + double mass = 63.55; + double lattice_constant = 3.615; + char lattice_type[] = "FCC"; + + int i; + + char LAMMPSFilename[] = "Cu_Zhou.eam"; + FILE *LAMMPSFile = fopen (LAMMPSFilename, "w"); + if (!LAMMPSFile) exit (ERROR); + + // Header for setfl format + fprintf (LAMMPSFile, \ + "#-> LAMMPS Potential File in DYNAMO 86 setfl Format <-#\n"\ + "# Zhou Cu Acta mater(2001)49:4005\n"\ + "# Implemented by G. Ziegenhain (2007) gerolf@ziegenhain.com\n"\ + "%d Cu\n"\ + "%d %20.20f %d %20.20f %20.20f\n"\ + "%d %20.20f %20.20f %s\n", + atomic_number, + Nrho, drho, Nr, dr, rmax, + atomic_number, mass, lattice_constant, lattice_type); + + // Embedding function + for (i = 0; i < Nrho; i++) + fprintf (LAMMPSFile, "%20.20f\n", F ((double)i*drho)); + // Density function + for (i = 0; i < Nr; i++) + fprintf (LAMMPSFile, "%20.20f\n", rho ((double)i*dr)); + // Pair potential + for (i = 0; i < Nr; i++) + fprintf (LAMMPSFile, "%20.20f\n", V ((double)i*dr) * (double)i*dr); + + fclose (LAMMPSFile); + return (EOK); +} diff --git a/tools/eam_generate/README b/tools/eam_generate/README new file mode 100644 index 0000000000..4ecce780dc --- /dev/null +++ b/tools/eam_generate/README @@ -0,0 +1,8 @@ +These are one-file C programs that generate tabulated embedded atom +method (EAM) potentials from analytic formulas. They each write out a +file in the DYNAMO setfl format which can be input by LAMMPS via the +pair_style eam/alloy command. + +The source files and potentials were provided by Gerolf Ziegenhain +(gerolf at ziegenhain.com). You could modify these programs to write +out your own EAM potentials in LAMMPS-compatible format. diff --git a/tools/eam_generate/W_Zhou.c b/tools/eam_generate/W_Zhou.c new file mode 100644 index 0000000000..c5bac924a7 --- /dev/null +++ b/tools/eam_generate/W_Zhou.c @@ -0,0 +1,111 @@ +#include <stdio.h> +#include <math.h> +#include <stdlib.h> +#include <unistd.h> +#include <stdarg.h> +#include <string.h> + +#define EOK 0 +#define ERROR -1 + +#define re 2.74084 +#define fe 3.48734 +#define rhoe 37.234847 +#define alpha 8.900114 +#define beta 4.746728 +#define A 0.882435 +#define B 1.394592 +#define kappa 0.139209 +#define lambda 0.278417 +#define Fn0 -4.946281 +#define Fn1 -0.148818 +#define Fn2 0.365057 +#define Fn3 -4.432406 +#define F0 -4.96 +#define F1 0. +#define F2 0.661935 +#define F3 0.348147 +#define eta -0.582714 +#define Fe -4.961306 + +double V (double r) { + return ( + ( A*exp (-alpha * (r/re-1.) ) ) / (1. + pow (r/re-kappa, 20.)) + - + ( B*exp (-beta * (r/re-1.) ) ) / (1. + pow (r/re-lambda, 20.)) + ); +} + +double rho (double r) { + return ( + (fe * exp (-beta * (r/re-1.))) / (1. + pow (r/re-lambda, 20.)) + ); +} + +double F (double rho_) { + double rhon = .85*rhoe, + rho0 = 1.15*rhoe; + if (rho_ < rhon) + return ( + Fn0 * pow (rho_/rhon-1., 0.) + + Fn1 * pow (rho_/rhon-1., 1.) + + Fn2 * pow (rho_/rhon-1., 2.) + + Fn3 * pow (rho_/rhon-1., 3.) + ); + else if (rhon <= rho_ && rho_ < rho0) + return ( + F0 * pow (rho_/rhoe-1., 0.) + + F1 * pow (rho_/rhoe-1., 1.) + + F2 * pow (rho_/rhoe-1., 2.) + + F3 * pow (rho_/rhoe-1., 3.) + ); + else if (rho0 <= rho_) + return ( + Fe*(1. - log( pow (rho_/rhoe, eta) ) ) * pow (rho_/rhoe, eta) + ); +} + +int main (void) { + int Nr = 10001; + double rmax = 2.5*3.157; + double dr = rmax/(double)Nr; + int Nrho = 10001; + double rhomax = rho (0.); + double drho = rhomax/(double)Nrho; + + int atomic_number = 1; + double mass = 183.84; + double lattice_constant = 3.157; + char lattice_type[] = "BCC"; + + int i; + + char LAMMPSFilename[] = "W_Zhou.eam"; + FILE *LAMMPSFile = fopen (LAMMPSFilename, "w"); + if (!LAMMPSFile) exit (ERROR); + + // Header for setfl format + fprintf (LAMMPSFile, \ + "#-> LAMMPS Potential File in DYNAMO 86 setfl Format <-#\n"\ + "# Zhou W Acta mater(2001)49:4005\n"\ + "# Implemented by G. Ziegenhain (2007) gerolf@ziegenhain.com\n"\ + "%d W\n"\ + "%d %20.20f %d %20.20f %20.20f\n"\ + "%d %20.20f %20.20f %s\n", + atomic_number, + Nrho, drho, Nr, dr, rmax, + atomic_number, mass, lattice_constant, lattice_type); + + // Embedding function + for (i = 0; i < Nrho; i++) + fprintf (LAMMPSFile, "%20.20f\n", F ((double)i*drho)); + // Density function + for (i = 0; i < Nr; i++) + fprintf (LAMMPSFile, "%20.20f\n", rho ((double)i*dr)); + // Pair potential + for (i = 0; i < Nr; i++) + fprintf (LAMMPSFile, "%20.20f\n", V ((double)i*dr) * (double)i*dr); + + fclose (LAMMPSFile); + return (EOK); +} -- GitLab