Sudden jumps in MLFF energies

[xiaolin_liu1]

I'm trying to run a long NVT simulation for an amorphous material at 300 K. I first train MLFF in the NPT ensemble in chunks of 20 ps to avoid large deformations of the cell for 60 ps. This is the last INCAR file for training:

Code: Select all

SYSTEM = alumina 
ISTART =0                                                                                                                   
ICHARG =2                                                                                                                  
ISYM = 0
ML_LMLFF         = .TRUE.
ML_MODE = train   

POTIM = 1                                                                                                             
GGA = PE                #   RP stands for RPBE, PE for PBE, 91 for PW91                                                       
 -----------------------------------------                                                                                  

                                                                             
Parameter for parallelization:                                                                                              
---------                                                                                                                   
--------------------------------                                                                                            
LPLANE = .TRUE.                                                                                                          
------------------------------------------                                                                                  
Electronic Relaxation:                                                                                                      
------------------------------------------                                                                                  
   LREAL  = Auto                                                        
   ENCUT  = 520                                                                                                      
                                                                                                                            
   PREC   = Normal          
   EDIFF  = 1e-5                                           
   ALGO  = Normal                          
   NELM=  60                                                                                    
   NELMIN = 6            
   ISPIN  = 1                                                      
   IVDW = 11                                                                                                               
-------------------------------------------                                                                                 
Ionic Relaxation:                                                                                                           
-------------------------------------------                                                                                                                    
   NSW    = 20000                                                                  
   #AIMD using Nose-Hoover                                                 
   IBRION = 0                 
   ISIF   = 3                  
   TEBEG = 300
   TEEND = 300 
   MDALGO = 3

   LANGEVIN_GAMMA=5  5  5 5
   LANGEVIN_GAMMA_L=5 
--------------------------------------------                                                                                
DOS Related:                                                                                                                
--------------------------------------------                                                                                
   ISMEAR = 0                                                                                                               
   SIGMA  = 0.05                                                                                                                                                                                                                
--------------------------------------------                                                                                
Output Control:                                                                                                             
--------------------------------------------                                                                                
   LWAVE  = .FALSE.                                                                                                         
   LCHARG= .FALSE.
   NWRITE = 1
   ML_MB = 3000   
   NCORE= 16

Then I did a refit to get the following errors:

Code: Select all

# ERR               nstep      rmse_energy       rmse_force      rmse_stress
# ERR                   2                3                4                5
# ERR ######################################################################
ERR                     0   1.50428578E-04   3.34169445E-02   3.21668029E-01

I averaged the errors from 50 snapshots that are equidistant in time from a 50-ps run using the MLFF to get test errors.
They are :

Code: Select all

 energy error per atom, rmse force,     rmse stress
3.325842696609328e-05, 0.03471839929114923, 0.35828575204040314

The errors seem fine. I then run an NVT simulation using the MLFF only. This is the INCAR:

Code: Select all

SYSTEM = alumina 
ISTART =0                                                                                                                   
ICHARG =2                                                                                                                  
ISYM = 0                                                                                                     
ML_LMLFF         = .TRUE.
ML_MODE = run
ML_OUTBLOCK = 20
GGA = PE                #   RP stands for RPBE, PE for PBE, 91 for PW91                                                    
 
 -----------------------------------------                                                                                  

                                                                             
Parameter for parallelization:                                                                                              
---------                                                                                                                   
--------------------------------                                                                                            
LPLANE = .TRUE.                                                                                                          
------------------------------------------                                                                                  
Electronic Relaxation:                                                                                                      
------------------------------------------                                                                                  
   LREAL  = Auto                                                              
   ENCUT  = 520                                                                                                      
                                                                                                                            
   PREC   = Normal        # precision:Medium/High/Low   default medium . the same time, ROPT accuracy of approximately     
                             # 1meV/atom (10e-3). recommend to use PREC Normal for everyday calculations in VASP.4.5        
   EDIFF  = 1e-5           # converge criterion: default=1e-4                                                             
   ALGO  = Normal             # algotithm (8-CG, 48-RMM), Fast means intially use DAV algo then use RMM.                     
   NELM=  60                  # default is 60                                                                               
   NELMIN = 6                # the minimum number of electronic SC steps  Default is 2                                      
                                                
   ISPIN  = 1                # spin polarized calculation (2-yes, 1-no)                                                     
   IVDW = 11                                                                                                               
-------------------------------------------                                                                                 
Ionic Relaxation:                                                                                                           
-------------------------------------------                                                                                                                    
   NSW    = 4000000           #     number of steps for IOM                                                                   
   #AIMD using Nose-Hoover                                                 
   IBRION = 0                   # choose molecular dynamics 
   MDALGO = 2                   # using Nose-Hoover thermostat    
   ISIF = 2                     # compute stress tensor but do not change box volume/shape 
   TEBEG = 300
   TEEND = 300 
   POTIM = 1.0                  # time step in femto seconds 
   SMASS = 1.0                  # setting the virtual mass for the Nose-Hoover thermostat                                                                                    
--------------------------------------------                                                                                
DOS Related:                                                                                                                
--------------------------------------------                                                                                
   ISMEAR = 0                                                                                                               
   SIGMA  = 0.05                                                                                                             
   LORBIT =11                                                                                                                
--------------------------------------------                                                                                
Output Control:                                                                                                             
--------------------------------------------                                                                                
   LWAVE  = .FALSE.                                                                                                         
   LCHARG= .FALSE.  
   NCORE= 16

A Nose-Hoover thermostat was used, and this is the energy profile:

test_energies.png

The green curve is the running average over 2ps. I see sudden jumps in potential energies. Is it unphysical, or should it be seen as a relaxation of the material?

To get better statistics, I also simulated with a 2 by 2 by 2 supercell containing 1424 atoms. The change in potential energy is more severe, and I'm also getting spikes in the energy profile:

2units_300_npt_ML_run_222_energies.png

Does this mean the MLFF is unable to handle the larger supercell?
Thanks for your insights!

[jonathan_lahnsteiner2]

Dear Xiaolin Liu,

What I saw from the files you provided is that you are training at 300K. Then you are running your simulation in production mode at the same temperature. It is always advisable to incorporate trainings data from higher temperatures into your training data. So I would recommend to include trainings data from maybe 400K into your trainings data by running training continuation in the NPT ensemble. Otherwise it is not guaranteed to get a stable machine learning force field for the desired temperature. The additional training will show more of the available phase space to the MLFF and therefore will make it more stable. What atomic species do you have in your system. If you have light atoms as hydrogen your time step might be too large. Please let me know which atoms you have.

In principle the MLFF should be able to handle larger boxes without issues. Please try to include trainings data from higher temperatures in your data set and check if the issues remain.

All the best Jonathan

[xiaolin_liu1]

Hi Jonathan,

Thanks for your suggestions. The lightest atom in my system is Oxygen, so a time step of 1 fs should be fine.

What I'm trying to do is run long NVT trajectories at 300 K to 500 K in an increment of 50 K. What's the best strategy to train VASP MLFF for this purpose? Should I train the MLFF at each temperature? Or should I use a single force field trained at 500 K and validate at each temperature?

[jonathan_lahnsteiner2]

Dear Xiaolin Liu,

I would suggest to do a heating run during training. To start with a temperature slightly below the value you want to run at and then heat up until the 500K which I suggested. To apply an NPT heating run you have to set the following tags in your INCAR file:

Code: Select all

 IBRION = 0                 

 ISIF   = 3                  

 TEBEG = 200
 TEEND = 500 
 MDALGO = 3

 LANGEVIN_GAMMA=5  5  5 5
 LANGEVIN_GAMMA_L=5 

After this I would probably do a short NPT run at 500K, just to be sure enough structures were picked up. I would recommend to use a temperature gradient of roughly 2K/ps faster.
The higher temperature structures are important because atoms do larger displacements at higher temperatures on average. So it is more likely to sample those and add to your data set. But also the step wise procedure you mention will be possible, but also there I would recommend constant temperature NPT simulations for each temperature step.

All the best Jonathan