{"id":"3628d56b-c730-452a-b2ce-aafa65f06285","arxiv_id":"2411.11988","paper_version":1,"verdict":"CONDITIONAL","confidence":"MODERATE","novelty_score":5.0,"correctness_risk":"medium","formal_verification":"none","parameter_count":3,"one_line_summary":"MixPI implements mixed-time-slicing path integral molecular dynamics in CP2K, allowing per-atom bead counts, and validates the implementation on water structure and a quantized-electron cobalt system.","lead":"This paper introduces MixPI, a CP2K-based program that lets each atom in a path integral molecular dynamics simulation use a different number of quantum beads. It benchmarks the program on liquid water and on a solvated cobalt ion with a 1024-bead electron ring polymer.","discovery_kind":"extension","skeptic_critique":{"model":"deepseek-v4-flash","headline":"Centroid-PME electrostatics makes the demonstrated benchmark agreement an internal consistency check rather than a validation of the exact mixTS Hamiltonian.","rationale":"The paper's central claim is that MixPI provides a general atomistic implementation of mixed-time-slicing PIMD. The theory in Section II is a plausible sequential Trotter factorization, and the implementation strategy of using CP2K as an external library is reasonable. However, the benchmarks do not yet establish that the code samples the exact mixTS Hamiltonian: the centroid-PME approximation replaces long-range electrostatic bead-bead interactions with centroid interactions, and both the all-replica reference and the mixTS runs use the same approximate force model. This makes Figs. 3-5 a test of algorithmic consistency rather than a test of the correctness of Eq. 13. The issue is especially pronounced for the 1024-bead electron: if electrostatic forces act only on the ring-polymer centroid, the external field cannot localize the electron distribution or couple its shape to the solvent, so the Co3+ + e- demonstration is weakened as evidence of a faithful quantized-electron simulation. Because the approximation is disclosed and drawn from prior ring-polymer contraction work, this is a validation gap rather than a demonstration of error in the mixTS formalism itself. The recommended verdict remains CONDITIONAL, pending an error estimate or a comparison with exact bead-resolved electrostatics.","tokens_in":15403,"tokens_out":10546,"duration_ms":113545,"concrete_test":"Re-run the bulk-water all-replica (NH=NO=32) simulation with exact bead-resolved PME electrostatics, for example using standard all-bead CP2K or an i-PI driver, and compare the O-O, O-H, and H-H radial distribution functions to the centroid-PME all-replica result reported in Fig. 3. If the rdfs differ by more than the stated statistical uncertainty, the paper's benchmark validates only the centroid-force model rather than Eq. 13, and the headline claim needs qualification. A secondary check is to compare the electron bead density in the Co3+ + e- system between centroid-PME and exact bead-resolved electrostatics to test whether the observed delocalization is an artifact of the centroid approximation.","verdict_should_be":"UNCHANGED","load_bearing_attack":"The most load-bearing assumption is introduced in Section II after Eq. 13: the smooth-particle mesh Ewald electrostatics between path-integral beads is replaced by forces between ring-polymer centroids, with no error estimate given. This approximation is used in both benchmark systems, and because the same centroid-PME force model supplies the all-replica reference and the mixTS runs, the agreement in Figs. 3-5 only shows that MixPI reproduces its own approximate reference; it does not validate the exact mixTS-PIMD Hamiltonian in Eq. 13 for Coulomb interactions. The concern is acute in the Co3+ + e- example: electrostatic forces act only on the electron ring-polymer centroid, so the external Coulomb field does not reshape the 1024-bead distribution. The displayed delocalized electron and the Co-O rdf agreement with Co2+ may therefore largely follow from placing a -1 point charge at the Co3+ position, rather than from a faithful quantum treatment of the electron. The paper does flag the approximation and cites prior ring-polymer contraction work, so this is not an internal inconsistency, but no error estimate or comparison to exact bead-resolved electrostatics is provided, leaving the central software claim insufficiently validated.","agreement_with_reader":"agree"},"referee_report":{"model":"deepseek-v4-flash","summary":"The manuscript introduces MixPI, a standalone driver that uses CP2K as a force-evaluation library to perform path integral molecular dynamics in the mixed-time-slicing (mixTS) regime, in which different atoms or particles can be assigned different ring-polymer bead numbers. After presenting a three-particle mixTS Hamiltonian, the paper describes the implementation workflow, including bead-exclusion lists, normal-mode propagation, and input/output conventions. Two applications are reported: q-SPC/Fw bulk water with several bead-number combinations (NH = 32/NO = 4 is shown to reproduce the all-replica NH = NO = 32 radial distribution functions), and aqueous Co2+ modeled as a classical Co3+ ion plus a 1024-bead electron ring polymer. The authors argue that MixPI fills a gap in available PIMD software and that mixTS reduces computational cost while enabling observable-specific convergence studies.","tokens_in":15642,"tokens_out":8774,"duration_ms":87311,"significance":"If the theoretical and numerical issues raised below are resolved, MixPI would be a useful open-source contribution: it targets a real need for PIMD simulations in large condensed-phase systems where only a subset of atoms require explicit quantization. The code is publicly available on GitHub, interfaces with a widely used electronic-structure/MD package, and the water benchmarks demonstrate the intended observable-specific convergence behavior in a clean, reproducible setting. The paper also gives explicit timing breakdowns, which help assess the practical cost of the approach. The central claim of being the first general open-source atomistic mixTS driver is plausible, although the current validation is weakened by the use of a centroid approximation for electrostatic forces without an error estimate, and by the deferral of the full derivation to the User Manual.","major_comments":[{"comment":"As printed, Eq. (13) appears to contain an extra prefactor N_i/N_j multiplying the two-body inter-bead sum. For the N1 = 4, N2 = 2 example described in the text and in Fig. 1d, the correct pair potential is (1/N1) sum over gamma=1..N2 of sum over alpha=1..N1/N2 of V12(q_{1,2(alpha-1)+gamma}, q_{2,alpha}); Eq. (13) as typeset instead yields an additional factor N1/N2 and hence double-counts the interaction. Please correct the equation or clarify the notation, since this is the central Hamiltonian on which the implementation is based.","section":"Section II, Eq. (13)"},{"comment":"The paper states that complete derivations of the all-replica and mixTS Hamiltonians are included in the MixPI User Manual, but the derivation of Eq. (13) is not shown in the manuscript. Because the correctness of the index structure in Eq. (13) is the central theoretical claim, the derivation should be presented in the paper or in an appendix, not only in external documentation. At minimum, show explicitly how the asymmetric Trotter splitting in Eq. (9) leads to the nested bead-index sums in Eq. (13).","section":"Section II, derivation of Eq. (13)"},{"comment":"The smooth-particle mesh Ewald electrostatics is approximated by forces between ring-polymer centroids, with no error estimate. Both the mixTS runs and the all-replica reference runs in Figs. 3-5 use this same centroid approximation, so the reported agreement validates that MixPI reproduces its own approximate reference; it does not validate the exact bead-resolved mixTS Hamiltonian for Coulomb interactions. Since the central software claim is general mixTS-PIMD, please provide a quantitative assessment of the centroid approximation, for example by comparing centroid-PME against bead-resolved PME in a small test system, or by reporting the electrostatic energy/force error in the water benchmark.","section":"Section II, paragraph after Eq. (13), and Section IV"},{"comment":"In the Co3+ + electron simulation, the centroid approximation means that the external electrostatic potential acts uniformly on all electron beads and does not reshape the electron's internal ring-polymer distribution. The delocalized electron cloud shown in the inset of Fig. 5b is therefore essentially the free ring-polymer width, with the centroid confined near the ion by the electrostatic field. Consequently, the observed agreement between the Co3+ + electron and classical Co2+ rdfs does not constitute evidence that a quantized electron with bead-resolved Coulomb interactions reproduces Co2+ solvation structure. Please either replace the centroid electrostatics for this test, or reframe the example as a demonstration of a centroid-driven free ring polymer and state the limitation explicitly.","section":"Section IV.B and Fig. 5"},{"comment":"Equation (14) defines the free ring-polymer Hamiltonian with a negative sign in front of the harmonic spring term, H0_RP = sum p^2/(2m) - (Nm/(2 beta^2 hbar^2))(q_alpha - q_{alpha+1})^2, whereas Eq. (12) has the correct positive sign. If this is a typographical error, it should be corrected; if it reflects the implemented integrator, the Hamiltonian would have an attractive (rather than restoring) spring potential and would not conserve the ring-polymer distribution. Please clarify and fix.","section":"Section III.C, Eq. (14)"}],"minor_comments":[{"comment":"In the sentence describing the quantum-classical case, 'N14' should read 'N1 = 4', and 'will interact will all the beads' should read 'will interact with all the beads'.","section":"Section II (text near Fig. 1)"},{"comment":"There are duplicate 'the' occurrences: 'discuss the the required input' and 'shown in the the MixPI User Manual'; please proofread.","section":"Section III.B"},{"comment":"The table header contains the typo 'In parantheses' (should be 'parentheses'), and the Fig. 3 caption contains 'fully coverged' (should be 'fully converged').","section":"Table I and Fig. 3 caption"},{"comment":"The symbol 'N0' in the caption should be 'NO' to match the notation used elsewhere for the oxygen bead number.","section":"Fig. 4 caption"},{"comment":"The particle count '1298' appears without a clear column header or explanatory text; aligning the table format with Table I would improve readability.","section":"Table II"},{"comment":"The Angstrom symbol and spacing in '14 ˙A box' are typeset inconsistently; please use a uniform notation for units.","section":"Section IV.A"}],"recommendation":"major_revision","confidential_remarks":"The manuscript makes a strong novelty claim ('to the best of our knowledge, a general PIMD software ... is not currently available'). I would ask the authors to double-check this against current capabilities of i-PI and other general PIMD drivers, since any counterexample would affect the framing. The centroid-PME issue and the apparent Eq. (13) prefactor are the main technical risks; both are fixable but need to be addressed before the paper can serve as a reliable reference for the mixTS implementation."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"Short version: this is a genuine software contribution—a CP2K driver for mixed-time slicing PIMD with per-atom bead counts—and the paper is honest about what is new and what is borrowed. The water benchmark is the strongest part. The main weakness is the centroid PME approximation, which makes the benchmarks consistency checks of an approximate force model rather than tests of the exact mixTS Hamiltonian. I'd still send it to review.\n\nThe paper does well: it ships code, documents the workflow, and demonstrates the method on two systems. The mixTS Hamiltonian is from Steele et al. 2011, and the authors don't pretend otherwise. What is new is the general implementation: a standalone driver that uses CP2K as a library, handles arbitrary bead assignments, and uses an exclusion list to manage the mixed-bead interactions. The water rdfs reproduce the all-replica result with NH=32, NO=4, and the comparison across NH=1, NO=4 and NH=4, NO=1 nicely shows observable-specific NQE convergence. That is a useful demonstration. The Co3+ plus 1024-bead electron run is a plausible proof of concept, though the delocalized electron density is the only real output and it's not quantified beyond a snapshot.\n\nThe soft spots are proportionate. The centroid-PME approximation is introduced without an error estimate. Since the same approximate force model is used for both the mixTS and the all-replica reference, the agreement in Figs 3–5 is only an internal consistency check. I agree with that reading. But I don't think it's fatal: the approximation is common, cited, and the paper flags it. The bigger practical issue is reproducibility: no commit hash, no input files in the paper, and no convergence study for the 1024-bead electron. For a software paper, those are the things a referee should push on.\n\nWho will use this: anyone doing PIMD on large systems where most atoms are classical, and people building on CP2K for path integral methods. The paper is worth a serious referee, and it should get one. My recommendation: send to review, with the expectation that the authors add a commit hash, input data, and either an error estimate for the centroid approximation or a clear statement that the benchmarks are intended to test the software implementation, not the exactness of the mixTS Hamiltonian for electrostatics.","headline":"Useful open-source mixTS PIMD driver; the benchmarks check internal consistency under an approximate centroid PME, not the exact Hamiltonian, but the software claim holds and it deserves review.","tokens_in":16164,"tokens_out":2924,"would_cite":true,"duration_ms":28920,"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":"The paper introduces MixPI, which it presents as the first general atomistic implementation of mixed-time-slicing path integral molecular dynamics, with per-atom bead counts validated on water and an aqueous cobalt system.","keywords":["mixed time slicing","path integral molecular dynamics","ring polymer beads","nuclear quantum effects","CP2K","radial distribution functions","quantized electron","aqueous Co2+"],"falsifier":"Re-run the q-SPC/Fw water benchmark at NH=32, NO=4 with exact bead-resolved Ewald electrostatics instead of centroid forces. If the resulting O-O or O-H radial distribution functions differ from the NH=NO=32 reference by more than the differences MixPI currently reports, the benchmark agreement would be an artifact of the centroid approximation rather than evidence that the mixTS Hamiltonian itself is correct.","tokens_in":15222,"feed_emoji":"💧","tokens_out":5717,"duration_ms":57525,"temperature":0.7,"pith_summary":"This paper introduces MixPI, an open-source driver built on CP2K that implements mixed time slicing (mixTS) path integral molecular dynamics, in which different atoms in the same simulation can be represented by different numbers of ring-polymer beads. The authors argue that no general atomistic mixTS software existed before: existing PIMD packages assume a uniform bead count for every particle, which becomes wasteful when only a few atoms carry most nuclear quantum effects. They derive the general mixTS ring-polymer Hamiltonian and validate the code by reproducing converged all-replica radial distribution functions for q-SPC/Fw water with hydrogen atoms at 32 beads and oxygen atoms at 4, and by simulating a Co3+ ion with a 1024-bead ring-polymer electron that recovers the first-solvation-shell structure of Co2+. The point is that mixTS is both a cost-saving route to converged PIMD and a diagnostic for which atoms contribute quantum effects to a given observable.","feed_headline":"MixPI runs path integral MD with a different bead count per atom","feed_subtitle":"A CP2K-based open-source driver matches full-quantum water and solvated-cobalt benchmarks with per-atom bead control.","key_machinery":"The load-bearing object is the mixTS ring-polymer potential Vmix in Eq. 13, which extends the standard all-replica potential by letting particle i have Ni beads, with all bead-number ratios Ni/Nj enforced as integers, and the electrostatic part approximated with centroid forces following ring-polymer contraction practice. It is paired with a split-operator integrator (Eq. 15) that evolves the free ring-polymer Hamiltonian exactly in normal modes and the external potential with half steps, while CP2K supplies the force evaluation. The book-keeping of which bead on one ring polymer interacts with which bead on another is the part that ordinary PIMD software lacks.","core_discovery":"The central claim is that a single code can correctly sample the quantum-classical mixTS Hamiltonian with a unique bead number per particle. Starting from a three-particle asymmetric Trotter factorization of the partition function, the paper derives H_mix = H0 + Vmix (Eqs. 11-13), where the free ring-polymer part H0 is unchanged and the potential Vmix assigns each pair and three-body interaction to the correct bead index with 1/Ni scaling. The implementation generates one system containing all beads, uses exclusion lists to control bead-bead interactions, and integrates the free ring polymers exactly through their normal modes with a split propagator. Benchmarks show that in q-SPC/Fw water the NH=32, NO=4 mixTS run reproduces the NH=NO=32 all-replica O-O, O-H, and H-H radial distribution functions, while the quantum-classical NH=32, NO=1 run misses oxygen-atom quantum effects; for aqueous Co3+ plus a 1024-bead electron, the resulting rdfs match classical Co2+ within line widths.","pith_inferences":["The observed observable-specific convergence suggests a practical adaptive protocol: start a production run with a low-bead calculation, identify which rdfs are converged, and raise beads only for the atom types that change the target observable; MixPI's per-atom bead input makes this a one-line change per atom type.","Because the paper's validation of the mixTS Hamiltonian is entangled with the centroid-Ewald approximation, a clean test of mixTS itself would require a system with no long-range electrostatics, or an exact bead-resolved treatment, before generalizing the 32/4 water result.","If per-atom bead counts are used for real-time RPMD or CMD rather than equilibrium properties, the different ring-polymer spring constants and normal-mode frequencies across atoms could alter how the centroid dynamics is thermostatted and how correlation functions are computed; this is an extension the current equilibrium-focused benchmarks do not address.","When only a few atoms are quantum, mixTS and ring polymer contraction may be complementary: contraction reduces beads for short-range forces, while mixTS reduces beads for atoms, suggesting combined schemes could be applied with the same exclusion-list bookkeeping."],"forward_implications":["Any additive force-field PIMD calculation can now assign different bead counts to different atom types in a single run, so the expensive all-quantum treatment is reserved for atoms that actually need it.","In q-SPC/Fw water, the O-O rdf converges at NH=NO=4 while the O-H and H-H rdfs require more hydrogen beads; this makes observable-specific convergence tests possible within one code.","For a classical Co3+ ion plus a 1024-bead ring-polymer electron in classical water, MixPI reproduces the first-solvation-shell structure of classical Co2+, providing a route to direct simulations of electron transfer in condensed phase.","Because mixTS reduces the number of force evaluations mainly when potentials are additive, its largest speedups come with force-field models; many-body potentials such as DFT will see cost similar to all-replica PIMD.","The current implementation requires bead counts that are powers of two, which constrains the achievable bead-ratio combinations."],"supporting_citations":[{"why":"Defines mixed time slicing and supplies the asymmetric Trotter factorization that the MixPI Hamiltonian is built from.","marker":"[62]"},{"why":"Provides the PILE thermostat and the symplectic split-operator integration scheme used to evolve the free ring polymers.","marker":"[38]"},{"why":"Establishes ring polymer contraction and the centroid force approximation that MixPI uses for Ewald electrostatics.","marker":"[48]"},{"why":"Supplies the q-SPC/Fw water model and the all-replica reference behavior for the bulk water benchmarks.","marker":"[79]"},{"why":"The prior in-house quantum-classical PIMD study of cobalt self-exchange whose aqueous Co2+ system MixPI reproduces.","marker":"[4]"},{"why":"CP2K provides the optimized force and energy evaluation routines that MixPI drives as an external library.","marker":"[45]"}],"fun_headline_variants":["MixPI: per-atom beads for path integral MD","Path integral MD with a different bead count per atom","MixPI lets each atom pick its own bead number","Per-atom beads: mixed-time slicing in CP2K","Quantum dynamics with adaptive bead counts"],"cache_read_input_tokens":3200,"weakest_assumption_plain":"The benchmarks assume that replacing the smooth-particle-mesh-Ewald electrostatic forces between individual beads by forces between ring-polymer centroids is accurate enough for these systems; the paper gives no error estimate for the 32/4 water run or the cobalt-electron system, so the agreement with all-replica results is only as strong as that approximation.","fun_headline_variants_meta":{"raw":{"variants":["MixPI: per-atom beads for path integral MD","Path integral MD with a different bead count per atom","MixPI lets each atom pick its own bead number","Per-atom beads: mixed-time slicing in CP2K","Quantum dynamics with adaptive bead counts"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000415,"raw_usage":{"total_tokens":2203,"prompt_tokens":1063,"completion_tokens":1140,"prompt_tokens_details":{"cached_tokens":384},"prompt_cache_hit_tokens":384,"prompt_cache_miss_tokens":679,"completion_tokens_details":{"reasoning_tokens":1065}},"tokens_in":679,"tokens_out":1140,"duration_ms":10434,"temperature":1.0,"reasoning_tokens":1065,"cache_read_input_tokens":384,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-12T18:02:29.980046+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Re-run the q-SPC/Fw water benchmark at NH=32, NO=4 with exact bead-resolved Ewald electrostatics instead of centroid forces. If the resulting O-O or O-H radial distribution functions differ from the NH=NO=32 reference by more than the differences MixPI currently reports, the benchmark agreement would be an artifact of the centroid approximation rather than evidence that the mixTS Hamiltonian itself is correct.","supporting_citations":[{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"Defines mixed time slicing and supplies the asymmetric Trotter factorization that the MixPI Hamiltonian is built from."},{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"The prior in-house quantum-classical PIMD study of cobalt self-exchange whose aqueous Co2+ system MixPI reproduces."}],"review_version":1}