# How to optimize zinc oxide bulk structure at certain pressure with LAMMPS

I am trying to simulate the Zinc oxide bulk structure (periodic) at constant pressure(0.0 ~ 20.0 GPa) to observe phase transition with LAMMPS.

However, I approached a problem that the fluctuating pressure in each step of optimisation and importantly the final step of the optimisation does not have the value of pressure that I set.

*5th-edit: the problem is solved *

The last step of the pressure is not the exact the pressure that I set for which presumably a form of showing residual force of minimisation. Basically, it went down to the basin but not exactly the lowest point of the basin with a certain tolerance.

Primarily the calculation went wrong because I have used wrongly formatted structure data in the read_data command. The wrong data file was produced by three different software... which was very unlucky...

I sincerely apologise for the question as I made many people confused.

Here is my input:

units       metal
atom_style  charge
dimension   3
boundary    p p p

#box tilt large                          # allow lammps to have large tilt facotr (xy yz xz)
#replicate       3 3 3

pair_style  comb3 polar_off
pair_coeff  * * ./lammps_ip Zn O

neighbor        2.0 bin
neigh_modify every 1 delay 0 check yes

set type 1    charge  2.0  # Zn
set type 2    charge -2.0  # O

##############################################################
#                   type of output data                      #
##############################################################

fix             1 all qeq/comb 2 0.00003 file fq.out
fix             2 all box/relax iso 0.0                 # Pressure = unit [bars], 1 GPa = 10,000 bars, Constant pressure

thermo          1
thermo_style    custom step press temp enthalpy etotal pe ke evdwl ecoul cella cellb cellc cellalpha cellbeta cellgamma

#min_style cg
minimize 0 0.0 10000 10000

write_data data.*        #produce file for "read_data" command
print " All done "


The corresponding output data is: Step Press c_1 v_press Temp v_etot v_relax TotEng PotEng KinEng

    1157   -186.71231   -186.71231   -186.71231            0  -2.9709191            0   -2.9709191   -2.9709191            0
1158    -128.7401    -128.7401    -128.7401            0   -2.9709192            0   -2.9709192   -2.9709192            0
1159   -124.22057   -124.22057   -124.22057            0   -2.9709192            0   -2.9709192   -2.9709192            0
Loop time of 464.543 on 1 procs for 1159 steps with 216 atoms

Minimization stats:
Stopping criterion = linesearch alpha is zero
Energy initial, next-to-last, final =
-2.25448773082     -2.97091924357     -2.97091924466
Force two-norm initial, final = 351.736 2.00988
Force max component initial, final = 342.312 0.805865
Final line search alpha, max atom move = 4.67711e-07 3.76912e-07
Iterations, force evaluations = 1159 1175

• What would you like to achieve, perfect non-fluctuating pressure? You write yourself that the pressure will fluctuate during an NPT simulations. This is usually acceptable and the fluctuations typically go down when the system is properly equilibrated. May 26 '20 at 13:13
• Nice first question. Welcome to the site and we hope to see more of you! May 26 '20 at 13:15
• @Fabian Thank you for the comment. Yes, I understand that NPT simulation will provides pressure fluctuation, but what I trying to do is just simple global optimisation with constant pressure. Is this considered as NPT simulations? May 26 '20 at 13:34
• @CodyAldaz I will take your advice in next posts. Apologies for the inconveniences. May 27 '20 at 15:36
• In your last edit you have said that you solved the problem. I think instead of putting the answer into the question, you should write an "answer" to your own question. Jun 12 '20 at 15:08

However, I am not clear which command have to be used to optimise the structure at a certain pressure to have the constant pressure at the final step.

For this you need to perform a time integration from the Isobaric Isothermal ($$NPT$$) ensemble. Essentially performing these dynamics, you equilibrate the systems at certain $$T$$ and $$P$$. In LAMMPS documentation you can find the necessary commands for $$NPT$$ simulation. A very generic example is given below.

fix 1 all npt temp 300.0 300.0 $$(100.0*dt) iso 0.0 0.0$$(1000.0*dt)
run 1000000


After equilibration, plot $$E(timestep)$$. If the slope becomes zero, your system is properly equilibrated at that $$T$$ and $$P$$, if you plot the $$P(timestep)$$ after equilibration, it should always fluctuate around Final $$P$$ as per your input.

• Thank you for the advice. However, what I want to do is geometry optimisation (H=E + pv) not MD simulation. (minimisation enthalpy) Jun 10 '20 at 9:45
• Is there any specific reason to run optimization only via LAMMPS "minimize" command? May I ask what is the purpose of your work? "Optimization" doesn't need to adhere to "minimize" command just because "minimize" command is designed to minimize something. Jun 12 '20 at 2:33
• Sometimes I also run very low-temperature (like 0.1K) NPT MD or NVT MD to optimize the system energy and geometry, or sometimes I use minimize command and very short NPT or NVT or NVE MD to stabilize the system. The choice of procedure usually depends on my goal or target of the simulation. Jun 12 '20 at 2:34
• Also, if you want to really "optimize" something only with "minimize" command, then you can test other minimize algorithms with different criteria conditions. Jun 12 '20 at 2:35
• You may know about this, but don't forget to monitor system PE, not only pressure, to judge if your system is optimized or not. Jun 12 '20 at 2:58