Automatic triclinic boxes in LAMMPS?

[t-brink]

When I run a LAMMPS simulation with a barostat and I also apply a barostat to shear components, the simulation can fail when the input simulation cell is orthogonal:

job2 = pr.create_job(pr.job_type.Lammps, "test2")
job2.structure = <something orthogonal>
job2.potential = ...
job2.calc_md(
    temperature=1000.0,
    temperature_damping_timescale=100.0,
    pressure=[0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
    pressure_damping_timescale=500.0,
    time_step=1.0,
    n_ionic_steps=100000,
)

The reason is that LAMMPS needs to be told if shear of a box is allowed. The command is change_box all triclinic. Shouldn't this command be added automatically when applying a barostat with shear components?

[t-brink]

Some notes so that the complications are not forgotten, I will have to dig deeper later:

• This could be solved by always writing a "triclinic" input file for LAMMPS in pyiron.lammps.structure.LammpsStructure.simulation_cell(). Since "triclinic" in LAMMPS does not really describe the shape of the simulation box (it really defines what features are available) it does not make much sense to force every orthorhombic box to not have the triclinic features. I sometimes want to start from a cubic box and shear it.
• But, currently, only fully periodic boxes are allowed to be triclinic (pyiron.lammps.structure.UnfoldingPrism.__init__()), which makes some sense, since there are limitations on what is allowed in LAMMPS. I would like to relax this a bit, since I have successfully used triclinic boxes with open boundaries before and there are valid use cases. I think the problematic part is mostly related to barostatting/deforming shear components in non-periodic directions. My feeling is currently that the checks should be done wherever the simulation is set up, not in the structure object. Could be wrong, though.

Tread carefully here.

[liamhuber]

@jan-janssen I don't believe this is code smell.

When I run a LAMMPS simulation with a barostat and I also apply a barostat to shear components, the simulation can fail when the input simulation cell is orthogonal

This sounds like a bug

[jan-janssen]

As far as I see this was addressed by yourself in pyiron/pyiron#1083 so we can close this issue. If this is not fixed feel free to open the issue again.

[t-brink]

This has nothing to do with pyiron/pyiron#1083, this bug still exists. If you look at the files in my second comment, they have not been touched by the referenced PR. There is currently a database problem on the cluster, so I cannot test it, but would you please run the provided minimal example in the original comment before assuming this is fixed?

[jan-janssen]

What about?

job2.calc_md(
    temperature=1000.0,
    temperature_damping_timescale=100.0,
    pressure=0.0,
    pressure_damping_timescale=500.0,
    time_step=1.0,
    n_ionic_steps=100000,
)

[t-brink]

I think you are misunderstanding the bug. What you posted might well work, but is completely irrelevant because it doesn't change the shear of the box (if I remember correctly). What about:

job2.structure = <something orthogonal>
job2.potential = ...
job2.calc_md(
    temperature=1000.0,
    temperature_damping_timescale=100.0,
    pressure=[0.0, 0.0, 0.0, 6.0, 7.0, 0.0],
    pressure_damping_timescale=500.0,
    time_step=1.0,
    n_ionic_steps=100000,
)

(Please also note the crucial part that the input structure has to be orthogonal with zero shear components! If the input box is already "triclinic", as it is called in LAMMPS, then this should works.)

[jan-janssen]

So you mean communicating the change box from pyiron to LAMMPS once the strain changes, like it is done for interactive jobs in #63

[samwaseda]

I'll take care of this one.