From 3286fcb13a99eb81a511d371e982eb5943fbc62a Mon Sep 17 00:00:00 2001 From: Aidan Thompson Date: Fri, 21 Jan 2022 19:29:23 -0700 Subject: [PATCH] Added correct units for pressure, still getting wrong answer --- examples/numdiff/in.numdiff | 5 ++--- src/EXTRA-FIX/fix_numdiff_stress.cpp | 5 ++++- 2 files changed, 6 insertions(+), 4 deletions(-) diff --git a/examples/numdiff/in.numdiff b/examples/numdiff/in.numdiff index f925aa8395..d451a56d1e 100644 --- a/examples/numdiff/in.numdiff +++ b/examples/numdiff/in.numdiff @@ -31,9 +31,8 @@ dump errors all custom 100 force_error.dump v_errx v_erry v_errz variable ferrsq atom (fx-f_numforce[1])^2+(fy-f_numforce[2])^2+(fz-f_numforce[3])^2 compute faverrsq all reduce ave v_ferrsq -fix avfaverrsq all ave/time 100 1 100 c_faverrsq ave running compute myvirial all pressure NULL virial - -thermo_style custom step temp pe press c_faverrsq f_avfaverrsq c_myvirial[1] f_numstress[1] +variable ratio11 equal f_numstress[1]/c_myvirial[1] +thermo_style custom step temp pe press c_faverrsq c_myvirial[*] f_numstress[*] v_ratio11 thermo 100 run 500 diff --git a/src/EXTRA-FIX/fix_numdiff_stress.cpp b/src/EXTRA-FIX/fix_numdiff_stress.cpp index 67ca091d61..a2dd5f9b8a 100644 --- a/src/EXTRA-FIX/fix_numdiff_stress.cpp +++ b/src/EXTRA-FIX/fix_numdiff_stress.cpp @@ -217,7 +217,10 @@ void FixNumDiffStress::calculate_stress() // compute a finite difference force in each dimension int flag,allflag; - double denominator = 0.5 / delta; + double nktv2p = force->nktv2p; + double inv_volume = 1.0 / (domain->xprd * domain->yprd * domain->zprd); + + double denominator = 0.5 / delta * inv_volume * nktv2p; for (int idir = 0; idir < NDIR_STRESS; idir++) { displace_atoms(nall, idir, 1.0);