Make omega_mass proportional to (N+1)kT
According to papers like Martyna, Tobias, Klein (JCP 1994, doi:10.1063/1.467468 section II.F) and Martyna, Tuckerman, Tobias, Klein (Mol. Phys. 1996, doi:10.1080/00268979600100761 section 2.5), the mass of the cell parameters should be proportional to (Ndof + dim) / dim, or in other words, Natoms + 1.
This commit is contained in:
@ -798,7 +798,7 @@ void FixNH::setup(int vflag)
|
|||||||
|
|
||||||
if (pstat_flag) {
|
if (pstat_flag) {
|
||||||
double kt = boltz * t_target;
|
double kt = boltz * t_target;
|
||||||
double nkt = atom->natoms * kt;
|
double nkt = (atom->natoms + 1) * kt;
|
||||||
|
|
||||||
for (int i = 0; i < 3; i++)
|
for (int i = 0; i < 3; i++)
|
||||||
if (p_flag[i])
|
if (p_flag[i])
|
||||||
@ -1826,7 +1826,7 @@ void FixNH::nhc_press_integrate()
|
|||||||
// Update masses, to preserve initial freq, if flag set
|
// Update masses, to preserve initial freq, if flag set
|
||||||
|
|
||||||
if (omega_mass_flag) {
|
if (omega_mass_flag) {
|
||||||
double nkt = atom->natoms * kt;
|
double nkt = (atom->natoms + 1) * kt;
|
||||||
for (int i = 0; i < 3; i++)
|
for (int i = 0; i < 3; i++)
|
||||||
if (p_flag[i])
|
if (p_flag[i])
|
||||||
omega_mass[i] = nkt/(p_freq[i]*p_freq[i]);
|
omega_mass[i] = nkt/(p_freq[i]*p_freq[i]);
|
||||||
|
|||||||
Reference in New Issue
Block a user