LAMMPS WWW Site - LAMMPS Documentation - LAMMPS Mailing List Archives
Re: [lammps-users] Different potential energy calculation results
[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

Re: [lammps-users] Different potential energy calculation results


From: Xiao Jia <xiaoj@...1508...>
Date: Tue, 19 Dec 2017 17:42:08 -0500

Dear Dr. Kohlmeyer,

In the code of pair_eam.cpp, the overall potential energy seems to be also summed from pairwise potential energy as 

recip = 1.0/r;
        phi = z2*recip;
        phip = z2p*recip - phi*recip;
        psip = fp[i]*rhojp + fp[j]*rhoip + phip;
        fpair = -psip*recip;
        f[i][0] += delx*fpair;
        f[i][1] += dely*fpair;
        f[i][2] += delz*fpair;
        if (newton_pair || j < nlocal) {
          f[j][0] -= delx*fpair;
          f[j][1] -= dely*fpair;
          f[j][2] -= delz*fpair;
        }
        if (eflag) evdwl = phi;
        if (evflag) ev_tally(i,j,nlocal,newton_pair,
                             evdwl,0.0,fpair,delx,dely,delz);

the value of phi should be the same output from single () function. 

This part seems to be similar with pairwise potential like Morse

        if (eflag) evdwl = phi;
        if (evflag) ev_tally(i,j,nlocal,newton_pair,
                             evdwl,0.0,fpair,delx,dely,delz);

So the I still have the curiosity about the difference in calculation of potential. 
Thanks!

On Tue, Dec 19, 2017 at 5:01 PM, Xiao Jia <xiaoj@...1508...> wrote:
Dear Dr. Kohlmeyer,

I am not sure about the results from compute group/group.
If that is the case, does it mean that the single() function cannot be used directly for the calculation the atom potential in EAM potential.
I feel by referring this function, the coeff of fpair is calculated correctly. 
Is there other way to refer the value of single pair value?

Thanks!


On Tue, Dec 19, 2017 at 4:55 PM, Axel Kohlmeyer <akohlmey@...24...> wrote:


On Tue, Dec 19, 2017 at 4:36 PM, Xiao Jia <xiaoj@...1508...> wrote:
Dear Dr. Kohlmeyer,

Sorry for the missing information. 
I am using the latest stable version of lammps, 11Aug2017. 
I tried another case with Morse potential, which generates consistent potential value to the log file.
While for EAM potential, it does not. 


​[...]​

I am not sure why different results can be generated for different potential, Morse and EAM,
and is there any difference to refer the same function single() by using different potentials.

​there is a very subtle difference between morse and eam. the former is a true pair-wise additive potential, the latter is not.
the PairEAM::single() function is dependent on position dependent properties in data structures, that are filled with suitable data during the original force computation.​

​do you get the same "wrong" number, when you are using compute group/group (which should be doing something similar to what you are doing with your modification?​

​axel.​

Thanks!


On Tue, Dec 19, 2017 at 9:45 AM, Axel Kohlmeyer <akohlmey@...29....> wrote:


On Tue, Dec 19, 2017 at 3:33 AM, Xiao Jia <xiaoj@...1508...> wrote:
Dear LAMMPS users,

I am referring the function single () defined in pair.cpp to calculate the pairwise force coefficient fpair and the potential energy.
eng=pair->single (i, j, itype, jtype, rsq, factor_coul, factor_lj, fpair);
The current potential I am using is EAM potential.

The weird observation I have is that my summation of the potential energy (eng value of each pair) from all the pairs does not match with the output from log file. The value I get is -21936114 and the value from log file is -63005273, the difference is close to 3 times. 

The calculation of pressure based on fpair seems to be consistent with the value from log file, but the difference in potential energy is hard to understand. 

​your e-mail is missing two important pieces of information.

1) which LAMMPS version you are using exactly. what you are reporting sounds a lot like a bug that was fixed not so long ago.
2) an example input deck that allows to quickly reproduce your issue (in case this is a real bug and still exists with the latest patch release).

axel.​


 

I will appreciate it so much if anyone can give me any idea about this issue.
Thanks!


------------------------------------------------------------------------------
Check out the vibrant tech community on one of the world's most
engaging tech sites, Slashdot.org! http://sdm.link/slashdot
_______________________________________________
lammps-users mailing list
lammps-users@...655....net
https://lists.sourceforge.net/lists/listinfo/lammps-users




--
Dr. Axel Kohlmeyer  akohlmey@...24...  http://goo.gl/1wk0
College of Science & Technology, Temple University, Philadelphia PA, USA
International Centre for Theoretical Physics, Trieste. Italy.




--
Dr. Axel Kohlmeyer  akohlmey@...24...  http://goo.gl/1wk0
College of Science & Technology, Temple University, Philadelphia PA, USA
International Centre for Theoretical Physics, Trieste. Italy.