Thanks for your instructions. I changed the pppm to ewald type and got the following errors:
switching the kspace style doesn’t do anything to solve your issue.
ERROR on proc 1: Bond atoms 1135 1136 missing on proc 1 at step 294 (…/ntopo_bond_all.cpp:63)
I attached the lammps input file for convenience. Could you please suggest a way to mitigate this problem during merging two boxes, Axel?
i already explained the reason for the problem. whenever you scale the box without taking precautions, you get bogus coordinates.
for individual atoms, you can explicitly set the image flags to zero first. but for molecules, you do not have that option.
when you do read_data add, the the total system will be scaled to a volume that encompasses both boxes. thus any coordinate in either data file with non-zero image flags will become invalid.
the most straightforward way to handle this would be to merge the coordinates with an external tool where you first write out unwrapped coordinates and then you can adjust the box any which way you like and later convert coordinates back to wrapped coordinates with image flags, but using the new box as reference.
the second option would be to use the change_box command where the group id is that of an empty group (so no atoms get displaced) and then change the box dimensions so that they are exactly the same for both parts. if done correctly, you can add the second data file without requiring any scaling of the box from the first data file and the image flag related problem becomes a non-issue.
there is no automated way of doing this, because it is not possible to tell what the intentions are and thus to determine how scaling/merging should be done.
that is why LAMMPS doesn’t even try.
Equilibrating a system under periodic boundaries and then “splitting” it to insert another (similarly equilibrated) system is not such a good idea in the first place.
you have 4 surfaces where there were none, so those atoms at or near the surface are suddenly in “unexpected” environment and far from equilibrium.
Also, especially for Molecular systems, those surfaces can be quite rugged and thus you need to leave extra space so that you have no overlaps or close contacts with unphysically high potential energy. and thus you will need a long time to Re-equilibrated the system.
if, however, you set up each subsystem with fixed boundaries and use a wall fix to keep atoms from leaving the box while already having them equilibrate at an interface, will make the interface smoother and the atoms are already rearranged for an interface, so you will save equilibration time. plus, you won’t need to leave as much extra vacuum, which means less effort to get rid of any vacuum bubbles.
Thank you so much.It is so insightful. It really helps.
I am a bit puzzled on the suggestion:
“if, however, you set up each subsystem with fixed boundaries and use a wall fix to keep atoms from leaving the box while already having them equilibrate at an interface, will make the interface smoother and the atoms are already rearranged for an interface, so you will save equilibration time. plus, you won’t need to leave as much extra vacuum, which means less effort to get rid of any vacuum bubbles.”
To use a wall fix, does it mean using the fix wall command do something like: “fix wallhi all wall/lj93 xlo -1.0 1.0 1.0 2.5 units box”, which imposes a force perpendicular to the wall?
In the suggestions, is it based on an assumption that only one simulation box is used? In this simulation box, the two materials are located in two parts. To equilibrate one part, the other part is assumed to be a wall.
Still use a p p p boundary condition for the wall fix case?
Please give me some instructions on this case, Axel. Thank you very much.
I used the fixed boundary and fix wall to constraint the simulation box. To meet the pppm calculation condition, I also used kspace_modify command. Please see the following scripts. However, when running the case, the program gave an error
"ERROR on proc 1: Out of range atoms - cannot compute PPPM (…/pppm.cpp:1940)
ERROR on proc 3: Out of range atoms - cannot compute PPPM (…/pppm.cpp:1940)
"
Could you please have a look at the simulation setup?
Regards,
Michael
units real boundary p p f
atom_style full
bond_style harmonic
angle_style harmonic
pair_style lj/charmm/coul/long 9.0 10.0 10.0
pair_modify mix arithmetic
kspace_style pppm 0.0001 kspace_modify slab nozforce …
you should be using a coul/long pair style for this. the kspace slab mechanism makes modifications to the box dimensions and re-enables periodic boundaries (kspace must be periodic) and that defeats the purpose of the walls.
In my script (dimension 2), I plan 1) keep all at 300 K using nvt for 10 ps, 2) heat all from 300 K to a variable/atom value as a function of atomic positon and time using fix angevin with nve within 100 ps, 3) keep all at the value for 10 ps, and then 4) cool all from the value to 300 K within 100 ps.
velocity all create 300.0 451547 rot yes dist gaussian
timestep 0.002
fix 1 all nvt temp 300 300 0.2 drag 0.5 # 1)
fix 11 all enforce2d
run 5000
unfix 11
unfix 1
####define varibale Temp
units real
boundary p p f
atom_style full
bond_style harmonic
angle_style harmonic #pair_style lj/long/coul/long 9.0 10.0 10.0 pair_style coul/long 10.0
pair_modify mix arithmetic
…
bond_coeff 1 1000.0 1.0
angle_coeff 1 1000.0 109.47
pair_coeff 1 1 0.1553 3.166
pair_coeff 2 2 0.0 2.058
group spce type 1 2
It reminded that “ERROR: Incorrect args for pair coefficients (…/pair_coul_long.cpp:211) Last command: pair_coeff 1 1 1.0 1.0 2.5”.
Then the pair_coeff changes to the following: pair_coeff 1 1 pair_coeff 2 2
when running the Lammps, it said “ERROR: Pair style requires a KSpace style (…/pair_coul_long.cpp:246) ”
So a pppm or ewald is still needed to complete the calculation?
what you did is not what I suggested. on the contrary.
perhaps this is the point where you need to talk to somebody local.
this isn’t much of a LAMMPS problem, but a problem of your approach to problem solving and understanding implications of the actual physics and requirements of the models you are using. this is beyond the scope of this mailing list.
Thank you, Axel. I got the point. I’ll find someone to guide me on this problem. I am indeed not familiar with the physics and requirements of the models.
I have to take back my previous statement.
I just noticed that my e-mail did not contain the text that I wanted to write.
It was supposed to mean you “should not be using coul/long” but rather something with coul/cut or similar.
Thank you, Axel. That’s really really helpful. I went through the examples distributed with the source code and saw many models using the coul/cut method. I’ll try the coul/cut on my problems.
just so there is no misunderstanding (again), please note that pair style coul/cut only implements coulomb. that is not what you want. you need to use a “coul/cut” variant, i.e. a pair style with both the lennard-jones part and coulomb or pair style hybrid overlay that combines a lennard-jones variant with a coulomb-only pair style. coul/cut is not the only option. there are other coulomb method that do not require a kspace style (wolf, dsf, etc.). details are in the manual.
Thanks so much for instructing me on both theory and technique. Really appreciate it for so nice, insightful and helpful emails especially professors are busy.
using pair_style lj/cut 2.5 is very bad. you are modeling something very different from what you had before, and something that I explicitly warned you about in a previous e-mail.
the cutoff of 2.5 angstrom is too short for the unit setting you have and you are missing out on the coulomb interaction. that means that you are equilibrating to a very different state.
bottom line. the fact that a simulation does not stop with an error does not mean that it is a correct simulation.
Thank you, Axel. I corrected it to pair_style lj/cut/coul/wolf 0.2 9.0 10.0. It works and gives a result as the attached figure (It is a box with 1000 spce water molecules). Really appreciate your help.
