How to allow angles after simulation_box is defined?

Dear all,

Months ago, I ran the simulation with this simulation box:

create_box  3 simulation_box bond/types ${n_bond_type} extra/bond/per/atom ${n_bond_per_atom}

Because the equilibration is too time-consuming (~10 million timesteps), I always read the restart file which has reach equilibrim to save time. Today, I find it necessary to incorporate angles besides bonds. But I don’t want to create the simulation box afresh because that will be time-consuming. I wonder how to allow angles after read_restart.

NOTE: I am using the RHEO and BPM packages.

I came up with two ideas.

The first is to create a new simulation box which allows angles, output it to a dump file. I only want to change the box but not the atom attributes, so I use the following command to replace the box:

read_dump all-10b.dump 10000000 vx box yes replace no add no trim no purge no 

However, this unexpectedly changed some attributes, but I cannot identify what it changed. Perhaps this is related to the fact that bond_style bpm/spring requires bond history while read_dump modifies the bond history.

The second idea is to use the change_box command, but it does not allow users to incorporate angles … Worse still, users cannot change_box after reading restart file with per-atom info. But I don’t think allowing angles will cause atoms to move to new procs.

Any advice or comment will be appreciated!

Hi @jasmine969,

If I understand correctly what you try to achieve, note that you can convert a restart file to a data file easily, for example with the restart2data flag or a very simple input script. See here if you don’t know about this command line option.

This would allow you to read your file as a data file and add the required angles to your simulation. See the doc explaining the extra/angle/types, extra/angle/per/atom and extra/special/per/atom options. Note that you would have to assign them to the related atom nonetheless.

bond_style bpm/spring command — LAMMPS documentation

The reference state is stored by each bond when it is first computed in the setup of a run. Data is then preserved across run commands and is written to binary restart files such that restarting the system will not reset the reference state of a bond.

Hi @Germain , bond bpm/spring needs the initial state of the bonds, which is stored in the restart file but not the data file. If I use restart2data, I’m afraid the initial bond states will be lost.

Then I think you should do the opposite of what you describe. That is

  1. read your initial box adding the angles
  2. convert your restart file to a dump containing the coordinates of your atoms and box and
  3. read the coordinates of your atoms and the box from your dump file.

If you have any other info you need from the binary file, this will not do.

The change_box command is not suited for that kind of modification.

The last idea I can come up with is to list all the triplets of atoms you want to have an angle with, and script a series of create/bonds single/angle commands. You would still need to check that the angles per atom and special list are large enough to add your new angles.

You are starting a new simulation and the added angle potential is a significant change to the potential, so you will need to re-equilibrate your system and preserving any details from the restart files outside of positions, velocity, and topology is not of any importance.

You cannot change type information after the box is created.
The only way to add (harmonic) angle interactions while keeping the rest of the system is to use the fix restrain command — LAMMPS documentation
This has its own set of limitations and inconveniences, of course.

But I really need the initial bond states from the binary file.

Yes, indeed. However, it will take ~10 million timesteps if running from the initial state, but will take only 20 thousand timesteps if running from the binary file after adding angle potentials.

This is a matter of expediency. But I need to create over 10 thousand angles, then should I write 10 thousand fix restrain (though it can be done automatically by for-loop)?

There is no solution that addresses all your requirements.

Without knowing more details it is impossible to judge whether those requirements are meaningful. I suspect they are not.

1 Like

Alright. I will run from the beginning. And I still hope that in future, users can modify whether bond/angle/dihedral/improper is allowed after the simulation box is defined.

This is a fundamental design decision in LAMMPS and changing it would be a massive undertaking to be changed and likely cause lots of issues with the (large amount of) existing code. It is thus very unlikely to happen for such corner-case applications. Again, first you have to be successful in convincing people that this effort is warranted and your approach has not at all helped.

You can easily avoid this issue by reserving sufficient space on spec for such interactions ahead of time and then not use them until you need them.

1 Like

Unless you have created/broken bonds during the original run, here’s a solution that I think should work:

  1. Load your equilibrated restart file w/o angles and write a dump file of atom coordinates/velocities
  2. Reinitialize the system with all angles + bonds. The bonds will save their reference state.
  3. In your new system, use read_dump with replace yes to load the dump file and overwrite atom coordinates to restore your equilibrated geometry

Maybe this is what you need? Depending on the precision of the dump file, you may (or may not) want to do a short relaxation afterwards.

1 Like

Wow, it works. Thanks a lot!

1 Like