How to restart a LMPOTF run?

[bastonero]

Hello Developers and everyone!

I recently started using the LMPOTF class along with the direct LAMMPS interface input. First, thanks a lot for this, it is truly excellent and useful! Second, from time to time the DFT calculator may fail and I would like to restart the run from that point.

What is the best practice to restart such simulation?

I tried the following:

1. Dump the model during the otf simulation by adding self.sparse_gp.to_json("model.json",self.sparse_gp) in the post_dft_callback method
2. When willing to restart, I load the model via: SparseGP.from_json("./model.json")
3. Add also at the end of the file lmpotf.dft_calls = 1 and lmpotf.save("aiida.flare")

Unfortunately, after the first DFT call the simulation stops with the following error:

terminate called after throwing an instance of 'std::bad_array_new_length'
terminate called recursively

After modifying the source code of the LMPOTF class, it seems the problem is located in this section of the code:

                    self.sparse_gp.add_training_structure(structure)
                    self.sparse_gp.add_specific_environments(
                        structure, atoms_to_be_added
                    )
                    self.sparse_gp.update_matrices_QR()

I honestly don't understand how to solve this, or whether my solution is correct.

Any pointer would be extremely appreciated.

Thanks a lot, and best wishes.

[YuuuXie]

@anjohan

[anjohan]

Hi Lorenzo,

Thank you for trying out the LMPOTF functionality! Sorry for the late reply.

In theory, the JSON approach should work and be the fastest, but I have no experience with the JSON I/O.

A simpler, more reliable approach is to just retrain on the existing training data (assuming this does not take too long). In your Python script, you can pretrain the GP with code along the lines of

from flare.learners.lmpotf import LMPOTF, transform_stress

...

type2number = [28,22]
single_atom_energies = np.array([-10.8,-10.1])

...

lmpotf = LMPOTF(sparse_gp, ...)

number2type = {t:i for i,t in enumerate(type2number)}

frames = ase.io.read(fname, index=":")
for i, frame in enumerate(frames):
    print(f"pretraining on frame {i}")
    types = np.array([number2type[s] for s in frame.get_atomic_numbers()])
    structure = Structure(
                    np.array(frame.get_cell()),
                    types,
                    frame.get_positions(),
                    rcut,
                    lmpotf.descriptors
                )
    E = frame.get_potential_energy() - np.sum(single_atom_energies[types])
    F = frame.get_forces().ravel()
    S = frame.get_stress(voigt=False)
    structure.energy = [E]
    structure.forces = F
    structure.stresses = transform_stress(S)
    lmpotf.sparse_gp.add_training_structure(structure)
    if i==0:
        lmpotf.sparse_gp.add_random_environments(structure, [16])
    else:
        lmpotf.sparse_gp.add_uncertain_environments(structure, [16]) # or something smarter
    lmpotf.sparse_gp.update_matrices_QR()
lmpotf.save(lmpotf.model_fname)
lmpotf.dft_calls = len(frames)

where fname points to the concatenated xyz file of DFT outputs from previous runs. Note that you probably want to remove the initial python otf invoke from the LAMMPS script.

Alternatively, if you are having regular DFT crashes, you can try to modify/override the run_dft method to be more robust and rerun DFT. Note that the run_dft method receives the Structure object and can modify it freely in case you want to rerun with microscopically perturbed positions or anything like that.

[bastonero]

Thanks so much for the kind help!