I am trying to implement a pair style for core/shell model for modelling ionic materials. My pair style consist of a electrostatic term (Wolf Summation), Buckingham potential for calculating the energy effects of shells, a harmonic spring potential and a harmonic three body term, (Lattice dynamics of Forsterite, Price et al., 1987)
I validated the implementation in comparison with the core/shell model existing in GULP for the material Magnesium Oxide(MgO) and forsterite.
I verified the implementation of two body terms with the material MgO where I could reproduce the elastic constants from the script given in LAMMPS examples which is using constant pressure minimisation.
Whereas during the validation with harmonic three body term, I could reproduce Bulk modulus and Generalised Stacking Fault curves for Forsterite, the force and energy minimization at constant volume matches between LAMMPS and GULP.

And when I calculate the elastic constant at constant volume(O.Beckstein et al., 2001) I could match the results, where as I could not use the ELASTIC script in LAMMPS example for three body term.

I am facing problem with pressure associated with three body term, can someone help with it.

I am trying to implement a pair style for core/shell model for modelling
ionic materials. My pair style consist of a electrostatic term (Wolf
Summation), Buckingham potential for calculating the energy effects of
shells, a harmonic spring potential and a harmonic three body term, (Lattice
dynamics of Forsterite, Price et al., 1987)
I validated the implementation in comparison with the core/shell model
existing in GULP for the material Magnesium Oxide(MgO) and forsterite.
I verified the implementation of two body terms with the material MgO where
I could reproduce the elastic constants from the script given in LAMMPS
examples which is using constant pressure minimisation.
Whereas during the validation with harmonic three body term, I could
reproduce Bulk modulus and Generalised Stacking Fault curves for Forsterite,
the force and energy minimization at constant volume matches between LAMMPS
and GULP.
And when I calculate the elastic constant at constant volume(O.Beckstein et
al., 2001) I could match the results, where as I could not use the ELASTIC
script in LAMMPS example for three body term.

I am facing problem with pressure associated with three body term, can
someone help with it.

it is not clear what you are asking here. there are a lot of details,
but when you get to your specific problem, you are extremely vague. it
is not even clear why somebody would want to help you, since your
problem seems to be entirely an issue of the feature that you are
implementing.

To be simple, I implemented a pairstyle core/shell for core/shell model propesd by Dick and Overhauser (1958),

Using the pair style I checked the energy and forces for materials MgO and forsterite which are in agreement with the core/shell model in GULP.

By using constant volume mimimisation method I calculated the bulk modulus, elastic costant and GSF excess energies for both of my materials, and I could accurately reproduce.

I tried to fit in my pair style to the ELASTIC script given in examples of LAMMPS to calculate elastic constants at constant pressure, which works for MgO which does not have a three body term, whereas when I want to do it for forsterite I could not use the ELASTIC script. I face problem with pressure of the system arising from the harmonic 3 body term I implemented.

I used full neighboring list in my implementation.

To be simple, I implemented a pairstyle core/shell for core/shell model
propesd by Dick and Overhauser (1958),

Using the pair style I checked the energy and forces for materials MgO and
forsterite which are in agreement with the core/shell model in GULP.

By using constant volume mimimisation method I calculated the bulk modulus,
elastic costant and GSF excess energies for both of my materials, and I
could accurately reproduce.

I tried to fit in my pair style to the ELASTIC script given in examples of
LAMMPS to calculate elastic constants at constant pressure, which works for
MgO which does not have a three body term, whereas when I want to do it for
forsterite I could not use the ELASTIC script. I face problem with pressure
of the system arising from the harmonic 3 body term I implemented.

so you have to correct your implementation. i still don't see how this
is of anybody's concern but yours.
the ELASTIC script works fine as it is with manybody terms. in fact,
the very example that is given uses a Stillinger-Weber potential that
has 3-body terms.

For any 3-body interaction, if you take
a set of 3 coords where the 3 particles
are close to each other (may require adjusting
by periodic image of 1 or 2 of the atoms),
then the virial contribution of the 3-body term
is simple sum-over-i (Fi dot Ri), where
i = 1,2,3, and the Fi are the forces
on each atom due to the 3-body term.

To see if periodic images are messing you up,
I would first see if you can match the virial/pressure
for a small system with no periodicity (cluster).
If you are getting the forces right, you can’t be
far off.

We’d like to add a core/shell model to LAMMPS,
so please contribute your code if you work everything
out.

Thank you, for your suggestion. I identified the problem which was during the virial stress calculation in fdotr. Now I could have closer convergence to the model exsisting in GULP during constant pressure minimisation process.