search

mir-group / flare

An open-source Python package for creating fast and accurate interatomic potentials.
https://mir-group.github.io/flare
MIT License
274 stars 64 forks source link

ASE-OTFMD Temperature #274

Closed fgimbert closed 1 year ago

fgimbert commented 3 years ago

Hi FLARE team !

I am trying to learn how to use FLARE with Li3PO4 structure for OTF MD using ASE and VASP calculator. I am encountering a problem with temperature and I would like to know if it is due to me or maybe problems with ASE/VASP.

As mentioned in tutorial, I create random perturbations (by using pymatgen function) in the crystal to break symmetry and set std parameter to 1. At each step, I add only 3 atoms to training set. I tried VelocityVerlet and NVTBerendsen as md_engine, temperature at 500K and 100 MD steps.

For all the 100 steps, DFT calculations is called which is probably ok (while for perfect crystal only for the first steps). Maybe I need to add more steps. The temperature is increasing a lot until more than 3000 K and never stabilizing. Is it fine ? I saw there was a problem of temperature exploding with ASE module last year.

For such materials (Li3PO4 and similar), do you have any advice how to use FLARE ? My objective will be to realize similar process as in this paper (https://pubs.acs.org/doi/10.1021/acs.chemmater.9b04663) and I would really appreciate any help !

Last step output :

-Frame: 99 Simulation Time: 0.099 ps El Position (A) GP Force (ev/A) Std. Dev (ev/A) Velocities (A/ps) Li 2.4748 4.5497 6.0929 0.0000 0.0000 0.0000 0.2841 0.0000 0.0000 0.0000 0.0000 0.0000 Li 5.0266 0.9863 3.3808 0.8744 -0.0087 0.7862 0.2806 0.0000 0.0000 3.1317 30.6178 2.8243 Li 2.5526 1.5093 1.6014 -0.2848 1.6967 0.4876 0.2782 0.0000 0.0000 37.0446 2.4263 -4.0196 Li 5.0088 3.3325 0.9203 0.0756 0.2771 -0.0006 0.2782 0.0000 0.0000 0.1365 4.5572 11.7418 Li 5.0417 3.4766 4.2889 -0.6327 -0.2450 1.7050 0.2834 0.0000 0.0000 -16.7394 -44.4389 10.2649 Li 2.5953 1.5675 4.5902 -0.4455 1.1949 1.3598 0.2844 0.0000 0.0000 -4.0518 -0.2787 12.7069 P 2.4661 4.4512 2.8837 1.9769 -0.9405 4.6768 0.5449 0.0000 0.0000 12.0459 9.8351 1.4449 P -0.0308 0.8379 0.2903 0.0290 1.7307 -5.5441 0.5448 0.0000 0.0000 -10.5808 0.7882 -9.9637 O 0.7294 1.7525 5.2533 -2.6329 -0.9533 1.9703 0.2692 0.0000 0.0000 21.8784 10.7281 -10.5757 O 3.2831 0.6255 3.3114 -4.7531 -3.3571 -3.2352 0.2685 0.0000 0.0000 4.3684 15.5378 1.5692 O 0.8216 4.5236 3.2821 5.6842 -0.9776 -3.1150 0.2672 0.0000 0.0000 -37.6938 -11.6845 12.3707 O 3.3003 0.9799 0.2576 2.5578 -0.5712 -1.8406 0.2690 0.0000 0.0000 -19.1150 -4.9548 15.1440 O 2.9415 3.7106 1.4480 -1.3553 1.4108 4.5181 0.2691 0.0000 0.0000 9.4394 7.5929 18.7382 O 0.4468 1.3999 1.5473 1.9605 1.9916 2.0085 0.2616 0.0000 0.0000 15.0434 -18.7137 7.4201 O 0.5493 4.7952 -0.0398 -0.6763 -1.2987 0.6075 0.2675 0.0000 0.0000 13.3918 -8.4952 -20.9826 O 3.1460 3.6237 4.3178 -2.1961 0.4756 -3.3084 0.2715 0.0000 0.0000 -2.3006 -17.5485 -5.7898

Periodic cell (A): [[4.922448 0. 0. ] [0. 5.296402 0. ] [0. 0. 6.176136]]

Stress tensor (GPa): xx yy zz yz xz xy 68.921 60.970 70.907 0.716 -4.650 0.166

Pressure (GPa): 66.932719 Temperature (K): 3650.78 Kinetic energy (eV): 7.078492 Potential energy (eV): -104.585162 Total energy (eV): -97.506670

Wall time from start: 17274.23 s

Calling DFT...

DFT run complete. Number of DFT calls: 100 Wall time from start: 17289.61 s

*-Frame: 99 Simulation Time: 0.099 ps El Position (A) DFT Force (ev/A) Std. Dev (ev/A) Velocities (A/ps) Li 2.4748 4.5497 6.0929 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Li 5.0266 0.9863 3.3808 0.4090 -0.3023 0.1696 0.0000 0.0000 0.0000 3.1317 30.6178 2.8243 Li 2.5526 1.5093 1.6014 -0.9858 1.3817 1.6472 0.0000 0.0000 0.0000 37.0446 2.4263 -4.0196 Li 5.0088 3.3325 0.9203 -0.1441 -0.1739 0.4286 0.0000 0.0000 0.0000 0.1365 4.5572 11.7418 Li 5.0417 3.4766 4.2889 -1.0258 -1.2475 2.1646 0.0000 0.0000 0.0000 -16.7394 -44.4389 10.2649 Li 2.5953 1.5675 4.5902 -0.4544 1.1597 1.3687 0.0000 0.0000 0.0000 -4.0518 -0.2787 12.7069 P 2.4661 4.4512 2.8837 1.7935 -1.0216 4.9934 0.0000 0.0000 0.0000 12.0459 9.8351 1.4449 P -0.0308 0.8379 0.2903 -0.2817 1.6701 -5.9248 0.0000 0.0000 0.0000 -10.5808 0.7882 -9.9637 O 0.7294 1.7525 5.2533 -2.4877 -0.5457 1.3442 0.0000 0.0000 0.0000 21.8784 10.7281 -10.5757 O 3.2831 0.6255 3.3114 -3.4490 -3.2854 -3.3157 0.0000 0.0000 0.0000 4.3684 15.5378 1.5692 O 0.8216 4.5236 3.2821 5.8599 0.9572 -3.8852 0.0000 0.0000 0.0000 -37.6938 -11.6845 12.3707 O 3.3003 0.9799 0.2576 3.2044 -1.1052 -2.2612 0.0000 0.0000 0.0000 -19.1150 -4.9548 15.1440 O 2.9415 3.7106 1.4480 -1.4688 1.5783 4.4966 0.0000 0.0000 0.0000 9.4394 7.5929 18.7382 O 0.4468 1.3999 1.5473 1.5225 2.0719 3.1801 0.0000 0.0000 0.0000 15.0434 -18.7137 7.4201 O 0.5493 4.7952 -0.0398 -0.3985 -1.8063 0.3554 0.0000 0.0000 0.0000 13.3918 -8.4952 -20.9826 O 3.1460 3.6237 4.3178 -1.9765 0.7240 -3.7429 0.0000 0.0000 0.0000 -2.3006 -17.5485 -5.7898

Periodic cell (A): [[4.922448 0. 0. ] [0. 5.296402 0. ] [0. 0. 6.176136]]

Stress tensor (GPa): xx yy zz yz xz xy 20.204 6.584 13.980 1.602 -1.378 0.422

Pressure (GPa): 13.589454 Temperature (K): 3650.78 Kinetic energy (eV): 7.078492 Potential energy (eV): -104.522050 Total energy (eV): -97.443558

Wall time from start: 17289.74 s mean absolute error: 0.3946 eV/A mean absolute dft component: 1.7452 eV/A Adding atom [5, 7, 6] to the training set. Uncertainty: [0.28436108 0. 0. ]

Run complete.

YuuuXie commented 3 years ago

Hi @fgimbert, thanks for reporting this issue.

Just a few questions (since recently there're some issues #269 #273 related to the ASE OTF training)

Also, if you want to call DFT fewer times, you can try

fgimbert commented 3 years ago

Hi ! I was able to do more tests and I found that with use_mapping = False, the temperature was not exploding anymore. For now I will use so the function to map the GP force field for LAMMPS after the run is finished. It should give the same result as use_mapping = True ?

The GP forces were looking fine (different for different frames).

YuuuXie commented 3 years ago

Hi ! I was able to do more tests and I found that with use_mapping = False, the temperature was not exploding anymore. For now I will use so the function to map the GP force field for LAMMPS after the run is finished. It should give the same result as use_mapping = True ?

The GP forces were looking fine (different for different frames).

The mapping is an approximation method of GP, and it is possible that the 2body grid number and 3body grid number of your settings are too small such that the mapping has a large difference compared to the original GP prediction. And that's probably why when you use mapping the system behaves strangely. If the grid numbers are set to large enough the mapping results should be almost the same as GP but that might be expensive then. For now you can set use_mapping=False for the otf training and then map them.

xig126 commented 3 years ago

Hi ! I was able to do more tests and I found that with use_mapping = False, the temperature was not exploding anymore. For now I will use so the function to map the GP force field for LAMMPS after the run is finished. It should give the same result as use_mapping = True ? The GP forces were looking fine (different for different frames).

The mapping is an approximation method of GP, and it is possible that the 2body grid number and 3body grid number of your settings are too small such that the mapping has a large difference compared to the original GP prediction. And that's probably why when you use mapping the system behaves strangely. If the grid numbers are set to large enough the mapping results should be almost the same as GP but that might be expensive then. For now you can set use_mapping=False for the otf training and then map them.

Hi @YuuuXie ,

Would you mind sharing your experience on choosing suitable grid numbers for MGP? I'm currently running otf md for a 32-atom unit cell. I tested the grid sizes of MGP from 16 to 128 for "two body" function and 16 to 64 for "three body", but I found the temperature of the system still exploded even with the largest size of the grid.

Thanks.

YuuuXie commented 3 years ago

Hi @xig126, for 2-body mapping, you can set a large number of grid (128 or larger), it is ok because 2-body mapping is always very cheap to build. For 3-body mapping, usually [32,32,32] to [128,128,128] are reasonable, really depending on your cutoff. I would assume that [64, 64, 64] is enough.

As for your temperature exploding issue, maybe try turning off mgp as I mentioned here https://github.com/mir-group/flare/issues/273#issuecomment-815922829

YuuuXie commented 1 year ago

Close stale issue