{"id":"41166fe2-094d-46f9-9537-759e6eccac7a","arxiv_id":"2505.20418","paper_version":2,"verdict":"CONDITIONAL","confidence":"HIGH","novelty_score":5.0,"correctness_risk":"medium","formal_verification":"none","parameter_count":0,"one_line_summary":"Treating the grid potential as a zero-mass variable and solving a constraint equation instead of the Poisson equation directly converges in fewer multigrid iterations and scales linearly.","lead":"Poisson MaZe computes electrostatic forces in molecular dynamics by treating the grid values of the electric potential as zero-mass auxiliary variables and enforcing the Poisson equation as a dynamic constraint. It converges in fewer multigrid iterations than a direct Poisson solve, with the same per-iteration cost, and reproduces molten NaCl properties with linear scaling.","discovery_kind":"new_method","skeptic_critique":{"model":"deepseek-v4-flash","headline":"Speedup claim for Poisson MaZe rests on inconsistent residual scaling: the direct Poisson residual is (4π/h)≈73× the MaZe residual for the same potential error, so equal absolute tolerances compare different accuracies.","rationale":"The reader's weakest assumption concerned the unbiasedness of the zero-inertia limit and the effect of finite MG tolerance on the dynamics. That is a legitimate statistical-mechanics concern, but the more decisive issue is numerical: the benchmark that supports the central novelty compares residuals of two different equations with a fixed absolute tolerance, while the residuals are related by a factor of 4π/h ≈ 73. This makes the comparison unfair and likely explains the reported iteration savings. The method may still be a correct reformulation of the Poisson solve, and the structural/dynamical validation on molten NaCl is encouraging, but the demonstrated 'fewer iterations' and 'overall reduction of computational time' are not trustworthy until the comparison is normalized. A fair comparison could rescue the claim if Poisson MaZe remains faster at equal Poisson-residual accuracy; hence the verdict should be CONDITIONAL rather than outright REJECT, with the condition being a corrected performance benchmark. The reader's focus on the MaZe theorem is not wrong, but it does not identify the load-bearing weakness of the performance argument.","tokens_in":19623,"tokens_out":12414,"duration_ms":127030,"concrete_test":"Re-run the convergence benchmark of Section III A 1 using a common physical stopping criterion: stop both solvers when the infinity-norm residual of the discretized Poisson equation for the actual potential, ||M φ + 4π q/h||_∞, is below 10^-7. For Poisson MaZe, compute φ after each iteration as φ = φ_p − (4π/h) y and evaluate this true Poisson residual, not the Eq 23 residual, for the same set of 10^3 configurations. If the iteration counts for time-propagated direct MG and Poisson MaZe become comparable, or MaZe is not faster, the reported speedup is an artifact of absolute-residual scaling.","verdict_should_be":"CONDITIONAL","load_bearing_attack":"The central performance claim—that solving the Poisson MaZe linear system (Eq 15) requires fewer multigrid iterations than solving Eq 4 directly, at equal per-iteration cost—is supported by a comparison with inconsistent stopping criteria. The direct Poisson residual (Eq 22) is r_direct = M φ^ν + 4π q/h, while the MaZe constraint residual (Eq 23) is r_MaZe = M y^ν − σ_p, with σ_p = h/(4π) M φ_p + q. For any candidate potential φ = φ_p − (4π/h) y, these residuals are related exactly by r_direct = (4π/h) r_MaZe. With h = 0.172 Å, the factor is 4π/h ≈ 73. Stopping both solvers at the same absolute residual 10^-7 therefore allows Poisson MaZe to terminate with a Poisson-equation residual about 73 times larger than the direct solver. Since MG iteration count to an absolute tolerance depends on the initial residual and the target, this scaling mismatch can account for most or all of the reported iteration reduction (4 vs 8 cycles for time-propagated direct MG; 4 vs 28 for static). The claim in Section III A 1 that the two methods yield potentials agreeing to ~10^-10–10^-14 is also inconsistent with this scaling: an effective residual difference of several×10^-6 should produce field differences far larger than 10^-10. Appendix B's attribution of the speedup to 'improved conditioning' is also incorrect because the matrix M is identical; the apparent acceleration is a consequence of the absolute-residual criterion and the small RHS σ_p, not of the linear system's conditioning. Thus the paper's headline advantage—'converges with fewer cycles' and 'overall reduction of computational time'—is not established under equivalent accuracy and may be an artifact of inconsistent residual normalization.","agreement_with_reader":"partial"},"referee_report":{"model":"deepseek-v4-flash","summary":"The manuscript introduces Poisson MaZe, an extended-Lagrangian method for computing electrostatic interactions in molecular dynamics. The values of the discretized electrostatic potential are treated as zero-inertia auxiliary degrees of freedom, and the discretized Poisson equation is imposed as a constraint. The constraint is solved via SHAKE, leading to a linear system (Eq. 15) that differs from the direct Poisson system (Eq. 4). The authors validate the method on molten NaCl, reporting agreement with reference structural and transport properties, and demonstrate linear scaling with a multigrid solver. They also report that Poisson MaZe requires fewer solver iterations than direct multigrid solutions of the Poisson equation, and they discuss energy, momentum, time-reversibility, and stationarity properties. Appendices provide a toy-model validation and a conjugate-gradient variant.","tokens_in":19864,"tokens_out":19562,"duration_ms":193078,"significance":"If the performance claims are substantiated, Poisson MaZe would be a useful real-space alternative to FFT-based methods for large-scale simulations. The method is parameter-free, built on a clearly stated Lagrangian, and validated against external benchmarks (CHARMM SOR and published molten-salt data). The convergence and stationarity tests are valuable, and the linear-scaling property inherited from the multigrid solver is a strength. However, the central acceleration claim is not currently established because the residual comparison in Section III A 1 is not accuracy-equivalent, and the conditioning explanation in Appendix B is incorrect.","major_comments":[{"comment":"The comparison of iteration counts is not accuracy-equivalent. For any candidate potential φ = φ_p - (4π/h)y, the direct Poisson residual and the MaZe residual are related by r_direct = Mφ + (4π/h)q = -(4π/h)(My - σ_p) = -(4π/h) r_MaZe. With h = 0.172 Å this is a factor of about 73. Stopping both solvers at the same absolute tolerance 10^-7 therefore allows Poisson MaZe to terminate with a Poisson-equation residual about 73 times larger than the direct solver. The reported differences in iteration counts (4 vs 8/28 in Fig. 4) and the speedups (3.36±0.21 and 1.88±0.15 in Section III A 2) are consequently not comparisons at equal accuracy. The accompanying statement that the final potentials agree to ~10^-10–10^-14 is also inconsistent with this scaling: with tol=10^-7, the implied potential error for the MaZe solution is O(||M^{-1}|| × 7×10^-6), which is orders of magnitude larger than 10^-10. Please repeat the benchmark with equal Poisson-equation residual (or equal potential error) and report the resulting iteration counts and timings.","section":"III A 1, Eqs. (22)-(23)"},{"comment":"The claim that the MaZe acceleration \"stems from an improved conditioning of the problem\" is incorrect. The linear system solved by Poisson MaZe, Eq. (15), uses exactly the same matrix M as the direct Poisson system Eq. (4), so the condition number is identical. The faster convergence observed in Fig. 15 is plausibly due to the smaller right-hand side σ_p and the warm start, combined with the absolute-residual criterion discussed above; the conditioning explanation should be removed or revised.","section":"Appendix B"},{"comment":"The force used in production, Eq. (18), is not the derivative of the extended-Lagrangian potential V(φ;{rβ}) appearing in Eq. (11). The derivative of V with respect to r_α at fixed φ is -Σ_n φ_n ∂q_n/∂r_α, whose discretized form is the alternative force Eq. (20); Eq. (18) differs by using the central difference of the potential rather than the gradient of the B-spline weights. Since the exact-sampling statement in the Conclusions cites Ref. 31 for the equations of motion (11), the authors should clarify whether the unbiasedness theorem applies to the implemented force Eq. (18) or only to the Lagrangian force Eq. (20). If it applies only to Eq. (20), the exact-sampling claim for the actual algorithm needs qualification.","section":"II A and Conclusions"}],"minor_comments":[{"comment":"The factor 1/2 in front of Δt^2 Mη is inconsistent with the equation of motion φ¨ = -Mη in Eq. (11), for which the Verlet update would contain -Δt^2 Mη. Since the algorithm solves for y directly, this does not change the final update rule, but the derivation should be corrected or the definition of η clarified.","section":"II, Eq. (12)"},{"comment":"The references to \"arrows in Figure 4\" and \"the inset of Figure 4\" should refer to Figure 15, which is the PLCG convergence figure.","section":"Appendix B"},{"comment":"The diffusion coefficients and g(r) peak parameters are quoted without statistical uncertainties; given the 25 ps trajectory, reporting error bars would allow a quantitative assessment of the agreement with Refs. 41 and 42.","section":"III and Table I/Figure 3"},{"comment":"The top-axis labels appear to have lost their formatting (e.g., \"1043\" rather than a readable grid-count notation); please reformat the axis labels.","section":"Figure 5"}],"recommendation":"major_revision","confidential_remarks":"The manuscript contains a promising method and a careful validation of structural and transport properties, but the performance claims need to be re-benchmarked at equal accuracy before publication. The force/Lagrangian mismatch is also worth addressing. I lean toward major revision rather than rejection because the issues are fixable within the manuscript's scope."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"Dear colleague,\n\nThe thing you should know: the paper's central speedup claim—that solving Eq 15 with multigrid needs fewer cycles than solving Eq 4—doesn't survive a close look at the residual definitions. The two residuals are proportional, not comparable at equal absolute tolerance.\n\nWhat is genuinely new: the authors apply the MaZe framework to Poisson electrostatics, deriving an extended Lagrangian where grid potentials are zero-inertia variables and the Poisson equation is a constraint. The resulting linear system (Eq 15) is different from the standard Poisson solve, and the SHAKE/MG coupling is a sensible implementation. The derivation is careful, and the molten NaCl validation (g(r), diffusion) matches the benchmarks. The work-based energy estimator in Eq 19 is a useful practical addition.\n\nThe soft spots, in order of importance. First, the performance comparison is unfair. For a given potential error, the direct residual (Eq 22) and the MaZe residual (Eq 23) are related exactly by r_direct = (4π/h) r_MaZe. With h = 0.172 Å, that's a factor of 73. Stopping both solvers at 10^-7 means the MaZe solution is allowed to have a Poisson residual 73 times larger. That alone can explain the 4-vs-8 iteration difference. The authors' claim that the final potentials agree to 10^-10–10^-14 is hard to square with this factor. They need to compare at equivalent physical accuracy—say, by setting the MaZe tolerance to 10^-7/(4π/h) or by reporting the direct residual achieved by each method. Second, Appendix B attributes the acceleration to 'improved conditioning,' which is wrong: the matrix M is identical in Eq 4 and Eq 15. The acceleration, if real, comes from the different right-hand side and initial guesses, not from conditioning. Third, the benchmark section lacks error bars on the physical observables, and there is no released code or data, so independent reproduction of the key number is not possible.\n\nNone of this kills the method. The MaZe formulation is coherent and the validation is solid. But the headline numerical advantage is not established, and the paper as written overstates it. This deserves a serious referee, but a conditional one: the performance section needs a corrected comparison and the conditioning statement needs to be fixed. If the speedup survives an apples-to-apples test, this is a useful contribution to the real-space electrostatics toolbox.\n\nBest.","headline":"The MaZe-for-Poisson formulation is a real idea, but the reported multigrid speedup looks like a residual-normalization artifact, so the paper needs a corrected performance comparison before the main claim can be trusted.","tokens_in":20507,"tokens_out":4416,"would_cite":false,"duration_ms":42642,"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":"Poisson MaZe computes electrostatic forces by evolving the grid potential as zero-inertia auxiliary variables constrained by the Poisson equation, reducing the multigrid iterations needed and giving linear-scaling real-space electrostatics.","keywords":["Poisson equation","electrostatic interactions","molecular dynamics","mass-zero constrained dynamics","multigrid","real-space Poisson solver","extended Lagrangian","molten NaCl"],"falsifier":"Run a moderate-size electrolyte or molten-salt simulation with multigrid tolerances of $10^{-4}$ and $10^{-10}$ and compare the resulting radial distribution functions and diffusion coefficients with a high-accuracy particle-mesh Ewald reference; any statistically significant difference between the two tolerance runs would show that the finite solver tolerance contaminates the sampled ensemble.","tokens_in":19364,"feed_emoji":"⚡","tokens_out":10335,"duration_ms":104022,"temperature":0.7,"pith_summary":"Poisson MaZe is a real-space method for computing electrostatic forces in molecular dynamics that avoids reciprocal-space Fast Fourier transforms. The paper argues that the value of the electrostatic potential at every grid point can be treated as a zero-mass auxiliary coordinate, with the discretized Poisson equation imposed as a dynamical constraint; the constraint update requires solving a linear system whose right-hand side is the mismatch between a Verlet-predicted potential and the ionic charge density. In a molten NaCl simulation the method reproduces the reference radial distribution functions and diffusion coefficients, and the multigrid solution of the constraint system converges in fewer cycles than solving the Poisson equation directly. The authors test the method as a complete electrostatic solver and demonstrate linear scaling with system size.","feed_headline":"Cut multigrid cycles by treating Poisson as a constraint","feed_subtitle":"A zero-inertia potential field solves electrostatics in real space with fewer solver iterations and linear scaling.","key_machinery":"The carrying object is the extended MaZe Lagrangian, in which the $N$ grid values of the electrostatic potential $\\phi_n$ enter as auxiliary variables with zero inertia and the constraint $\\sigma=\\frac{h}{4\\pi}M\\phi+q=0$ is enforced by Lagrange multipliers. A Verlet step proposes a trial potential $\\phi_p^{k+1}=2\\phi^k-\\phi^{k-1}$; requiring the corrected potential to satisfy the constraint turns the SHAKE step into the linear system $M y=\\sigma_p^{k+1}$ (Equation 15), where $y$ is the rescaled multiplier vector and $\\sigma_p$ is the constraint evaluated at the trial potential. This system shares the Laplacian matrix $M$ with the standard Poisson equation but has a different unknown and a right-hand side that measures how far the trial potential is from satisfying the Poisson equation, and it is solved iteratively with multigrid, initialized with the previous step's solution, which is the source of the faster convergence.","core_discovery":"The central claim is that embedding the Poisson problem in a mass-zero constrained extended dynamics yields a new linear system, $M y=\\sigma_p$, whose multigrid solution reaches the electrostatic potential faster than solving the standard discretized Poisson equation $M\\phi=-4\\pi q/h$. At a tolerance of $10^{-7}$, Poisson MaZe converges in about four multigrid cycles, whereas a time-propagated direct multigrid solve needs roughly twice as many and a static-start direct solve roughly seven times as many, at essentially the same cost per cycle. The resulting forces, computed with a central-difference gradient of the grid potential and B-spline charge assignment, conserve total momentum and, when combined with a work-based energy estimator, keep the total energy stable at a relative fluctuation of $3.2\\times 10^{-7}$. The validation on 250-ion molten NaCl matches the reference structural and transport observables, and the measured scaling of total CPU time with particle number is linear for all methods, with Poisson MaZe lowest on every tested size.","pith_inferences":["Beyond the paper, the same constraint-embedding trick should transfer to other elliptic equations solved inside a molecular dynamics loop, such as generalized Poisson or linearized Poisson-Boltzmann problems, whenever the previous time step supplies a good initial guess for the iterative solve.","If the per-cycle cost remains flat at high core counts, replacing the global FFT communication of particle-mesh Ewald with local multigrid smoothing could make real-space electrostatics competitive at very large processor counts; the paper only measures scaling on 8 processors, so this is an extrapolation rather than a demonstrated result.","Because the trial potential is a second-order time prediction, the convergence gain should grow for slowly moving ions and small time steps; the paper's solid-versus-molten comparison (slowdown factor about 1.8 versus about 2.9 to 3.7 on melting) is consistent with that expectation."],"forward_implications":["If the claim is right, long-range electrostatics can be computed entirely in real space with linear scaling, removing the global communication pattern of FFT-based particle-mesh methods.","At a multigrid tolerance of $10^{-7}$, Poisson MaZe converges in about four cycles per step, versus roughly seven cycles for a time-propagated direct solve and roughly 28 for a static-start direct solve, cutting total CPU time by factors of about 1.9 and 3.4 in the tested systems.","The molten NaCl runs reproduce the reference first-peak positions and heights for all three radial distribution functions and give diffusion coefficients $D_{\\text{Na}}=0.16\\times 10^{-3}$ cm$^2$/s and $D_{\\text{Cl}}=0.14\\times 10^{-3}$ cm$^2$/s, so the method yields correct structural and transport properties in a realistic liquid.","The same convergence advantage appears with a preconditioned conjugate-gradient solver, indicating that the benefit comes from the reformulated linear system rather than from a property specific to multigrid.","The work-based electrostatic energy estimator provides a stable conserved quantity for monitoring simulations, with relative total-energy fluctuations of $3.2\\times 10^{-7}$, while the naive discretized energy oscillates badly due to grid self-energy terms."],"supporting_citations":[{"why":"It proves that the zero-inertia limit of MaZe extended dynamics samples the exact probability density for the physical degrees of freedom, which is the premise that makes the ionic dynamics unbiased.","marker":"[31]"},{"why":"It supplies the SHAKE algorithm used to solve the constraint and to initialize the iterative cycle with the previous solution.","marker":"[26]"},{"why":"It supplies the multigrid operators (red-black Gauss-Seidel smoothing, 27-point restriction, trilinear prolongation) and the real-space Poisson solver baseline used for the comparison.","marker":"[23]"},{"why":"It supplies the central-difference force estimator and B-spline charge assignment that make the forces momentum-conserving.","marker":"[21]"},{"why":"It supplies the grid discretization of the Poisson equation and the comparison Poisson solver used for the two-charge validation test.","marker":"[33]"},{"why":"It provides the reference molten NaCl pair distribution functions and diffusion coefficients against which the method is validated.","marker":"[41]"}],"fun_headline_variants":["Mass-zero constraint cuts multigrid cycles for Poisson","Zero-inertia auxiliary field solves electrostatics in real space","Fewer solver cycles with mass-zero Poisson dynamics","Linear scaling electrostatics via zero-mass constraint","Real-space Poisson via zero inertia: fewer iterations"],"cache_read_input_tokens":3200,"weakest_assumption_plain":"The result depends on the zero-inertia limit of the constrained extended dynamics being exactly unbiased for the physical degrees of freedom, and on the finite multigrid tolerance not biasing the sampled configurations; if either fails, the dynamics will depart from the intended Boltzmann statistics even though the potential is correct.","fun_headline_variants_meta":{"raw":{"variants":["Mass-zero constraint cuts multigrid cycles for Poisson","Zero-inertia auxiliary field solves electrostatics in real space","Fewer solver cycles with mass-zero Poisson dynamics","Linear scaling electrostatics via zero-mass constraint","Real-space Poisson via zero inertia: fewer iterations"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.001117,"raw_usage":{"total_tokens":4683,"prompt_tokens":1009,"completion_tokens":3674,"prompt_tokens_details":{"cached_tokens":384},"prompt_cache_hit_tokens":384,"prompt_cache_miss_tokens":625,"completion_tokens_details":{"reasoning_tokens":3598}},"tokens_in":625,"tokens_out":3674,"duration_ms":28319,"temperature":1.0,"reasoning_tokens":3598,"cache_read_input_tokens":384,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-07T13:55:37.379854+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Run a moderate-size electrolyte or molten-salt simulation with multigrid tolerances of $10^{-4}$ and $10^{-10}$ and compare the resulting radial distribution functions and diffusion coefficients with a high-accuracy particle-mesh Ewald reference; any statistically significant difference between the two tolerance runs would show that the finite solver tolerance contaminates the sampled ensemble.","supporting_citations":[{"cited_title":"Coretti , author T","cited_arxiv_id":null,"evidence_quote":"It proves that the zero-inertia limit of MaZe extended dynamics samples the exact probability density for the physical degrees of freedom, which is the premise that makes the ionic dynamics unbiased."},{"cited_title":"\\ Ryckaert , author A","cited_arxiv_id":null,"evidence_quote":"It supplies the SHAKE algorithm used to solve the constraint and to initialize the iterative cycle with the previous solution."},{"cited_title":"Sagui \\ and\\ author T","cited_arxiv_id":null,"evidence_quote":"It supplies the multigrid operators (red-black Gauss-Seidel smoothing, 27-point restriction, trilinear prolongation) and the real-space Poisson solver baseline used for the comparison."},{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"It supplies the central-difference force estimator and B-spline charge assignment that make the forces momentum-conserving."},{"cited_title":"\\ Ryckaert , author G","cited_arxiv_id":null,"evidence_quote":"It supplies the grid discretization of the Poisson equation and the comparison Poisson solver used for the two-charge validation test."},{"cited_title":"Klapper , author R","cited_arxiv_id":null,"evidence_quote":"It provides the reference molten NaCl pair distribution functions and diffusion coefficients against which the method is validated."}],"review_version":1}