I am trying to make biomolecular condensates. I have a cubic box (500x500x500 A3) with supposed 100 chains of polymer-type chains. Now I want to make a condensate, so first I use the NPT ensemble (fix npt all npt temp 100.0 100.0 100 iso 1.0 1.0 100), and when density converges with box dimension (let’s say 150x150x150 A3), I use NVT ensemble to equilibrate the system. Everything is going well, but now I want to increase the box dimension to 500 A, which gives me an error because of the PBC effect.
The second way, I tried to run with unwrap coordinate, but it still is giving me errors like
ERROR on proc 5: Bond atom missing in image check (…/domain.cpp:792)
ERROR on proc 6: Bond atom missing in image check (…/domain.cpp:792)
One thing you can try to do is, after unwrapping the coordinates like you did, set up an energy minimization. It could be that when the coordinates were unwrapped across boundaries, some atoms ended up being way too close than they should depending on the conformation of the chains. And then, as they will experience an ultra large force as a consequence of the high proximity, they might end up too far away from one another, causing the error to occur. If this is the case, after the miniization, you should be able to run the dynamics with no problem.
Another thing, when you unwrapped the coordinates of the polymer, check if you zero’ed back the image flags columns. If you didnt, your atoms will still end up too far from one another when LAMMPS reads the data file.
thermo_style custom step temp pe ke etotal lx ly lz pxx pyy pzz
thermo 100
#--------------- Dump configurations at regular intervals -----------------------------
timestep 10
dump 1 all custom 100 input1.lammpstrj id mol type q xu yu zu ## unwrap coordinate
dump_modify 1 sort id
fix npt all npt temp 100.0 100.0 100 iso 1.0 1.0 1000
run 5000 ## run for 1ns
write_data final-structure_npt_100.dat nocoeff
undump 1
unfix npt
dump 2 all custom 100 input1.lammpstrj id mol type q xu yu zu ## unwrap coordinate
dump_modify 2 sort id
fix fxnve all nve
fix fxlange all langevin 100 300 10.0 548669595 ## run for 50ns
run 1000
write_data final-structure_nve_100.dat nocoeff
Just a small remark (not the solution to your problem, as I maintain the solution I mentioned previously, but rather a simple remark): if you want to simulate the NVT ensemble as you said in the initial message, it’s better to use fix nvt instead of fix nve + fix langevin. The way the Langevin thermostat affects the equations of motion is apparently not such that allows you to truly sample the NVT ensemble (see Difference between nve+langevin and nvt dynamics! - #2 by akohlmey).
thermo_style custom step temp pe ke etotal lx ly lz pxx pyy pzz
thermo 100
#--------------- Dump configurations at regular intervals -----------------------------
timestep 10
dump 1 all custom 100 input1.lammpstrj id mol type q xu yu zu ## unwrap coordinate
dump_modify 1 sort id
fix npt all npt temp 100.0 100.0 100 iso 1.0 1.0 1000
run 5000 ## run for 1ns
write_data final-structure_npt_100.dat nocoeff
undump 1
unfix npt
change_box all x final 0 500 boundary p p p
change_box all y final 0 500 boundary p p p
change_box all z final 0 500 boundary p p p
dump 2 all custom 100 input1.lammpstrj id mol type q xu yu zu ## unwrap coordinate
dump_modify 2 sort id
fix fxnve all nve
fix fxlange all langevin 100 300 10.0 548669595 ## run for 50ns
run 1000
write_data final-structure_nve_100.dat nocoeff
initial, I just deleted that line and sent it, but now this is the final one.
You should try visualising your system and the results should be immediately clear. You can’t just change_box when there are bonds stretched across box boundaries, because LAMMPS stores atom coordinates as “internal to” the box (keeping track of wrapping effects via image flags), so when the box gets expanded, atoms on opposite sides of a box boundary simply get that distance added to their existing interatomic distance (which therefore makes their bond lengths completely unphysical).
The easiest thing to do is to dump unwrapped coordinates, start a new LAMMPS simulation run with the correct box size, and then read those unwrapped coordinates back in.
Ah, I thought you had already tried creating a LAMMPS .dat file with the unwrapped coordinates and with your desired box dimension when you said:
Everything is going well, but now I want to increase the box dimension to 500 A, which gives me an error because of the PBC effect.
The second way, I tried to run with unwrap coordinate, but it still is giving me errors like
In any case, my argument of high proximity of atoms after unwrapping the coordinates makes no sense now that I stopped to think more carefully about it… so if the problem persists even after unwrapping the coordinates, I would be out of ideas. The idea of image flags would also not make sense since your system is still periodic…
In a scenario where the problem persists even after unwrapping the coordinates, I suppose that it is in the realm of possibilities that the conformation of one or more chains sitting at the edge are not stable without the other surrounding chains (that will no longer exist) once you enlarge the box. And then maybe it is possible that you have an instability issue leading to an error message if eventually the dynamics leads to atoms flying too far from one another. Maybe an energy minimization in a scenario like this could still help. But then the underlying issue is not really because of “atoms in the initial configuration being too close”.