virial adjustments for MDI stress

This commit is contained in:
Steve Plimpton
2022-07-27 09:38:23 -06:00
parent 10c8a41ad8
commit 925196a24f
3 changed files with 23 additions and 40 deletions

View File

@ -371,7 +371,8 @@ void FixMDIQM::post_force(int vflag)
}
}
// optionally request stress tensor from MDI engine, convert to virial
// optionally request stress tensor from MDI engine, convert to 6-value virial
// MDI defines virial tensor as intensive (divided by volume), LAMMPS does not
// qm_virial = fix output for global QM virial
if (virialflag) {
@ -607,18 +608,18 @@ void FixMDIQM::unit_conversions()
mdi2lmp_force = angstrom_to_bohr / ev_to_hartree;
}
// pressure or stress units = force/area = energy/volume
// stress units = force/area = energy/volume
mdi2lmp_pressure = 1.0;
lmp2mdi_pressure = 1.0;
if (lmpunits == REAL) {
lmp2mdi_pressure = (kelvin_to_hartree / force->boltz) /
(angstrom_to_bohr * angstrom_to_bohr * angstrom_to_bohr) / force->nktv2p;
(angstrom_to_bohr * angstrom_to_bohr * angstrom_to_bohr);
mdi2lmp_pressure = 1.0 / lmp2mdi_pressure;
} else if (lmpunits == METAL) {
lmp2mdi_pressure =
ev_to_hartree / (angstrom_to_bohr * angstrom_to_bohr * angstrom_to_bohr) / force->nktv2p;
ev_to_hartree / (angstrom_to_bohr * angstrom_to_bohr * angstrom_to_bohr);
mdi2lmp_pressure = 1.0 / lmp2mdi_pressure;
}
}

View File

@ -1493,6 +1493,9 @@ void MDIEngine::send_pe()
/* ----------------------------------------------------------------------
<STRESS command
send 9-component stress tensor (no kinetic energy term)
should be intensive quantity (divided by volume in pressure compute)
MDI stress tensor units are energy/volume,
so conversion factor includes nktv2p to convert pressure back to virial
---------------------------------------------------------------------- */
void MDIEngine::send_stress()
@ -1837,6 +1840,8 @@ void MDIEngine::unit_conversions()
}
// pressure or stress units = force/area = energy/volume
// MDI energy/volume = Hartree/Bohr^3,
// so need to remove LAMMPS nktv2p from pressure
mdi2lmp_pressure = 1.0;
lmp2mdi_pressure = 1.0;