pyiron / pyiron_atomistics

pyiron_atomistics - an integrated development environment (IDE) for atomistic simulation in computational materials science.
https://pyiron-atomistics.readthedocs.io
BSD 3-Clause "New" or "Revised" License
44 stars 15 forks source link

Quasi Harmonic Approximation does not perform structure optimisation ? #467

Closed samwaseda closed 2 years ago

samwaseda commented 2 years ago

I can see that in the Quasi-Harmonic Approximation job, there is no structure minimisation done for each of the strained structures, which gives wrong displacement fields. I did a small calculation in a Ni-4x4x4-cell containing 1 H atom at 2% strain, which is frankly very small. You can see the error on the free energy in the figure below. We can say that the error per atom is relatively small, but since we're mostly interested in the defect energy, the total energy matters, and I have to say an error of 0.1 eV is devastating.

I opened a PR here, in which I raise a ValueError if the structure is not perfectly symmetric, meaning QHA does calculations only for things like fcc unit cell and so on.

There is actually a way to correct this error even without performing a structure optimization, and that's something I also implemented in my personal QHA class, so if there's someone who really needs to do QHA containing defects, they can talk to me.

Unknown
jan-janssen commented 2 years ago

Does not this depend on the reference job? If you assign a reference job with calc_minimize() it should optimize the atomistic positions of the structure.

samwaseda commented 2 years ago

In QHA we create an artificial displacement field, for which the force field is calculated. If you set calc_minimize, i.e. a structure optimisation is performed, then the displacement field is lost. Of course, for Murnaghan it's a different story.

jan-janssen commented 2 years ago

Can you post an example?

samwaseda commented 2 years ago

Of what?

samwaseda commented 2 years ago

You mean something like this?

lmp = pr.create.job.Lammps('lmp_phono')
lmp.structure = my_structure
lmp.calc_minimize()
phono = lmp.create_job('PhonopyJob', 'phono')
qha = phono.create_job('QuasiHarmonicJob', 'qha')
qha.run()

This won't work because we explicitly (and above all: correctly) forbid calc_minimize for Phonopy

jan-janssen commented 2 years ago

Ok, so you basically want the relaxation for each volume before calculating the phonons? The risk with this is that you could end up with different symmetries for different strains, resulting in different phonon modes.

samwaseda commented 2 years ago

Yes, but that's not what I want. That's what I must do if I want to calculate the free energy for a given volume.

The risk with this is that you could end up with different symmetries for different strains, resulting in different phonon modes.

Yes, I will end up with the correct phonon modes. If there's a symmetry breaking (which is usually beyond QHA though), then what the current implementation returns is so badly wrong that it is not even an approximation.

jan-janssen commented 2 years ago

... (which is usually beyond QHA though) ...

But that is the name of the job, right? When I read Quasi Harmonic Approximation then I assume the symmetry is the same for all strains.

samwaseda commented 2 years ago

Yes, so either the code should fail or the user knows what they are doing. Either way is fine by me but currently our class simply meddles through.

jan-janssen commented 2 years ago

I am fine with the error message you add in https://github.com/pyiron/pyiron_atomistics/pull/466 I am just not sure how to add the structure optimization in the workflow. How do you do it currently?

samwaseda commented 2 years ago

I have my own force fitting class, in which I estimate the optimised structure based on the forces and the energies (which works perfectly if there's no structure transformation). I don't think we can twist our Phonopy class to include this though. I'll try to put it in pyiron as soon as I have the feeling that it works always reliably.