PyXtal_FF stress predict the negative value of DFT stress

[yanxon]

PyXtal_FF stress is trained based on the negative value of the DFT stress: https://github.com/qzhu2017/PyXtal_FF/blob/669c8a46b87d8b89686af8af2a14ba8b77d8bd6b/pyxtal_ff/utilities/__init__.py#L378

We need to fix the PyXtal_FF stress to predict the true DFT stress.

To fix this, we need to replace stress in two places:

1. descriptors, i.e. change this line to positive https://github.com/qzhu2017/PyXtal_FF/blob/669c8a46b87d8b89686af8af2a14ba8b77d8bd6b/pyxtal_ff/descriptors/EAD.py#L171
2. utilities

[yanxon]

@qzhu2017

This is probably why ASE calculator will predict the negative of the stress as well.

If you look at it here, PyXtal_FF stress has the opposite sign with LAMMPS stress: https://github.com/qzhu2017/PyXtal_FF_LAMMPS/tree/main/Cu_snap/Cu_NN

Previously, we impose PyXtal_FF to predict the negative of DFT stress.

Now, I make the change here: 5908e9eecae69cae2cb0615c4fbebf08aed7f46a.

[yanxon]

@qzhu2017

I believe I completed this task. I check with LAMMPS and the energy, force, and stress computations are consistent, including the hybrid ZBL potential & NN potential.

I'm going to close this issue.

[yanxon]

Train with Si dataset.

from pyxtal_ff import PyXtal_FF
import os

TrainData = "data/training.json"
TestData  = "data/test.json"

if not os.path.exists(TrainData):
    if not os.path.exists('data'):
        os.mkdir('data')
    os.chdir('data')
    url = 'https://raw.githubusercontent.com/materialsvirtuallab/mlearn/master/data/Si/'
    print('Downloading the training and test data')
    os.system('wget ' + url + TrainData.split('/')[-1])
    os.system('wget ' + url + TestData.split('/')[-1])
    os.chdir('..')

descriptor = {'type': 'SNAP',
              'weights': {'Si': 1.0},
              'Rc': 5.0,
              'parameters': {'lmax': 3},
              'base_potential': {'inner': 2.0,
                                 'outer': 4.0}
             }

model = {'system' : ['Si'],
         'hiddenlayers': [16, 16],
         'path': 'Si-snap/',
         'optimizer': {'method': 'lbfgs'},
         'force_coefficient': 2e-2,
         'stress_coefficient': 2e-3,
         'alpha': 1e-6,
         'epoch': 1000,
         }

ff = PyXtal_FF(descriptors=descriptor, model=model)
ff.run(mode='train', TrainData=TrainData, TestData=TestData)

Scripts

Si.cif:

# generated using pymatgen
data_Si
_symmetry_space_group_name_H-M   'P 1'
_cell_length_a   5.46872800
_cell_length_b   5.46872800
_cell_length_c   5.46872800
_cell_angle_alpha   90.00000000
_cell_angle_beta   90.00000000
_cell_angle_gamma   90.00000000
_symmetry_Int_Tables_number   1
_chemical_formula_structural   Si
_chemical_formula_sum   Si8
_cell_volume   163.55317139
_cell_formula_units_Z   8
loop_
 _symmetry_equiv_pos_site_id
 _symmetry_equiv_pos_as_xyz
  1  'x, y, z'
loop_
 _atom_site_type_symbol
 _atom_site_label
 _atom_site_symmetry_multiplicity
 _atom_site_fract_x
 _atom_site_fract_y
 _atom_site_fract_z
 _atom_site_occupancy
  Si  Si0  1  0.25000000  0.75000000  0.25000000  1.0
  Si  Si1  1  0.00000000  0.00000000  0.50000000  1.0
  Si  Si2  1  0.25000000  0.25000000  0.75000000  1.0
  Si  Si3  1  0.00000000  0.50000000  0.00000000  1.0
  Si  Si4  1  0.75000000  0.75000000  0.75000000  1.0
  Si  Si5  1  0.50000000  0.00000000  0.00000000  1.0
  Si  Si6  1  0.75000000  0.25000000  0.25000000  1.0
  Si  Si7  1  0.45000000  0.45000000  0.45000000  1.0

Si.lammps.struc:

Si.data (written by ASE) 

8    atoms 
1  atom types
0.0      5.4687279999999996  xlo xhi
0.0      5.4687279999999996  ylo yhi
0.0      5.4687279999999996  zlo zhi

Atoms 

     1   1      1.3671819999999999      4.1015459999999999      1.3671819999999999
     2   1                       0                       0      2.7343639999999998
     3   1      1.3671819999999999      1.3671819999999999      4.1015459999999999
     4   1                       0      2.7343639999999998                       0
     5   1      4.1015459999999999      4.1015459999999999      4.1015459999999999
     6   1      2.7343639999999998                       0                       0
     7   1      4.1015459999999999      1.3671819999999999      1.3671819999999999
     8   1      2.4609275999999998      2.4609275999999998      2.4609275999999998

lammps input file:

# ---------- Initialize Simulation --------------------- 
clear 
units metal 
dimension 3 
boundary p p p
atom_style atomic
read_data Si.data

mass 1 1. # place holder

# Complete pair style and coeff
pair_style hybrid/overlay &
mliap model nn Si-snap/NN_weights.txt &
descriptor sna Si-snap/DescriptorParam.txt &
zbl 2.0 4.0
pair_coeff 1 1 zbl 14.0 14.0
pair_coeff * * mliap Si

variable        e equal pe
variable N equal count(all)

thermo_style custom pe pxx pyy pzz pxy pxz pyz
thermo 1

dump 1 all custom 1 dump.myforce.* id type fx fy fz

run             0
variable        emin equal $e
variable EperAtom equal (${emin}/${N})

Prediction Output

PyXtal_FF via ASE calculator:

Energy (eV/atom):  -5.313921748054366

Forces (eV/A): 
[[-0.597469346  0.599759447 -0.597469346]
 [-0.224402594 -0.224402594  0.219674215]
 [-0.597469346 -0.597469346  0.599759447]
 [-0.224402594  0.219674215 -0.224402594]
 [-0.803306282 -0.803306282 -0.803306282]
 [ 0.219674215 -0.224402594 -0.224402594]
 [ 0.599759447 -0.597469346 -0.597469346]
 [ 1.6276165    1.6276165    1.6276165  ]]

Stress (eV/A^3): 
[ 0.007724698  0.007724698  0.007724698 -0.007505016 -0.007505016
 -0.007505016]

LAMMPS output:

Energy (eV/atom):  -5.31392175

Forces (eV/A): 
1 1 -0.597469 0.599759 -0.597469
2 1 -0.224403 -0.224403 0.219674
3 1 -0.597469 -0.597469 0.599759
4 1 -0.224403 0.219674 -0.224403
5 1 -0.803306 -0.803306 -0.803306
6 1 0.219674 -0.224403 -0.224403
7 1 0.599759 -0.597469 -0.597469
8 1 1.62762 1.62762 1.62762

Stress (bar):
Pxx Pyy Pzz Pxy Pxz Pyz 
12376.32     12376.32     12376.32   -12024.359   -12024.359   -12024.359