{"id":"e6be185c-3b36-4008-981e-5995e64ead07","arxiv_id":"1908.08617","paper_version":1,"verdict":"CONDITIONAL","confidence":"MODERATE","novelty_score":5.0,"correctness_risk":"medium","formal_verification":"none","parameter_count":4,"one_line_summary":"Classical stochastic MD with explicit gold layers reproduces ~20 pW/K thermal conductance in alkanedithiol single-molecule junctions and shows harmonic Landauer calculations miss anharmonic effects at room temperature.","lead":"The authors simulate heat flow through single molecules with a stochastic molecular dynamics framework that filters thermal noise through explicit gold layers. They find this reproduces measured alkane conductances around 20 pW/K and shows conjugated polyynes conduct less heat than saturated alkanes.","discovery_kind":"new_application","skeptic_critique":{"model":"deepseek-v4-flash","headline":"Untested Langevin coupling gamma may set the computed conductance; the reported 20 pW/K match with experiment and the MD-vs-Landauer anharmonicity claim are not robust without a gamma-scan.","rationale":"The reader's stated weakest assumption is that Jtot is independent of the partition weights C_tau,j and the imaginary plane. This is a natural thing to worry about, but it is actually protected by the structure of Eqs. (1)-(3): because sum_j C_tau,j = 1, the pairwise fluxes J_tau,ij are antisymmetric and satisfy dE_i/dt = sum_j J_ij + bath power. In a steady state, summing over one side of any dividing plane gives the current equal to the power injected by the baths on that side, independent of C and plane. So I do not think the partition/plane issue is the load-bearing risk; a bug in implementing Eq. (3) would show up as a violation of this balance, but the manuscript does not provide code or a self-check. The more serious gap is the untested bath coupling gamma. The MD conductance is computed with gamma = 1 ps^-1, but no sensitivity study is shown. In a finite model the bath-cluster interface is in series with the molecule, and the conductance can be limited by that interface; the four-layer convergence test of Fig. 4 is consistent with either the explicit gold layers correctly filtering the bath or the bath contact dominating for every layer count. The Landauer calculations in the SI likewise depend on an unreported gamma_L/R. Since the paper's physical conclusion is that classical MD matches experiment better than the harmonic Landauer result (and hence that anharmonicity is relevant), a missing gamma in either channel makes that comparison untestable. This supports the conditional verdict: addressable by a gamma scan and by reporting the Landauer gamma, but not a demonstrated error.","tokens_in":14884,"tokens_out":20827,"duration_ms":209064,"concrete_test":"Perform a gamma-scan for HS(CH2)6SH with the reported three-layer gold geometry: run SNEMD with gamma = 0.1, 0.3, 1, 3, and 10 ps^-1 on the outer gold layer, keeping force field, geometry, bias, and time step fixed, and compute the steady-state conductance with the paper's plane-current definition. If the conductance changes by more than the published SE (about 0.74 pW/K) across this range, the central 20 pW/K value is not robust to the bath-coupling parameter. In parallel, recompute the Landauer transmission of the same junction with the same gamma values to determine whether the MD-vs-Landauer gap, and the anharmonicity inference, survives.","verdict_should_be":"UNCHANGED","load_bearing_attack":"The SNEMD method has a free parameter that is never tested: the Langevin friction coefficient gamma of the Markovian baths attached to the outermost gold layer. The paper fixes gamma = 1 ps^-1 (Section 2) but reports no scan over gamma. In a finite junction, gamma controls the bath-electrode contact resistance; for small gamma the measured current is set by the bath coupling rather than by the molecule, and the layer-convergence test in Fig. 4 cannot distinguish 'the gold layers filter the bath' from 'the bath contact dominates regardless of layer count.' If the 20 pW/K value moves with gamma, the agreement with Cui et al. is an artifact of the chosen coupling. The same gap affects the Landauer comparison: SI Eqs. (S1)-(S5) use white baths with constant gamma_L/R, but the value is not reported. The conclusion that classical MD fits experiment better than harmonic Landauer, and hence that anharmonicity matters, is unverifiable unless both calculations use the same, physically justified gamma. The reader's plane-independence concern is less likely to be the soft spot: Eqs. (1)-(3) satisfy local energy conservation for any partition with sum_j C_tau,j = 1, so in steady state the current across a plane equals the bath power on one side and is independent of C and plane, provided Eq. (3) is implemented correctly. A numerical check is still useful, but the untested bath coupling is the load-bearing gap.","agreement_with_reader":"partial"},"referee_report":{"model":"deepseek-v4-flash","summary":"The authors present a stochastic nonequilibrium molecular dynamics (SNEMD) framework for phononic heat transport in single-molecule junctions, in which Markovian Langevin baths are attached to the outermost of several explicit gold layers that filter the bath spectrum. The heat current is computed by partitioning interaction energies among atoms and summing pairwise fluxes across an imaginary plane. The method is applied to alkanedithiols and polyynes with UFF force fields, and the computed thermal conductances (~20 pW/K for alkanes) are compared with recent experimental values (Cui et al., Nature 2019) and with harmonic Landauer calculations based on the same force field. The paper reports that classical MD agrees better with experiment than the harmonic Landauer calculation, suggests that anharmonicity plays a role, and finds that polyynes have lower thermal conductance than alkanes despite their higher electronic conductance.","tokens_in":15219,"tokens_out":7686,"duration_ms":75016,"significance":"If the reported values are robust, the SNEMD approach provides a practical, generally applicable tool for simulating heat transport in molecular junctions, with the explicit-gold-layer filtering being a useful methodological refinement over direct Markovian baths. The paper is careful in reporting statistical errors and in checking convergence with the number of gold layers, and it compares against both independent ab initio Landauer results and experimental data. The polyyne result is an interesting falsifiable prediction about conjugation and thermal transport. However, the central quantitative claims rest on untested choices of the bath coupling and are partly confounded by classical/quantum statistical differences, so the significance is conditional on those points being resolved.","major_comments":[{"comment":"The Langevin friction coefficient γ is a free parameter that is fixed to a single value (apparently 1 ps⁻¹, stated in a garbled sentence in Section 2) in the MD simulations, but no scan over γ is reported. In a finite junction with explicit gold layers, γ controls the bath–electrode contact resistance; for small γ the measured current is set by the bath coupling rather than by the molecule, and the layer-convergence test in Fig. 4 cannot distinguish 'the gold layers filter the bath' from 'the bath contact dominates regardless of layer count'. The value of γ used in the Landauer calculation (SI Eqs. S1–S5) is also not reported, so the comparison in Fig. 8 and the conclusion that classical MD captures anharmonicity better than harmonic Landauer are not reproducible. A γ-scan (e.g. 0.1–10 ps⁻¹) with a statement of the plateau range is required before the 20 pW/K agreement with experiment can be considered robust.","section":"Section 2 (Model and Calculations) and SI Eqs. (S4)-(S5)"},{"comment":"The claim that the better agreement of MD with experiment 'indicat[es] anharmonicity plays a tunning role' is confounded by the different statistics of the two methods. The MD simulation is classical (equipartition), while the Landauer calculation uses Bose–Einstein occupancies (SI Eq. S1). For modes with ħω ≲ kBT the classical and quantum populations differ, and the high-frequency modes that are suppressed in the quantum calculation may still carry some current in the MD if the gold spectral density has weight there. The observed MD–Landauer difference could therefore reflect classical vs. quantum statistics rather than anharmonicity. To support the anharmonicity conclusion, the authors should compare MD with a classical Landauer calculation (using the classical limit of the Bose–Einstein distribution, kBT/ħω), or provide a frequency-resolved decomposition of the MD current, or perform an MD run with the harmonic part of the force field only.","section":"Section 3 (Results and Discussion) and Conclusion"},{"comment":"The statement that 'the heat current through the molecule will be measured the same, regardless of where we chose to draw this imaginary plane' is asserted but never numerically verified. While the local energy-conservation argument with Σ_j C_{τ,j}=1 makes plane independence plausible in the continuum limit, the implementation with multi-body force-field terms and discrete atoms warrants a check. A simple test is to compute the total current across planes located at different positions along the molecule (e.g., at each backbone atom) and show that the values agree within statistical error. If the current varies with plane position, the reported conductance is an artifact of the chosen flux decomposition rather than a physical observable.","section":"Section 2, paragraph after Eq. (4)"}],"minor_comments":[{"comment":"There are numerous typographical errors, including 'Nonequlibrium', 'alkanethiol' (SI), 'molecualr', 'polyyenes', 'utlizing', 'tunning', and 'calcualted'; these should be corrected throughout.","section":"Abstract and throughout"},{"comment":"The expression for Welch's t-test, η(x,y) = |E[x]-E[y]|/√(σ_x²+σ_y²), is a t-statistic, not a p-value; reporting '13%' as the result of a t-test is ambiguous and should be clarified as a p-value with the degrees of freedom.","section":"Table 1 caption"},{"comment":"The sentence 'the coupling strength between the Markovian bath and outermost layer of explicit bulk (region is 1) is ps−1' is garbled; the value and units of γ should be stated clearly in a complete sentence.","section":"Section 2 (Model and Calculations)"},{"comment":"The statement that 'the similarity of the results indicate that heat transport is dominated by modes in the lower frequency range' applies cleanly only to the polyyne series; for the alkanes the MD and Landauer results differ markedly in Fig. 8, so the claim should be qualified by chain type.","section":"Abstract and Section 3"},{"comment":"No data or code availability statement is provided; given that the method is implemented in a customized GROMACS package, providing the input topologies and analysis scripts would enhance reproducibility.","section":"Supporting Information"},{"comment":"The legend text 'the ones with and without hydrogen atoms are alkanes and conjugated polyyenes respectively' is confusing; it should say 'saturated alkanedithiols and polyynes'.","section":"Figure 7 caption"}],"recommendation":"major_revision","confidential_remarks":"The paper is from an established group and the method is potentially useful, but the main quantitative claims are conditional on the missing γ-scan and on disentangling anharmonicity from classical-statistics effects. Both issues are fixable within the scope of a revision. I would also encourage the editor to ask for a clear data/code availability statement, as the customized GROMACS implementation is central to the method's reproducibility."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"Colleague,\n\nThis paper is worth reading, but the main quantitative claim has a gap. The SNEMD setup—white-noise Langevin baths filtered through explicit gold layers—is a sensible way to model molecular heat conduction, and the computed alkane conductances clustering around 20 pW/K are pleasing. But the bath friction gamma is fixed at 1 ps^-1 with no scan, and the Landauer calculation in the SI uses the same type of white bath without reporting its gamma. If the computed current depends on gamma, the agreement with the Cui et al. experiment is an artifact of that choice, and the conclusion that classical MD outperforms harmonic Landauer at room temperature—attributed to anharmonicity—does not yet hold.\n\nThe paper does several things well. The extended bath model is a genuine improvement over bare Nose-Hoover or simple thermostats. Using the same UFF force field for the classical MD and for the Hessian in the Landauer calculation is a fair consistency check. The layer-convergence test in Fig. 4 is a good sanity check, as is the all-atom vs. unified-atom comparison. The polyyne result—lower thermal conductance than saturated alkanes despite higher electronic conductance—is a useful datapoint for thermoelectric design. And the paper does not fit its conductance to the experimental target, which deserves credit.\n\nThe stress-test note is on target: the untested gamma is a bigger problem than the arbitrary partition of potential energy among atoms. The partition weights sum to one, so local energy conservation should make the steady-state current across a plane independent of those weights, provided the implementation is consistent. Still, the paper asserts plane-independence without testing it, and a numerical check would be cheap. More important, the paper ships no code, does not state the restraint force constant, and does not report the Landauer gamma_L/R. Those are ordinary omissions, but they matter here because the central number may be sensitive to them.\n\nNone of this suggests the results are wrong. It means they are not yet verifiable. The framework is likely useful and the comparison between classical anharmonic MD and harmonic Landauer is a question people in the field care about. I would send this to referees with a clear request: show a gamma scan, report every parameter, and test plane-independence numerically.\n\nBest,\n[You]","headline":"A useful SNEMD framework for molecular junctions, but the missing gamma-scan and unreported Landauer coupling leave the headline MD-vs-Landauer claim unverified.","tokens_in":15695,"tokens_out":3547,"would_cite":false,"duration_ms":38547,"reading_group":"maybe","serious_thinker":"yes","would_accept_peer_review":true},"rs_alignment":null,"lean_confirmation":null,"pith_extraction":{"msc":[],"pacs":[],"model":"deepseek-v4-flash","headline":"Classical molecular dynamics with explicit gold layers reproduces the measured ~20 pW/K heat conductance of single-molecule alkane junctions.","keywords":["stochastic nonequilibrium molecular dynamics","thermal conductance","single-molecule junctions","phonon heat transport","Langevin baths","alkanedithiols","polyynes","anharmonicity"],"falsifier":"Recompute the same alkanedithiol junctions with strongly unequal partition weights, for example putting 90 percent of each stretching energy on one atom, or draw the counting plane through the gold layers instead of at the molecule-substrate interface; if the steady-state current and the 20 pW/K plateau change by more than the statistical error, the central claim collapses.","tokens_in":14701,"feed_emoji":"🔥","tokens_out":9400,"duration_ms":91771,"temperature":0.7,"pith_summary":"This paper tries to establish that phononic heat conduction through a single molecule between two gold electrodes can be captured by classical stochastic nonequilibrium molecular dynamics, without tuning a thermostat directly on the molecule. The key move is to couple white-noise Langevin baths only to the outermost atoms of several explicit gold layers, letting the metal layers filter the noise into a realistic vibrational spectrum before it reaches the molecule. For alkanedithiol single-molecule junctions the simulations give a thermal conductance that saturates near 20 pW/K as more gold layers are added, in line with recent experiments, and at room temperature the classical result tracks the measured behavior better than harmonic quantum Landauer calculations built from the same force field. The paper reads this as evidence that low-frequency modes dominate molecular heat conduction and that anharmonicity plays a real role at ambient conditions, not just a quantum-mechanical correction.","feed_headline":"Simulation reproduces single-molecule heat flow at 20 pW/K","feed_subtitle":"Stochastic molecular dynamics with explicit gold leads lands on the experimental thermal conductance for alkanedithiol junctions.","key_machinery":"The machinery is the atom-resolved heat-current decomposition combined with an extended-bath thermostat. Each potential-energy term $V_\\tau$ is partitioned among the atoms it connects by weights $C_{\\tau,j}$ (equal weights here), so the flux from atom $j$ to atom $i$ through interaction $\\tau$ is $J_{\\tau,ij}=C_{\\tau,j}\\mathbf{f}_{\\tau,i}\\cdot\\mathbf{v}_i - C_{\\tau,i}\\mathbf{f}_{\\tau,j}\\cdot\\mathbf{v}_j$; summing these across an imaginary plane perpendicular to the chain gives the molecular heat current. The reservoirs enter as Langevin white noise applied only to the outermost gold layer, while the explicit gold layers in between provide a realistic spectral density that the molecular vibrations actually see. The assumption that the current is independent of the partition weights and of the plane position lets the authors report a single conductance $\\kappa = J/(T_{\\rm hot}-T_{\\rm cold})$.","core_discovery":"The central claim is that a stochastic nonequilibrium molecular dynamics (SNEMD) scheme with an explicitly modeled metal-molecule interface can serve as a general, atomistic tool for heat conduction in single-molecule junctions in and beyond linear response. Using the Universal Force Field and several layers of gold atoms as a filter between the molecule and Markovian reservoirs, the paper obtains steady-state heat currents and defines the junction conductance as the current divided by the 50 K bath bias. For alkanedithiols, the conductance saturates at about 20 pW/K once three gold layers are included, matching the most recent single-molecule experiments; for conjugated polyynes, the simulations find lower thermal conductance than for saturated alkanes of the same length, opposite to their electronic conductance ordering. A side-by-side comparison with harmonic Landauer-type calculations from the same force field shows that the classical MD results agree with experiment at room temperature, while the quantum harmonic results show non-monotonic length dependence, suggesting that anharmonicity and finite-bias effects are important in the experimentally accessed regime.","pith_inferences":["If the plane-independence assumption holds, the same flux decomposition could be applied per interaction type to quantify how much of the junction resistance is interfacial versus intramolecular, which the paper only sketches qualitatively through temperature profiles.","A decisive test of the anharmonicity interpretation would be to run the same SNEMD simulations at several bath biases, such as 10 K, 50 K, and 100 K; if conductance is bias-independent, the MD-versus-Landauer difference must come from anharmonicity rather than nonlinear response.","The framework's prediction that polyynes are low-heat-conducting but high-electron-conducting could be converted into a concrete thermoelectric figure-of-merit estimate only by also computing the electronic conductance and Seebeck coefficient for the same junction geometries."],"forward_implications":["Alkanedithiol single-molecule junctions should have a room-temperature thermal conductance near 20 pW/K that depends only weakly on chain length, rather than falling steeply.","Adding more than three explicit gold layers does not change the computed conductance, so the bath-filtering approximation is converged for gold leads.","At room temperature, classical full-force-field simulations are a more reliable route to junction thermal conductance than harmonic quantum Landauer calculations, which miss anharmonic contributions.","Conjugated polyyne molecular wires, despite conducting electrons well, conduct less heat than saturated alkanes, which would make them favorable for thermoelectric applications.","The same SNEMD framework can be applied to other molecules and metal leads to map structural dependence of interfacial heat conduction beyond linear response."],"supporting_citations":[{"why":"Supplies the recent experimental single-molecule thermal conductance values (about 20 pW/K for alkanedithiols) that the simulations are matched against.","marker":"37"},{"why":"Provides the ab initio Landauer-type calculations whose non-monotonic length dependence is compared with the present MD and Landauer results.","marker":"33"},{"why":"Establishes the expectation that low-frequency modes dominate molecular heat conduction, the theoretical baseline for the classical approximation.","marker":"32"},{"why":"Supplies the scheme for partitioning interaction potential energy among atoms, on which the heat-current decomposition is built.","marker":"40"},{"why":"Provides the Langevin-dynamics algorithm used to thermostat the outer bath layers.","marker":"31"},{"why":"Defines the Universal Force Field used for the molecular junction and gold lead interactions.","marker":"44"},{"why":"Supports the finding that only the first few gold layers affect the molecular heat transport, motivating the layer-convergence test.","marker":"45"},{"why":"Provides the molecular dynamics engine the custom simulation code is built around.","marker":"41"}],"fun_headline_variants":["Stochastic MD matches measured heat flow in single-molecule junctions","Simulation reproduces 20 pW/K thermal conductance for alkanedithiols","New molecular dynamics framework predicts nanoscale heat transport","Extending heat conduction simulations beyond linear response","Molecular junction heat flow reproduced by stochastic dynamics"],"cache_read_input_tokens":3200,"weakest_assumption_plain":"In Section 2 the paper assumes that the heat current through the molecule is the same regardless of where the imaginary counting plane is drawn, and the whole calculation also depends on the arbitrary choice of splitting each interaction's energy equally among the atoms involved; if the measured current changes with that bookkeeping, the reported conductance values are artifacts of the flux definition.","fun_headline_variants_meta":{"raw":{"variants":["Stochastic MD matches measured heat flow in single-molecule junctions","Simulation reproduces 20 pW/K thermal conductance for alkanedithiols","New molecular dynamics framework predicts nanoscale heat transport","Extending heat conduction simulations beyond linear response","Molecular junction heat flow reproduced by stochastic dynamics"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000289,"raw_usage":{"total_tokens":1691,"prompt_tokens":944,"completion_tokens":747,"prompt_tokens_details":{"cached_tokens":384},"prompt_cache_hit_tokens":384,"prompt_cache_miss_tokens":560,"completion_tokens_details":{"reasoning_tokens":666}},"tokens_in":560,"tokens_out":747,"duration_ms":6722,"temperature":1.0,"reasoning_tokens":666,"cache_read_input_tokens":384,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-14T11:34:29.937880+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Recompute the same alkanedithiol junctions with strongly unequal partition weights, for example putting 90 percent of each stretching energy on one atom, or draw the counting plane through the gold layers instead of at the molecule-substrate interface; if the steady-state current and the 20 pW/K plateau change by more than the statistical error, the central claim collapses.","supporting_citations":[{"cited_title":"A.; Kl \\\" o ckner, J","cited_arxiv_id":null,"evidence_quote":"Supplies the recent experimental single-molecule thermal conductance values (about 20 pW/K for alkanedithiols) that the simulations are matched against."},{"cited_title":"o ckner, J. C.; B \\","cited_arxiv_id":null,"evidence_quote":"Provides the ab initio Landauer-type calculations whose non-monotonic length dependence is compared with the present MD and Landauer results."},{"cited_title":"Thermal conductance through molecular wires","cited_arxiv_id":null,"evidence_quote":"Establishes the expectation that low-frequency modes dominate molecular heat conduction, the theoretical baseline for the classical approximation."},{"cited_title":"Contribution of inter- and intramolecular energy transfers to heat conduction in liquids","cited_arxiv_id":null,"evidence_quote":"Supplies the scheme for partitioning interaction potential energy among atoms, on which the heat-current decomposition is built."},{"cited_title":"J.; de Vries, A","cited_arxiv_id":null,"evidence_quote":"Provides the Langevin-dynamics algorithm used to thermostat the outer bath layers."},{"cited_title":"K.; Casewit, C","cited_arxiv_id":null,"evidence_quote":"Defines the Universal Force Field used for the molecular junction and gold lead interactions."},{"cited_title":"L.; Yan, T.; Hase, W","cited_arxiv_id":null,"evidence_quote":"Supports the finding that only the first few gold layers affect the molecular heat transport, motivating the layer-convergence test."},{"cited_title":"R.; Smith, J","cited_arxiv_id":null,"evidence_quote":"Provides the molecular dynamics engine the custom simulation code is built around."}],"review_version":1}