{"id":"314fe3ac-706d-4c12-8068-8ce422e0f2dd","arxiv_id":"1908.07341","paper_version":1,"verdict":"CONDITIONAL","confidence":"MODERATE","novelty_score":4.0,"correctness_risk":"medium","formal_verification":"none","parameter_count":0,"one_line_summary":"Numerical tests confirm the shaken dynamics reproduces the predicted Ising phase transition curve, while mixing-time and GPU benchmark results support its use as a fast parallel sampler.","lead":"This paper numerically tests the 'shaken dynamics', a parallel Monte Carlo method for Ising spin models, and verifies its predicted phase transition curve across a range of lattice geometries. It also measures how fast the method mixes and shows a GPU implementation that runs about 500 times faster than a single CPU core.","discovery_kind":"new_application","skeptic_critique":{"model":"deepseek-v4-flash","headline":"Equilibration of the 200x200 critical-curve runs is unchecked; 300k warm-up may be below the mixing time near Jc(q), so the variance ridge in Figs. 6-7 may not be an equilibrium signal.","rationale":"The paper's central claim has two parts: a phase transition along Jc(q), and closeness of the shaken equilibrium measure to the square-lattice Gibbs measure for large q. The second part is supported by exact Propp-Wilson sampling, so the warm-up issue does not affect it. The first part, however, is validated only by the time-averaged magnetization variance in Section 3.1, which assumes the chain has reached equilibrium after 300,000 steps. Near criticality, mixing times are long, and the paper's own coalescence-time data for L=32 show values comparable to or exceeding the warm-up length at several points. Without an equilibration diagnostic or finite-size scaling, the observed variance ridge cannot be unambiguously attributed to the equilibrium phase transition. This is the most load-bearing concern because it directly affects whether Figs. 6-7 test eq. (11). A secondary, smaller issue is that the displayed formula for the marginal transition probability P_box in Section 2 appears to have Z_tau in both denominators; the correct marginal from eqs. (15)-(16) should have Z_sigma in the first denominator. Algorithm 3 appears to implement the correct update, so this is likely a typo rather than a substantive error, but it should be corrected. Neither concern changes the overall CONDITIONAL verdict: the qualitative conclusions remain plausible, while the quantitative support is not yet fully solid.","tokens_in":12494,"tokens_out":16897,"duration_ms":172855,"concrete_test":"Rerun the Section 3.1 protocol for a fixed set of q values (e.g., q=0.2, 0.65, 1.0, 2.5) with warm-up lengths of 3e5, 1e6, and 3e6 steps, and with both all-plus and all-minus initial conditions, on L=200 and also L=64; record the J-location of the maximum magnetization variance for each q. If the peak position changes by more than the grid spacing (0.025) or depends on the initial condition, the 300,000-step warm-up is insufficient and the Fig. 6-7 validation of eq. (11) is unsupported. If the peak is stable across warm-up lengths and initial conditions, the concern is resolved.","verdict_should_be":"UNCHANGED","load_bearing_attack":"Section 3.1 locates the phase transition by the variance of the magnetization of the shaken dynamics on a single 200x200 torus after 300,000 warm-up steps, with no mixing diagnostic. This is the only evidence in the paper that the numerically observed ridge validates eq. (11). But near Jc(q) the chain is expected to be slow; the paper's own coalescence measurements in Fig. 9 for L=32 reach 10^5 to 10^6 steps at criticality, and mixing on L=200 should be slower. If 300,000 steps is below the mixing time for a substantial part of the (q,J) grid, the time-averaged magnetization and its variance are not equilibrium quantities, and the ridge in Fig. 7 could be displaced or broadened by relaxation from the all-minus initial condition. The text states that the choice 'turned out to be good enough' but gives no quantitative check, such as agreement between runs from all-plus and all-minus starts, or stability of the variance peak as the warm-up length increases. Because the claimed validation of the critical curve rests on this measurement, the numerical support for the central claim is weaker than the figures suggest. This concern does not affect the exact-sampling comparisons in Figs. 10-13, which use Propp-Wilson and are not subject to the warm-up assumption.","agreement_with_reader":"partial"},"referee_report":{"model":"deepseek-v4-flash","summary":"The manuscript presents numerical simulations of the 'shaken dynamics', a parallel probabilistic cellular automaton for the 2D Ising model on a bipartite lattice with parameters J and q. The authors estimate the critical curve Jc(q) in the (q,J) plane via magnetization-variance measurements, study coalescence times and mixing behavior, compute spin-spin correlations, compare the shaken dynamics' equilibrium measure with the Gibbs measure of the square-lattice Ising model for large q, and report CPU/GPU implementations with a benchmark. The paper is positioned as numerical support for the companion theoretical papers [1] and [2], whose predictions it aims to verify.","tokens_in":12767,"tokens_out":6425,"duration_ms":59850,"significance":"If the numerical evidence is quantitatively sound, the paper would provide useful independent validation of the analytical critical curve and of the claim that the shaken dynamics approximates the square-lattice Ising Gibbs measure for q≥2.5. Its strengths include the use of Propp-Wilson exact sampling for the equilibrium comparisons in Figs. 10-13, which avoids warm-up bias in those comparisons, and the concrete algorithmic description with a GPU benchmark. The central limitation is that the critical-curve identification rests on an unchecked equilibration assumption for a single lattice size, and several key quantities are reported without error bars or quantitative comparisons to the analytic predictions.","major_comments":[{"comment":"The only numerical evidence for the critical curve is the variance ridge of the magnetization on a single 200×200 torus after a fixed 300,000-step warm-up from the all-minus state. This does not establish equilibration: the paper's own coalescence measurements in Fig. 9 for L=32 reach 10^5-10^6 steps near criticality, and mixing on L=200 should be slower. Provide a mixing diagnostic (e.g., agreement between all-plus and all-minus starts, or a plot of the variance peak versus warm-up length) and a quantitative measure of the distance between the observed ridge and Jc(q). Without these, the validation of Eq. (12) is weaker than Fig. 7 suggests.","section":"§3.1, Eq. (11)/(12), Fig. 7"},{"comment":"The threshold variance ≥0.03 used to center the bars in Fig. 7 is arbitrary, and no justification is given for why this threshold identifies the transition rather than merely a fixed fluctuation level. Since the height and width of the variance peak depend on lattice size, run length, and correlation time, the threshold should be justified or replaced by a more robust estimator (for example, a Binder-cumulant or finite-size-scaling analysis).","section":"§3.1, Fig. 7"},{"comment":"Coalescence times are central to the mixing-time comparisons, yet they are reported as sample averages with no standard errors, no number of independent repetitions, and no indication of how many chains were coupled. This makes it impossible to judge whether the apparent speedups of the shaken or alternate dynamics over single-spin-flip dynamics are statistically meaningful, especially near criticality where coalescence times fluctuate strongly.","section":"§3.2, Figs. 8-9 and 14"},{"comment":"The statement that 'for q≥2.5 the approximation provided by the shaken dynamics is quite good' is not backed by a quantitative distance between the shaken-dynamics observables and the Gibbs-reference values. Several panels, e.g. the energy standard deviation panels in Fig. 13, show discrepancies that appear larger than the between-algorithm scatter, and no error bars are given. Report a quantitative measure (e.g., the maximum observable bias or an estimated total-variation distance) over the plotted J values.","section":"§3.2, Figs. 10-13"},{"comment":"The spin-spin correlation table is presented without error bars or a statement of the number of independent samples used, so the directional asymmetry conclusion and the claim that correlations decay rapidly below Jc rest on uncharacterized single-run estimates. For example, the q=0.05 supercritical row shows the SW-NE correlation at l=16 (0.739) exceeding the NW-SE value at the same distance (0.726), which is not discussed and may indicate large statistical uncertainty.","section":"§3.3, Table 1"}],"minor_comments":[{"comment":"The paragraph beginning 'In this framework, a new PCA parameterized by J and q...' is duplicated verbatim within the Introduction.","section":"Section 1"},{"comment":"The sentence 'The critical value of βc separates the ordered phase where all the spin have the same probability... from the ordered phase where the measure is polarized' should read 'disordered phase' in the first instance.","section":"Section 2, after Eq. (6)"},{"comment":"The text says the simulations are run 'for (J,q)∈{(0,2)×(0,2)} on a 80×80 grid', while Section 4 says 80 couples of (q,J) values were simulated; please clarify whether the parameter grid has 80 points or 6400 points.","section":"Section 3.1"},{"comment":"Figure 8 lacks a color scale label; the reader has to infer that the color encodes the logarithm of the average coalescence time, and no units or error information are provided.","section":"Section 3.2, Fig. 8"},{"comment":"The phrase 'we can not go beyond105 for this GPU' appears to have a missing exponent or line break; it should read something like '10^5'.","section":"Section 4.0.1"},{"comment":"The benchmark plot shows single timings with no error bars or repeated measurements; since the text claims a speed-up factor of approximately 500, a statement of measurement variability would strengthen the claim.","section":"Section 4.0.2, Fig. 20"}],"recommendation":"major_revision","confidential_remarks":"The paper is a numerical companion to two theoretical papers by the same group and would be a reasonable fit for physics.comp-ph once the numerical evidence is made quantitative. The main revision burden is to replace the unchecked 300,000-step warm-up assumption in §3.1 with explicit mixing diagnostics and to add error bars and quantitative comparisons throughout. I do not see evidence of circularity: the equilibrium comparisons use Propp-Wilson sampling, and the critical-curve measurement is an independent numerical probe of Eq. (12), provided the equilibration concern is addressed."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"Colleague,\n\nWhat you should know first: this is a numerical support paper for the shaken dynamics PCA introduced in the companion papers [1,2]. The genuinely new content is the numerical verification of the critical curve, the mixing-time comparisons, the correlation table, and the GPU benchmark. None of that is revolutionary, but it is useful and mostly done carefully.\n\nThe strong parts are the Propp-Wilson exact-sampling comparisons in Figs. 10-13. Those are not subject to a warm-up assumption, and they give real evidence that for q >= 2.5 the shaken dynamics equilibrium is close to the square-lattice Gibbs measure. The coalescence-time comparison with heat bath, renormalized by volume, is a fair and informative way to show that the parallel dynamics is not cheating. The correlation table qualitatively supports the anisotropy theorem from [1]. Credit where due: the numerical experiments are independent checks of the companion theorems, not fits to them.\n\nNow the soft spots, in proportion. The stress-test note is right. Section 3.1 locates the phase transition by looking at the magnetization variance on a single 200x200 torus after a fixed 300k warm-up, with no mixing diagnostic. The paper's own Fig. 8 shows coalescence times on L=32 reaching 10^5-10^6 steps near criticality; L=200 should be slower. The text admits the warm-up 'turned out to be good enough' without any quantitative check, such as agreement from all-plus vs all-minus starts or a variance-peak stability scan. So the variance ridge in Figs. 6-7 may not be an equilibrium signal, and the numerical support for eq. (11) is weaker than the figures suggest. This concern does not affect the exact-sampling sections, so the central q>=2.5 claim survives.\n\nMinor but real: no error bars on coalescence times, correlations, or the benchmark. The GPU speedup factor of 500 comes from one measurement, no runs shown for variance. Code is not provided, despite the promise of a Julia library. The writing has some rough edges (e.g., duplicated paragraph in the intro, 'stime stesps' typo). None of that changes the scientific content.\n\nBottom line: this paper deserves a serious referee, not a desk reject. The numerical work is worth publishing after the authors add mixing diagnostics for the critical-curve runs, report error bars, and ideally release code. If they address the equilibration issue, the paper is a solid contribution to computational statistical mechanics.\n\nFor the reading group I'd say maybe; it is not urgent but it is relevant if anyone cares about parallel samplers for Ising-type models. I would not cite it in my own next work, but I would point to [1,2] instead for the dynamics and to this paper as a reference for implementation details.","headline":"Useful numerical companion to the shaken-dynamics papers, with a real soft spot in the equilibration check for the critical-curve estimate; the exact-sampling parts hold up.","tokens_in":13254,"tokens_out":2121,"would_cite":false,"duration_ms":22488,"reading_group":"maybe","serious_thinker":"yes","would_accept_peer_review":true},"rs_alignment":null,"lean_confirmation":null,"pith_extraction":{"msc":["82B20","82B80","60J22"],"pacs":["05.10.-a","05.50.+q"],"model":"deepseek-v4-flash","headline":"The shaken dynamics, a fully parallel probabilistic cellular automaton, reproduces the Ising critical curve $J_c(q)=\\tanh^{-1}(\\sqrt{\\tanh^2 q+1}-\\tanh q)$ and, for $q\\ge2.5$, samples close to the square-lattice Gibbs measure.","keywords":["shaken dynamics","probabilistic cellular automaton","Ising model","parallel Markov chain Monte Carlo","critical curve","mixing time","spin-spin correlations","GPU simulation"],"falsifier":"Directly estimate the total-variation distance $\\|\\pi_s - \\pi_G\\|_{TV}$ from Propp-Wilson samples for $q=2.5$ on increasing lattice sizes; if the distance does not tend to zero as $|\\Lambda|\\to\\infty$, or if the magnetization variance peaks away from $J_c(q)$ on lattices larger than $200\\times200$, the central claims would be refuted.","tokens_in":12325,"feed_emoji":"🧲","tokens_out":8581,"duration_ms":66800,"temperature":0.7,"pith_summary":"The paper numerically investigates the \"shaken dynamics\", a parallel Markov chain for Ising spin systems whose transition rule depends on two parameters, $q$ and $J$, that control the effective lattice geometry. It aims to establish that the same dynamics, with no change of algorithm, can simulate Ising models ranging from one-dimensional chains (small $q$) to the honeycomb lattice ($J=q$) to the square lattice (large $q$). The central numerical findings are that the magnetization variance marks a sharp phase-transition curve matching the predicted $J_c(q) = \\tanh^{-1}(\\sqrt{\\tanh^2 q + 1} - \\tanh q)$, and that for $q \\ge 2.5$ the equilibrium measure of the shaken dynamics is close to the square-lattice Gibbs measure. This matters because it would validate the companion theoretical papers and offer a single parallel sampler for a whole class of lattice models.","feed_headline":"Parallel 'shaken' spins hit the Ising critical curve","feed_subtitle":"A two-parameter parallel sampler interpolates between 1D, honeycomb, and square-lattice Ising models, validated numerically near…","key_machinery":"The central object is the shaken dynamics, a probabilistic cellular automaton defined as the marginal on one sublattice of an alternate parallel heat-bath dynamics on a bipartite graph. Each half-step updates all spins simultaneously with probabilities derived from a two-parameter Hamiltonian $H(\\sigma,\\tau) = -\\sum_x [J\\sigma_x(\\tau_{x^\\uparrow}+\\tau_{x^\\to})+q\\sigma_x\\tau_x]$, where $q$ is a self-interaction that tunes the effective lattice geometry. The load-bearing identity is the critical-curve equation $1 = 2\\tanh J\\tanh q + \\tanh^2 J$, whose solution is $J_c(q)$; it comes from the even-subgraph expansion for the honeycomb-lattice Ising partition function. The paper uses this identity to predict where the magnetization variance should peak, and uses the partial-order-preserving update to estimate mixing via coalescence of two extremal chains.","core_discovery":"The paper claims that the shaken dynamics—a factorized, fully parallel update rule on a bipartite lattice—undergoes an order-disorder phase transition along the explicit curve $J_c(q) = \\tanh^{-1}(\\sqrt{\\tanh^2 q + 1} - \\tanh q)$, which limits to the square-lattice critical value $\\tanh^{-1}(\\sqrt{2}-1)\\approx 0.4407$ as $q\\to\\infty$ and passes through the honeycomb critical point $J=q\\approx 0.6585$. Using the variance of the magnetization as a phase indicator on a $200\\times200$ torus, the authors locate the transition points over a grid in $(q,J)$ and find they fall on the theoretical curve. They also measure coalescence times to estimate mixing, and compare the shaken dynamics' equilibrium samples against the square-lattice Gibbs measure for magnetization and energy, concluding that for $q\\ge 2.5$ the approximation is good. Two parallel implementations (multicore CPU and GPU) are benchmarked, with the GPU version running roughly 500 times faster than a single CPU core for large lattices.","pith_inferences":["If the $q\\ge2.5$ approximation holds in total variation and not just for the observables tested, the shaken dynamics could be a drop-in parallel replacement for single-spin-flip samplers wherever the exact Gibbs measure is not essential.","The pattern of coalescence times suggests tuning $q$ could trade bias against mixing speed; testing intermediate $q$ values (e.g., $q\\in[1.5,2.5]$) by direct total-variation estimates would map that trade-off.","The even-subgraph critical-curve method might extend to other doubly periodic planar lattices, potentially producing a family of tunable parallel samplers for anisotropic Ising-type models.","A natural follow-up is to measure autocorrelation times rather than coalescence times, since the latter may be pessimistic for the variance of estimators in practice."],"forward_implications":["The shaken dynamics serves as one parallel MCMC algorithm that can simulate the whole family of Ising models across lattice geometries by tuning only $q$ and $J$.","The critical curve is confirmed numerically, so the phase diagram of the two-parameter Hamiltonian is established by both theory and simulation.","For $q\\ge2.5$, shaken-dynamics samples approximate square-lattice Gibbs samples well enough to estimate magnetization and energy at the accuracies shown.","Coalescence-time comparisons indicate the parallel alternate dynamics mixes faster per attempted spin flip than the single-spin-flip heat bath, even after volume renormalization.","The GPU implementation is roughly 500 times faster than a single CPU core for large lattices, making large-scale Monte Carlo runs practical."],"supporting_citations":[{"why":"Supplies the Hamiltonian, the critical-curve theorem expressed by eq. (11), and the theorems on magnetization and correlation asymmetry used in the numerical tests.","marker":"[1]"},{"why":"Introduces the shaken dynamics, its equilibrium measure $\\pi_s = Z_\\sigma/Z$, and the theorem that this measure approaches the square-lattice Gibbs measure for large $q$.","marker":"[2]"},{"why":"Gives the even-subgraph critical-temperature formula for Ising models on doubly periodic planar graphs, from which the critical-curve equation (10) is taken.","marker":"[3]"},{"why":"Provides the PCA construction for sampling pair-interaction Gibbs measures that the parallel update scheme is built on.","marker":"[4]"},{"why":"Defines coalescence times and the sandwiching technique for monotone chains, used to estimate mixing times in the simulations.","marker":"[8]"},{"why":"Supplies the Propp-Wilson coupling-from-the-past algorithm used to draw exact equilibrium samples in the comparisons.","marker":"[13]"}],"fun_headline_variants":["Shaken dynamics maps Ising phase boundaries on a GPU","Parallel Ising sampler finds critical curve with GPU speed","Two-parameter shaken Ising dynamics matches theory","GPU shaken Ising: 500x faster parallel sampling","Shaken Ising dynamics: phase curve validated on GPU"],"cache_read_input_tokens":3200,"weakest_assumption_plain":"The paper assumes the companion-paper theorems are correct—the critical-curve equation and the large-$q$ convergence of the shaken equilibrium measure to the square-lattice Gibbs measure—and assumes that 300,000 warm-up steps bring the chain to equilibrium.","fun_headline_variants_meta":{"raw":{"variants":["Shaken dynamics maps Ising phase boundaries on a GPU","Parallel Ising sampler finds critical curve with GPU speed","Two-parameter shaken Ising dynamics matches theory","GPU shaken Ising: 500x faster parallel sampling","Shaken Ising dynamics: phase curve validated on GPU"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000714,"raw_usage":{"total_tokens":3198,"prompt_tokens":921,"completion_tokens":2277,"prompt_tokens_details":{"cached_tokens":384},"prompt_cache_hit_tokens":384,"prompt_cache_miss_tokens":537,"completion_tokens_details":{"reasoning_tokens":2199}},"tokens_in":537,"tokens_out":2277,"duration_ms":17206,"temperature":1.0,"reasoning_tokens":2199,"cache_read_input_tokens":384,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-14T15:01:26.739360+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Directly estimate the total-variation distance $\\|\\pi_s - \\pi_G\\|_{TV}$ from Propp-Wilson samples for $q=2.5$ on increasing lattice sizes; if the distance does not tend to zero as $|\\Lambda|\\to\\infty$, or if the magnetization variance peaks away from $J_c(q)$ on lattices larger than $200\\times200$, the central claims would be refuted.","supporting_citations":[{"cited_title":"Criticality of measures on 2-d Ising configurations: from square to hexagonal graphs","cited_arxiv_id":"1906.02546","evidence_quote":"Supplies the Hamiltonian, the critical-curve theorem expressed by eq. (11), and the theorems on magnetization and correlation asymmetry used in the numerical tests."},{"cited_title":"Shaken dynamics: an easy way to parallel Markov Chain Monte Carlo","cited_arxiv_id":"1904.06257","evidence_quote":"Introduces the shaken dynamics, its equilibrium measure $\\pi_s = Z_\\sigma/Z$, and the theorem that this measure approaches the square-lattice Gibbs measure for large $q$."},{"cited_title":"The critical temperature for the Ising model on planar doubly periodic graphs","cited_arxiv_id":null,"evidence_quote":"Gives the even-subgraph critical-temperature formula for Ising models on doubly periodic planar graphs, from which the critical-curve equation (10) is taken."},{"cited_title":"Häggström and O.H.G","cited_arxiv_id":null,"evidence_quote":"Defines coalescence times and the sandwiching technique for monotone chains, used to estimate mixing times in the simulations."}],"review_version":1}