REVIEW 2 major objections 4 minor 32 references
Defining velocities for accurate kinetic statistics in the GJF thermostat
T0 review · 2 major / 4 minor · reviewed 2026-08-14 · deepseek-v4-flash
Pith's one-line read For the GJF thermostat, no on-site velocity can give time-step-independent kinetic temperature, but a one-parameter family of two-point velocities can.
desk verdict A solid extension of the GJF kinetic-velocity program; the one-parameter family is the real result, and the no-on-site claim needs a qualifier. read the letter →
The pith
A machine-rendered reading of the paper's core claim, the machinery that carries it, and where it could break.
The reading
What carries the argument
The engine of the argument is the general finite-difference velocity ansatz, Eq. (18), with five unitless coefficients ($\gamma_1$ through $\gamma_5$) spanning all three-point position differences plus the two noise increments that bracket a time step. Plugging the GJF trajectory correlations into the mean-square velocity produces Eq. (27), whose requirement of time-step independence imposes algebraic constraints. For the two-point subclass the condition collapses to a single equation whose solution parameterizes the whole family by $\gamma_1$. The analysis of on-site velocities turns on the limiting behavior of $\gamma_3$ as $\alpha dt\to 0$, where the divergence in the noise coefficient $\gamma_4$ eliminates the possibility.
What would settle it
Take the on-site ansatz $w=(r_{n+1}-r_{n-1})/(2dt)+\lambda\beta_n/m$ with a finite, fixed $\lambda$, and insert the GJF trajectory correlations into the mean-square velocity; if any finite $\lambda$ makes $\langle w^2\rangle=k_B T/m$ for all stable $\Omega_0 dt$, the no-go claim for on-site velocities is wrong. Conversely, testing any two-point member of the family on a strongly anharmonic oscillator and seeing a systematic time-step trend in $\langle w^2\rangle$ larger than statistical error would bound how far the linear-system proof extends.
Extended reading notes
Core claim
The central claim is that kinetic statistics can be made exactly consistent with the GJF trajectory by choosing the right velocity definition. For the linear harmonic oscillator $f=-\kappa r$, the paper derives the general three-point finite-difference velocity $w=(\gamma_1 r_{n+1}+\gamma_2 r_n+\gamma_3 r_{n-1})/dt+(\gamma_4 \beta_n+\gamma_5 \beta_{n+1})/m$ and imposes that $\langle w w\rangle=k_B T/m$ for every stable time step. This forces the coefficients of $(\Omega_0 dt)^{-2}$ and $(\Omega_0 dt)^2$ in the mean-square velocity to vanish. The inverse-squared term rules out on-site velocities: in the frictionless limit any on-site velocity must approach the central difference $\gamma_1=1/2$, $\gamma_3=-1/2$, which drives the required noise coefficient to infinity. For two-point velocities ($\gamma_3=\gamma_4=0$), the remaining condition is $b\gamma_1^2+4(1-b)\gamma_5^2/b+4(1-b)\gamma_1\gamma_5=1$, yielding a one-parameter family $\gamma_5(\gamma_1)$ with $|b\gamma_1|\le 1$. Three highlighted members are Case A ($\gamma_1=1/\sqrt{b}$, $\gamma_5=0$, the 2GJ velocity), Case B ($b\gamma_1=1$, $\gamma_5=-1/2$), and Case C ($\gamma_1=1$). The paper also provides a leap-frog algorithm update that realizes any member and notes that configurational sampling remains unchanged.
Load-bearing premise
The results assume that a kinetically correct velocity must match the standard central-difference or leap-frog form when friction vanishes, and that matching equipartition for the harmonic oscillator at every stable time step is the defining test of correct kinetic statistics.
Editorial extensions
If this is right
- Users of the GJF thermostat can now choose a velocity definition matched to their measurement: Case A for equilibrium diffusive transport, Case C for drift or ballistic motion, and Case B as a maximal-amplitude option.
- All members of the family give exactly the same configurational sampling as GJF, since they are built on the same position trajectory; kinetic and configurational temperatures agree across the full stability range for linear systems.
- The family subsumes the two previously reported kinetically correct velocities, showing they are special cases of one continuous parameter.
- Lennard-Jones simulations across two thermodynamic states and three friction values confirm that the highlighted definitions measure kinetic temperature with time-step-independent averages, while the on-site GJF velocity departs increasingly with time step.
- Green-Kubo diffusion evaluation changes with the chosen velocity: Case A requires a right-Riemann sum and Case B a trapezoidal sum to reproduce the Einstein diffusion coefficient.
Reading between the lines
- The one-parameter freedom is likely a resource for tuning higher-order statistics: members of the family share the same mean kinetic energy but differ in velocity autocorrelation and fluctuations, so future work could pick $\gamma_1$ to minimize variance of kinetic temperature estimates.
- A similar no-go argument may apply to any integrator whose natural velocity is on-site and must match the central-difference limit, not just GJF; the divergence mechanism is generic.
- The linear-system result rigorously covers harmonic potentials, but the paper's Lennard-Jones tests suggest the family may work for nonlinear systems too; a systematic nonlinear proof or a counterexample remains open.
- For non-equilibrium simulations with a drift velocity, Case C's property of giving the correct ballistic velocity may make it preferable even though it is not a true half-step velocity; this distinction is a practically testable prediction.
Signed reviews
Editorial analysis
A structured set of objections, weighed in public.
Referee Report
Summary. The paper studies possible velocity definitions that can accompany the Grønbech-Jensen-Farago (GJF) thermostat. Starting from the GJF configurational trajectory, the authors write a general finite-difference velocity involving three consecutive positions and the two adjacent noise increments [Eq. (18)], impose the requirement that the mean-square velocity equal the equipartition value k_B T/m for every stable time step in a harmonic oscillator, and solve the resulting constraints. They conclude that no on-site velocity in this class can have time-step-independent kinetic temperature, while a one-parameter family of two-point velocities can; they give the explicit parametrization [Eqs. (31)-(34)], express the resulting leap-frog algorithms [Eqs. (35)-(41)], identify Cases A, B, and C, and test them on Lennard-Jones solids and liquids. The central derivation for the two-point family is explicit and internally consistent, and the simulations corroborate the predicted time-step independence for the highlighted choices.
Significance. The central constructive result is a useful unification: one parameter γ1 generates infinitely many two-point velocities with correct, time-step-independent mean kinetic energy for linear systems, including the previously known 2GJ half-step velocity (Case A) and the Farago velocity (Case B). The analytic derivation is not fitted to data, the algorithmic forms in Eqs. (35)-(41) are directly implementable, and the paper is careful to note that correct kinetic energy alone is not sufficient (Case D) and to give Green-Kubo consistency criteria. If the construction holds, the paper provides practical guidance for choosing a velocity estimator in GJF-based Langevin simulations. The main weakness is the advertised no-go statement for on-site velocities, which is narrower than the abstract suggests and needs qualification.
major comments (2)
- [Abstract and §II.B (after Eq. (29))] The no-go result for on-site velocities is proved only under the convention, stated in the text after Eq. (29), that a 'reasonable' on-site velocity must approach the central-difference stencil γ1=1/2, γ3=−1/2 as αdt→0. The abstract and the first sentence of Section IV state the impossibility without this qualifier. The convention is load-bearing: with only first-order consistency required, the factor bγ1 in Eq. (29) can vanish. For example, γ1=0, γ2=1, γ3=−1, γ5=0 makes Eq. (27) reduce to b + 4(1−b)γ4^2/b + 4(1−b)γ4 = 1, which has the finite root γ4=(−b+√(b^2+b))/2 for every b∈(0,1). This is a backward two-point velocity, a time-shifted member of the family used for the paper's positive result, with finite coefficients and no diverging noise term. The conclusion should therefore be reworded as a conditional statement, or proved without the central-difference limiting assumption.
- [Abstract and §II.C (Eq. (18))] The claim of a 'complete investigation of all possible finite difference approximations' exceeds what is shown. The analysis covers only the three-position, two-noise ansatz of Eq. (18); velocities built from more than three consecutive positions or from additional noise increments (for example β_{n−1}) are not treated. Both main conclusions should be stated as applying within the family Eq. (18), unless a supplementary argument establishes that no wider ansatz can work.
minor comments (4)
- [Eq. (33)] The notation for the square root is ambiguous; the expression should be written explicitly as sqrt{b(1−b^2 γ1^2)/(1−b)} so that the restriction (bγ1)^2≤1 follows immediately.
- [References] Refs. [27]-[29] are cited with DOIs and arXiv identifiers instead of complete bibliographic entries; please complete them in the journal's reference format.
- [Global] The text contains typographical artifacts, for example an extra space in 'G JF' in the title line and an awkward line break inside 'kinetic averages' in Section I; a careful proofread would improve readability.
- [Section III] The sentence describing the Lennard-Jones kinetic results as 'as expected from the analysis above' extrapolates an exact harmonic-oscillator result to a nonlinear system; the agreement is numerical validation rather than a consequence of the analytic derivation, and the text should say so explicitly.
Circularity Check
No significant circularity: the two-point velocity family is derived analytically from the stated equipartition constraint and the GJF trajectory, with prior GJF results forming independently validated inputs.
full rationale
After walking the derivation chain, I find no circular step. The paper's new content is a classification exercise: with the GJF trajectory (Eq. (16)) and the general three-point velocity ansatz (Eq. (18)) as inputs, the equipartition condition <ww>=kBT/m (Eq. (21)) is inserted into the correlation expansion (Eq. (27)); requiring the result to hold for every stable time step forces the coefficient conditions (Eqs. (28)-(29)), which in the two-point sector gamma3=gamma4=0 yields the explicit one-parameter family (Eqs. (33)-(34)). This is an analytic constraint-solution, not a fit: no parameter is tuned to the kinetic-energy data that the family is later said to produce. Cases A, B, and C are specializations of the solved family, and their Lennard-Jones simulations test the nonlinear regime not covered by the linear derivation, so those comparisons are not forced by construction. The prior GJF configurational-sampling results are invoked as inputs and are independently published and reproduced, including by non-overlapping authors (Ref. [28]); the paper also reruns configurational benchmarks in Section III. The one caveat worth flagging is Section II.B's on-site no-go: it is proved only under the stated central-difference limiting condition gamma1=1/2, gamma3=-1/2 as alpha dt -> 0, while the abstract and Section IV state it unconditionally. That is a precision/assumption issue, not circularity, because the proof does not assume the impossibility it concludes. Score 0.
Assumptions & free parameters
free parameters (1)
- gamma_1
assumptions (5)
- domain assumption Equilibrium correlations of the linear GJF oscillator, Eqs. (14) and (23)-(26), are correct and taken from prior work.
- ad hoc to paper A candidate velocity has the finite-difference form Eq. (18): three consecutive positions plus the two noise terms beta_n and beta_{n+1}.
- standard math Kinetic correctness means <w^2> = k_B T/m for every stable time step.
- domain assumption In the zero-friction limit the velocity must reduce to the standard central-difference or forward-difference form.
- domain assumption Exact kinetic correctness for linear systems approximately carries over to nonlinear systems.
Cite this review
Pith. "Pith review of Defining velocities for accurate kinetic statistics in the GJF thermostat." pith.science (2026). https://pith.science/paper/X7BWNS5O
@misc{pith2026190809225,
author = {Pith},
title = {Pith review of: Defining velocities for accurate kinetic statistics in the GJF thermostat},
year = {2026},
howpublished = {\url{https://pith.science/paper/X7BWNS5O}},
note = {Machine review of arXiv:1908.09225}
}
read the original abstract
We expand on two previous developments in the modeling of discrete-time Langevin systems. One is the well-documented Gr{\o}nbech-Jensen Farago (GJF) thermostat, which has been demonstrated to give robust and accurate configurational sampling of the phase space. Another is the recent discovery that also kinetics can be accurately sampled for the GJF method. Through a complete investigation of all possible finite difference approximations to the velocity, we arrive at two main conclusions:~1) It is not possible to define a so-called on-site velocity such that kinetic temperature will be correct and independent of the time step, and~2) there exists a set of infinitely many possibilities for defining a two-point (leap-frog) velocity that measures kinetic energy correctly for linear systems in addition to the correct configurational statistics obtained from the GJF algorithm. We give explicit expressions for the possible definitions, and we incorporate these into convenient and practical algorithmic forms of the normal Verlet-type algorithms along with a set of suggested criteria for selecting a useful definition of velocity.
Figures
Reference graph
Works this paper leans on
-
[1]
Moreover, in the limit αdt → 0 ( a, b → 1, β n = β n+1 = 0) the coefficient γ1 should become either γ1 → 1 2 or γ1 → 1, such that w becomes one of the two known velocities given in Eqs. (19) and (20) in that limit. Under these conditions, Eq. (29) yields the following noise term associated with β n: γ4β n = b 1 − b γ3 √ 2αk BT dt σ n = 2mγ3 √ 2 kBT αdt σ n,...
-
[2]
3E0, which results in a stable fcc (face centered cubic) crystal at a volume of V = 617
Two characteristically different tem- peratures are tested; kBT = 0 . 3E0, which results in a stable fcc (face centered cubic) crystal at a volume of V = 617 . 2558r3 0, and kBT = 0 . 7E0, which results in a liquid state at volume V = 824 . 9801r3
-
[3]
1: Statistical averages of potential energy ⟨Ep⟩ (from Eq
We show re- sults for three different friction coefficients: α = 1 mω 0, 6 FIG. 1: Statistical averages of potential energy ⟨Ep⟩ (from Eq. (52)) (a) and (b), and its standard deviation σ p (from Eq. (53)) (c) and (d) as a function of reduced time step ω 0dt for α = 1 mω 0. (a) and (c) show results for a crystalline fcc state at kBT = 0 . 3E0; (b) and (d) sho...
-
[4]
M. P . Allen, D. J. Tildesley, Computer Simulation of Liquids, Oxford University Press, Inc., New York, 1989
work page 1989
-
[5]
D. Frenkel and B. Smit, Understanding Molecular Sim- ulations: From Algorithms to Applications , (Academic Press, San Diego, 2002)
work page 2002
-
[6]
D. C. Rapaport, The Art of Molecular Dynamics Simu- lations, (Cambridge University Press, Cambridge, 2004)
work page 2004
- [7]
- [8]
Show all 32 references
-
[9]
W. G. Hoover, Phys. Rev. A 31, 1695 (1985)
1985
-
[10]
Schneider and E
T. Schneider and E. Stoll, Phys. Rev. B 17, 1302 (1978)
1978
-
[11]
Br¨ unger, C
A. Br¨ unger, C. L. Brooks, and M. Karplus, Chem. Phys. Lett. 105, 495 (1984)
1984
-
[12]
R. W. Pastor, B. R. Brooks, and A. Szabo, Mol. Phys. 65, 1409 (1988)
1988
-
[13]
W. F. van Gunsteren, H. J. C. Berendsen, Mol. Simul. 1, 173, (1988)
1988
-
[14]
R. J. Loncharich, B. R. Brooks, and R. W. Pastor, Biopolymers 32, 523 (1992)
1992
-
[15]
Vanden-Eijnden, G
E. Vanden-Eijnden, G. Ciccotti, Chem. Phys. Lett. 429, 310 (2006)
2006
-
[16]
Melchionna, J
S. Melchionna, J. Chem. Phys. 127, 044108 (2007)
2007
-
[17]
Leimkuhler, C
B. Leimkuhler, C. Matthews, Appl. Math. Res. Express 2013, 34 (2012)
2012
-
[18]
Grønbech-Jensen and O
N. Grønbech-Jensen and O. Farago, Mol. Phys. 111, 983 (2013)
2013
-
[19]
Paquet, H
E. Paquet, H. L. Viktor, BioMed Res. Int., 183918 (2015)
2015
-
[20]
Langevin, C
P. Langevin, C. R. Acad. Sci. Paris 146, 530 (1908)
1908
-
[21]
Parisi, Statistical Field Theory , (Addison-Wesley, Menlo Park, 1988)
G. Parisi, Statistical Field Theory , (Addison-Wesley, Menlo Park, 1988)
1988
-
[22]
L. F. G. Jensen and N. Grønbech-Jensen, Mol. Phys. 117, 2511 (2019)
2019
-
[23]
O. G. Jepps, G. Ayton, and D. J. Evans, Phys. Rev. E 62, 4757 (2000)
2000
-
[24]
Rickayzena and J
G. Rickayzena and J. G. Powles, J. Chem. Phys. 114, 4333 (2001)
2001
-
[25]
P. J. Davis, B.A. Dalton, and T. Morishita, Phys. Rev. E 86, 056707 (2012)
2012
-
[26]
Jackson, J
N. Jackson, J. Miguel Rubi, and Fernando Bresme, Mol. Sim. 42, 1214 (2016)
2016
-
[27]
Grønbech-Jensen, N
N. Grønbech-Jensen, N. R. Hayre, and O. Farago, Com- put. Phys. Commun. 185, 524 (2014)
2014
-
[28]
E. Arad, O. Farago, and N. Grønbech-Jensen, Isr. J. Chem. 56, 629 (2016)
2016
-
[29]
Farago, Physica A 534, 122210 (2019)
O. Farago, Physica A 534, 122210 (2019)
2019
-
[30]
L. F. G. Jensen and N. Grønbech-Jensen, Comput. Phys. Commun. (2019), DOI: 10.1016/j.cpc.2019.107011. arXiv:1902.02338
2019
-
[32]
Grønbech-Jensen, Mol
N. Grønbech-Jensen, Mol. Phys. (2019), DOI: 10.1080/00268976.2019.1662506. arXiv:1909.04380
2019
-
[33]
M. S. Green, Chem. Phys. 22, 398 (1954); R. Kubo, J. Phys. Soc. Jpn. 12, 570 (1957). Appendix A: Green-Kubo diffusion using the GJF on-site velocity vn for f = 0 The Green-Kubo equivalent of the Einstein expres- sion for diffusion Eq. (12) is in discrete-time calculated from the...
1954
Reviewed August 14, 2026 · model on record in the stance chip above.
Discussion (0). Continue with ORCID to comment.