Unphysical result in MD using trained ACE potential

[WJiangH]

Dear Authors,

I intend to train an ace potential for a high-entropy alloy and collected different types of data from DFT, including bulk and defective system under NVT and NPT simulations, as well as NEB barrier in order to capture the transition states. The dataset is large and therefore I used D-optimiality (by setting up # of function = 1500 and 2000) to select representative data used for training and the rest of data for testing. The log file show a reasonable accuracy with MAE~7mev (I think increase iter steps will improve accuracy further). So I made some test using ace potential in LAMMPS. By varying the lattice constant and calculate the instant energy (see the plot
), the energy with small lattice constant did not go to positive.

Also, when I perform static MD simulation for calculating the formation energy, it seems unstable. for example:

`LAMMPS (27 Jun 2024)
WARNING: Using I/O redirection is unreliable with parallel runs. Better to use the -in switch to read input files. (../lammps.cpp:551)
Reading data file ...
triclinic box = (0 0 0) to (10.616668 10.651563 10.644024) with tilt (-0.0018335981 0.01816288 -0.0074016512)
2 by 4 by 2 MPI processor grid
reading atoms ...
108 atoms
read_data CPU = 0.071 seconds
21 atoms in group Co
22 atoms in group Cr
21 atoms in group Fe
22 atoms in group Mn
22 atoms in group Ni
ACE version: 2023.11.25
Product evaluator is used
Loading /home/wenjiang0716/ace_training/model_rep_400/model_rep_700_11/output_potential.yaml
Total number of basis functions
Co: 37 (r=1) 947 (r>1)
Cr: 37 (r=1) 948 (r>1)
Fe: 37 (r=1) 948 (r>1)
Mn: 37 (r=1) 948 (r>1)
Ni: 37 (r=1) 948 (r>1)
Mapping LAMMPS atom type #1(Co) -> ACE species type #0
Mapping LAMMPS atom type #2(Cr) -> ACE species type #1
Mapping LAMMPS atom type #3(Fe) -> ACE species type #2
Mapping LAMMPS atom type #4(Mn) -> ACE species type #3
Mapping LAMMPS atom type #5(Ni) -> ACE species type #4
ID_X_Y_Z: 60 3 7.12198792685739 5.35045671717239 5.32475633463281
1 atoms in group delete_id
Deleted 1 atoms, new total = 107

CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE

Your simulation uses code contributions which should be cited:

• Type Label Framework: https://doi.org/10.1021/acs.jpcb.3c08419
  The log file lists these citations in BibTeX format.

CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE

Neighbor list info ...
update: every = 1 steps, delay = 0 steps, check = yes
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 8
ghost atom cutoff = 8
binsize = 4, bins = 3 3 3
1 neighbor lists, perpetual/occasional/extra = 1 0 0
(1) pair pace, perpetual
attributes: full, newton on
pair build: full/bin/atomonly
stencil: full/bin/3d
bin: standard
Setting up cg style minimization ...
Unit style : metal
Current step : 0
Per MPI rank memory allocation (min/avg/max) = 4.116 | 4.117 | 4.118 Mbytes
Step Temp PotEng Press Pxx Pyy Pzz Pxy Pxz Pyz Lx Ly Lz Volume
0 0 -829.98668 -29937.384 -30624.581 -29623.33 -29564.241 -60.670947 -1560.9965 617.31179 10.616668 10.651563 10.644024 1203.6699
100 0 -1058.5868 -285371.35 -576126.13 -104777.52 -175210.39 119501.27 132287.6 309227.97 10.616668 10.651563 10.644024 1203.6699
200 0 -1740.9004 -135014.24 139499.78 1508.7726 -546051.29 -443657.84 436187.88 -737319.72 10.616668 10.651563 10.644024 1203.6699
300 0 -5845.5873 -5392059 -7863757.4 -2708723.5 -5603695.9 4957366.2 717140.7 3789853.5 10.616668 10.651563 10.644024 1203.6699
400 0 -39906.5 -39874290 -33202197 -36982422 -49438252 18452448 -6770210 -7947180.1 10.616668 10.651563 10.644024 1203.6699
500 0 -140583.1 -56486193 -6210923.7 -27732394 -1.3551526e+08 45348861 -58859942 18272217 10.616668 10.651563 10.644024 1203.6699
600 0 -420448.31 16859876 1.0042968e+08 1526931.3 -51376983 69001160 -24659668 30066903 10.616668 10.651563 10.644024 1203.6699
700 0 -551839.73 1.6523462e+08 1.6434509e+08 1.7167652e+08 1.5968226e+08 23916357 -21419026 9369950.3 10.616668 10.651563 10.644024 1203.6699
800 0 -602356.21 2.3959591e+08 2.5683149e+08 2.3053058e+08 2.3142567e+08 -10160399 -16568885 -23248133 10.616668 10.651563 10.644024 1203.6699
900 0 -935898.24 1.2749972e+08 2.0040675e+08 1.8245796e+08 -365563.87 53649536 -67188333 -30131160 10.616668 10.651563 10.644024 1203.6699
1000 0 -1675289.9 -9.8137099e+08 -1.2578969e+09 -2.1879812e+08 -1.467418e+09 4.5245684e+08 2.2896784e+08 -65724170 10.616668 10.651563 10.644024 1203.6699
1089 0 -16907212 -9.6741758e+09 -1.5588405e+10 -1.8225208e+09 -1.1611601e+10 -2.0705924e+09 -7.2735849e+08 6.0771053e+08 10.616668 10.651563 10.644024 1203.6699
Loop time of 8.95029 on 16 procs for 1089 steps with 107 atoms

100.0% CPU use with 16 MPI tasks x no OpenMP threads`

The energy drop dramatically. When performing relaxation in LAMMPS, should the energy and force convergence tolerance keep consistency with the setting in DFT? which is relatively coarse.

I also attached my input and log files. I am also wondering if I correctly set the core-repulsion? since the system is a regular fcc and the lattice constant from DFT is around 3.55, there is rare case that atoms will be very close to each others.

input.txt
log.txt

I notice it is very hard to train a potential for high-entropy alloy considering its various local configurations especially for defects. Any suggestions help. Thanks
Best,
Wenjiang

[yury-lysogorskiy]

You are lacking core-repulsion - 99% typical problem of MLIP.
Thre options to fix it:
1a) add the DFT data that have short bond distances to the training set and upfit
1b) generate augmented ZBL data with pace_augment and upfit potential (use aug_factor param)
2) add core repulsion post-factum by using pace_corerep.

Also, I would suggest to uncomment repulsion:auto line in the pacemaker input file and continue fit, since you loss function is still going down quite well (you can simply put max iter to 3000 or 5000 and wait until early stopping will work)