[lammps-users] Z(r) for EAM potentials?

Dear All,

My name is Alberto Fraile and I am a PhD student in the University of Madrid, in the Nuclear Fusion Department. I am working mainly on MD simulation with LAMMPS.

We are trying to implement an EAM potential in our LAMMPS scripts but we found a couple of problems.

We are following this paper [Zhou et al, Acta mater. 49 (2001) 4005?4015] in order to create a potential file for Pb in EAM format.
We have modified the EAM example Al_Zhou.c, (in LAMMPS/tools/eam_generate/) just changing the parameters for Al for those of Pb given in the previous reference.

Now the program aparently creates a EAM potential in setfl format. How can we made this file to be read by LAMMPS?

I have tried to change this to funcfl format, changing the header etc ... but here is the problem.
Changing the header is easy: 1 line comment instead of 3 etc.
But after that the program creates three arrays: functions F(r), rho(r) and phi(r).
According to LAMMPS manual (page 439) funcfl format expects three functions, F(rho), Z(r), being Z the effective charge density, and rho(r).
Then I have to change too the order in the ?print section?

This is a piece of the program Al_Zhou.c to explain the problem better:

// 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) [email protected]...\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);

(here V=phi(r), the pair potential)

The question now is simple:
How can I calculate Z(r)? Or should I write Z^2(r)?

According with the file ?pair_EAM.c? (in LAMMPS/src/MANYBODY) I think again that LAMMPS expects the aforementioned functions F, rho and Z (if format is funcfl), (and here is a difference with the mannual in the order!) and F, rho and Z2 (if format is setfl).

Below, (in pair_EAM.c) I can see (in a comment) that Z2 is phi*r, but I can not see what is Z. Is it Z the same as Z2? I think that Z2 is Z^2.

I have read for example in this paper [S. M. Foiles PRB Vol 32 N6 (1985)] that the effective electron charge function, Z(r), is related with the pair potential by the equation phi=Z^2(r)/r, then seems that Z(r) is just sqrt(phi*r). But it is not possible to do the sqrt of phi*r because our phi (calculated as explained before) is negative sometimes (and so phi*r).

In short:

Which is the correct order to F, rho and Z in the funcfl format?
How can I calculate Z(r)?

Thank you very much for your help
Alberto Fraile

The LAMMPS doc page for doc/pair_eam.html explains
the format of these files in great detail. If you think
there is an error in one of those pages, please let us
know.

You are correct that some of the files require Z(r), and others
Z^2, as I recall. Also, I think it is correct that one of
the formats is not general enough to allow for negative
functions due to the sqrt issue. But the other formats are,
so you can always use one of those. E.g. the Finnis/Sinclair
formulation is the most general, as I recall.

Steve