OpenMP parallelization of multiple precision Taylor series method

I. Hristov; I. Puzynin; N. Shegunov; R. Hristova; S. Dimova; T. Puzynina; Z. Sharipov; Z. Tukhliev

REVIEW 2 major objections 4 minor 10 references

OpenMP parallelization of multiple precision Taylor series method

T0 review · 2 major / 4 minor · reviewed 2026-08-14 · deepseek-v4-flash

Challenge this review Re-run · record.json Download PDF Read on arXiv ↗ Review JSON

Pith's one-line read A simple OpenMP parallelization of the multiple-precision Taylor series method achieves about 75 percent parallel efficiency on 28 cores, a 21x speedup, and reproduces a published Lorenz benchmark over [0,5000] on a single CPU node.

desk verdict Useful short engineering report on OpenMP parallel reduction for GMP-based Taylor integration; scaling numbers are believable but the 5-step benchmark is too thin to support the production-time claim. read the letter →

arxiv 1908.09301 v1 pith:DGLA6ZTN submitted 2019-08-25 cs.MS math.DS

S. Dimova , I. Hristov , R. Hristova , I. Puzynin , T. Puzynina , Z. Sharipov , N. Shegunov , Z. Tukhliev This is my paper · ORCID

classification cs.MSmath.DS MSC 65L0565Y05

keywords OpenMPTaylorseriesmethodmultipleprecisionarithmeticLorenzsystemparallelreductioncleannumericalsimulationchaoticdynamicalsystemsGMP

verification ladder T0 review T1 audit T2 compute T3 formal

The pith

A machine-rendered reading of the paper's core claim, the machinery that carries it, and where it could break.

The reading

This paper shows that the computationally heavy part of the multiple precision Taylor series method for chaotic ODEs—the recursive calculation of Taylor coefficients, which costs $O(N^2)$ operations—can be parallelized with a few OpenMP pragmas and a manual reduction. On the Lorenz system with $N=2000$ and about 2412 decimal digits, the authors measure about 75 percent parallel efficiency on 28 cores and a 21x speedup, and use this to integrate the equations over $[0,5000]$ in about nine days and fourteen hours on a single CPU node. The importance is that clean numerical simulation (the approach for mathematically reliable chaotic solutions proposed in [1]) becomes feasible on moderate shared-memory hardware rather than requiring large clusters. The same parallelization pattern applies to any dynamical system whose right-hand side is automatically differentiable.

What carries the argument

The central object is the Taylor coefficient recursion for the Lorenz system, where each new coefficient $\alpha_{i+1}$, $\beta_{i+1}$, $\gamma_{i+1}$ depends on all previous coefficients through two convolution sums $s_1 = \sum_{k=0}^i \alpha_{i-k}\gamma_k$ and $s_2 = \sum_{k=0}^i \alpha_{i-k}\beta_k$. These sums are the $O(N^2)$ bottleneck and are a classical parallel reduction. The paper's mechanism is a manual OpenMP reduction: because GMP's heap-allocated numbers cannot be privatized, each thread accumulates its partial sums in shared containers `sum[2*tid]` and `sum[2*tid+1]` with a private temporary variable, and after the parallel loop the master thread sums the containers. This turns the serial $O(N^2)$ recursion into a parallel reduction while leaving all other parts sequential.

What would settle it

Run the same $N=2000$, about 2412-digit method over a long interval (for example, $[0,1000]$) in both serial and 28-thread OpenMP modes, and compare the observed speedup with the 21x measured on the five-step test; if the long-run speedup falls well short of 21x or the per-step time grows with step count, the short-test scalability does not predict production performance.

Watch

Extended reading notes

Core claim

The central claim is that a straightforward OpenMP parallelization of the multiple precision Taylor series method for the Lorenz system achieves near-linear speedup on one shared-memory node, with about 75 percent parallel efficiency on 28 cores (21x speedup) for a short five-step test, and that this suffices for a full production integration over $[0,5000]$ in about 9 days and 14 hours, reproducing benchmark table S1 of [5]. The key obstacle is that GMP multiple-precision numbers are heap-allocated and cannot be privatized, so the standard OpenMP reduction clause does not apply; the authors instead store thread-local partial sums in a shared array of containers and manually unload them after each parallel loop. Because the parallelized part is only the convolution-like sums in the Taylor coefficient recursion, the method transfers directly to any system whose right-hand side is obtained by sums and products of elementary functions.

Load-bearing premise

The load-bearing premise is that the five-step timing test used for the scalability table is representative of the full-length production run over $[0,5000]$; if memory bandwidth, cache, or load-balance behavior differ over thousands of steps, the reported 75 percent parallel efficiency and 21x speedup would not transfer to the long integration.

Editorial extensions

If this is right

The same OpenMP pattern should transfer directly to any dynamical system whose right-hand side is built from sums, products, and compositions of elementary functions, because the Taylor coefficient recursion always produces convolution-type sums.
Long-term, mathematically reliable simulations (clean numerical simulation) become tractable on a single multicore node rather than requiring large MPI clusters.
The parallelized code reproduces benchmark table S1 of [5] on the interval $[0,5000]$ in about 9 days and 14 hours on 28 cores, offering an independent check of previously published reliable Lorenz solutions.
Even modest thread counts give substantial speedups: 4 cores yield 3.63x, 8 cores 6.92x, 14 cores 11.4x, and 28 cores 21.0x.

Reading between the lines

Editorial extensions of the paper, not claims the author makes directly.

The manual container-based reduction could be reused for other user-defined types that OpenMP cannot privatize, such as interval arithmetic libraries, extending the approach beyond GMP.
The five-step timing test used for scaling likely underestimates memory-bandwidth effects; a production-oriented study would measure per-step times over the entire $[0,5000]$ run.
Because Taylor steps are serially dependent, the next untested performance gain would be pipelining or overlapping consecutive steps across threads, which this paper does not attempt.

Share X Bluesky LinkedIn Reddit HN

Signed reviews

No signed human review yet.

Editorial analysis

A structured set of objections, weighed in public.

Desk editor's note, referee report, and a circularity audit.

Desk Editor's Note

read the letter

Punchline: the paper contributes a clean OpenMP implementation of the multiple-precision Taylor method with manual thread-local partial sums for GMP heap numbers, and reports 21x speedup on 28 cores. That is a real, if modest, extension of the MPI-based work of Wang et al. The serial baseline matches ref [5] within about 9%, so the measurements are plausible.

What's good: the recurrence (3) is standard and clearly written. The explanation of why built-in reduction can't be used (user-defined types and heap allocation) is correct, and the per-thread container scheme in Fig. 3 is straightforwardly described. The paper is honest that its goal is a moderate-resource node, not the 1200-core record. The short-test speedups in Table 1 are consistent with the algorithm: 90% efficiency at 4 cores decreasing to 75% at 28 is what you'd expect.

Soft spots, in order of importance:

1. The scalability table is measured on only 5 integration steps. For the [0,5000] production run, that's 500,000 steps, and memory behavior can change over sustained execution—cache/TLB, heap fragmentation, load balance. The paper asserts it integrated the full interval in ~9d14h and repeated table S1 from [5], but does not show timing curves, a log, or the reproduced table. As written, the production claim is an assertion, not a demonstrated result. The stress-test note is right to flag this.

2. Timings are single runs, no error bars. For parallel speedup numbers, run-to-run variance can be a few percent; that's not fatal here, but it would be cheap to repeat.

3. Comparing to ref [5] is not on matched hardware; the authors acknowledge this, and the comparison is 'our 28 cores beat their 50 cores,' which is suggestive but not rigorous. They should soften that language.

4. GMP thread-safety is cited to an OpenMP textbook rather than GMP docs; minor.

5. 'Large class of chaotic systems' is plausible because the method is automatic differentiation, but demonstrated on only Lorenz.

Overall: the central idea is sound and the engineering contribution is real. The gap between the short benchmark and the long-run claim is the main weakness. If the authors publish code/data and show the actual production run, this becomes a useful reference. As is, it's an honest but incomplete report.

I'd send it to peer review in a computational-science venue, with the expectation that the production claim be backed by more evidence. A serious referee would be worth the time.

Referee Report

2 major / 4 minor

Summary. The paper proposes and describes an OpenMP parallelization of a multiple-precision Taylor series integrator for the Lorenz system, in which the convolution-like sums appearing in the Taylor coefficient recurrence are computed by a manually implemented parallel reduction. The authors explain why OpenMP's built-in reduction clause cannot be used directly for GMP user-defined types, and they present container arrays for per-thread partial sums. They report scalability on one CPU node (28 cores), with a 21.0x speedup and 75.1% parallel efficiency on a short 5-integration-step test, and they state that this enabled integrating the Lorenz equations on [0,5000] in about 9 days and 14 hours, also claiming to have repeated benchmark table S1 from reference [5].

Significance. The paper provides a clear and potentially useful implementation recipe: using per-thread containers for partial sums and manual unloading, combined with the observation that GMP is thread-safe, is a simple way to parallelize the dominant reduction loops in Taylor-series integration on a shared-memory node. If the reported performance numbers are reproducible, the approach offers a moderate-resource route to long reliable chaotic simulations for systems with polynomial right-hand sides, complementing the large-scale MPI efforts in the cited literature. The derivation of the recurrence is standard and correctly stated, and the pseudocode in Figs. 1-3 is instructive. However, the headline performance and reproduction claims rest on weak evidence, as detailed below.

major comments (2)

[Section 5, Table 1] The central performance claims—75.1% parallel efficiency, 21.0x speedup, and the estimated ~9 days and 14 hours for the [0,5000] integration—are based solely on a short simulation test of 5 integration steps. The text states that the equations were integrated on [0,5000] in that time, but no wall-clock measurement of that long run is given; the reported time is consistent with a linear extrapolation of the 5-step timing (174.05 s / 5 steps * 500,000 steps / 21.0). Memory bandwidth, cache behavior, heap allocation patterns, and load balance can all change over 500,000 steps, so the efficiency measured on a cache-warm 5-step run cannot be assumed to transfer to the production run. Please either provide a direct timing for the full [0,5000] integration with details of the run, or soften the claims to what was actually measured.
[Section 5, comparison with reference [5]] The claim that the authors "repeated the benchmark table S1 from paper [5]" is not substantiated: no table, no numerical comparison of trajectories or computed quantities, and no discussion of verification or tolerances is provided. In addition, the comparison of parallel efficiency with reference [5] is performed on different hardware (2 x Intel Xeon E5-2695v3 here versus the platform used in [5]), so the assertion of better parallel efficiency is not a controlled comparison. Please show the actual reproduction data and qualify the cross-hardware comparison explicitly.

minor comments (4)

[Abstract and Section 1] The system in question is the Lorenz system; "Lorentz equations" should be corrected to "Lorenz equations" in the abstract and introduction.
[Section 2] "Salztman's parameter values" should be "Saltzman's parameter values" (or "Lorenz's standard parameter values").
[Table 1] Each row of Table 1 appears to be a single timing measurement. Reporting multiple runs and their variance (or at least stating that timings were stable) would make the scalability claim more robust.
[Section 5] The paper does not state whether the parallel implementation reproduces the serial solution bit-for-bit or to within the working precision; a brief statement on numerical reproducibility of the computed solution would be valuable.

Circularity Check

0 steps flagged · score 0.0 of 10

No circularity: the paper is a direct empirical measurement of OpenMP speedups, with no fitted input renamed as prediction and no load-bearing self-citation.

full rationale

The paper's central claims are performance measurements, not derived predictions. The speedup and parallel efficiency in Table 1 are obtained by timing the same parallel code against its own serial baseline under fixed parameters (N=2000, precision about 2412 decimal digits, 5 integration steps). No parameter is fitted to the reported result, and no quantity used as input is defined in terms of the output. The comparison with reference [5] is an external benchmark: the authors explicitly use the initial conditions, step size, and Table 1/S1 from Wang, Liu, and Li, who are not authors of this paper. The paper does not invoke any uniqueness theorem, nor does it import an ansatz from its own prior work; the parallel reduction scheme is explained directly from the serial recurrence (3) and the GMP code in Fig. 2. The skeptical concern that a 5-step timing test may not extrapolate to the full [0,5000] run is a legitimate robustness or correctness-risk issue about measurement representativeness, but it is not a circularity: the paper reports what was timed (5 steps) and separately reports the full production run, without dressing an extrapolation up as an independent prediction. Under the stated review rules, 'this is not standard consensus' and 'the benchmark may not generalize' are not circularity arguments. Accordingly, no step satisfies the requirement of reducing by construction to its own input, and a score of 0 is appropriate.

Assumptions & free parameters 0 free parameters · 3 assumptions · 0 invented entities

The central claim rests on the standard Taylor recurrence, the convergence behavior of high-order Taylor integration, and the thread-safety of GMP. No new entities or fitted constants are introduced; the Lorenz parameters and the order and precision choices are carried over from prior literature [1,5].

assumptions (3)

standard math Leibniz rule applied to the polynomial right-hand side of the Lorenz system yields the coefficient recurrence (3).
Section 2, equation (3). The recurrence is the standard foundation for Taylor series integration of polynomial ODEs.
domain assumption The solution is analytic on the integration interval so that a sufficiently high-order truncated Taylor series approximates it to the chosen precision for the estimated order N.
Section 1, following the clean numerical simulation framework of reference [1], which assumes the Taylor method converges with estimated N and precision.
domain assumption GMP library functions are thread-safe when called concurrently from distinct OpenMP threads.
Section 4 states 'An important feature of the GMP library is its thread-safety [9]' and the parallel reduction relies on concurrent mpf_mul and mpf_add calls. The cited reference [9] is an OpenMP textbook, not the GMP manual, which is a citation error.

how reviews work

0 comments

Cite this review

Pith. "Pith review of OpenMP parallelization of multiple precision Taylor series method." pith.science (2026). https://pith.science/paper/DGLA6ZTN

@misc{pith2026190809301,
  author       = {Pith},
  title        = {Pith review of: OpenMP parallelization of multiple precision Taylor series method},
  year         = {2026},
  howpublished = {\url{https://pith.science/paper/DGLA6ZTN}},
  note         = {Machine review of arXiv:1908.09301}
}

read the original abstract

OpenMP parallelization of multiple precision Taylor series method is proposed. A very good parallel performance scalability and parallel efficiency inside one computation node of a CPU-cluster is observed. We explain the details of the parallelization on the classical example of the Lorentz equations. The same approach can be applied straightforwardly to a large class of chaotic dynamical systems.

Discussion (0). Continue with ORCID to comment.

Reference graph

Works this paper leans on

10 extracted references · 10 canonical work pages

[5]

”Clean numerical simulation for some chaotic systems using the parallel multiple-precision Taylor scheme.” Chinese science bulletin 59.33 (2014): 4465-4472

Wang, Pengfei, Yong Liu, and Jianping Li. ”Clean numerical simulation for some chaotic systems using the parallel multiple-precision Taylor scheme.” Chinese science bulletin 59.33 (2014): 4465-4472

work page 2014
[1]

”On the reliability of computed chaotic solutions of non-linear differential equa- tions.” Tellus A: Dynamic Meteorology and Oceanography 61.4 (2008): 550-564

Liao, Shijun. ”On the reliability of computed chaotic solutions of non-linear differential equa- tions.” Tellus A: Dynamic Meteorology and Oceanography 61.4 (2008): 550-564

work page 2008
[2]

”A software package for the numerical integration of ODEs by means of high-order Taylor methods.” Experimental Mathematics 14.1 (2005): 99-117

Jorba, Angel, and Maorong Zou. ”A software package for the numerical integration of ODEs by means of high-order Taylor methods.” Experimental Mathematics 14.1 (2005): 99-117

work page 2005
[3]

”Breaking the limits: the Taylor series method.” Applied mathematics and computation 217.20 (2011): 7940-7954

Barrio, Roberto, et al. ”Breaking the limits: the Taylor series method.” Applied mathematics and computation 217.20 (2011): 7940-7954. 8 S. Dimova, I. Hristov, et al

work page 2011
[4]

Wang, Pengfei, Jianping Li, and Qian Li. ”Computational uncertainty and the application of a high-performance multiple precision scheme to obtaining the correct reference solution of Lorenz equations.” Numerical Algorithms 59.1 (2012): 147-159

work page 2012
[6]

Liao, ShiJun, and PengFei Wang. ”On the mathematically reliable long-term simulation of chaotic solutions of Lorenz equation in the interval [0, 10000].” Science China Physics, Me- chanics and Astronomy 57.2 (2014): 330-335

work page 2014
[7]

”Deterministic nonperiodic ﬂow.” Journal of the atmospheric sciences 20.2 (1963): 130-141

Lorenz, Edward N. ”Deterministic nonperiodic ﬂow.” Journal of the atmospheric sciences 20.2 (1963): 130-141

work page 1963
[8]

Methods and applications of interval analysis

Moore, Ramon E. Methods and applications of interval analysis. Society for Industrial and Applied Mathematics, 1979

work page 1979

Show all 10 references

[9]

Using OpenMP: portable shared memory parallel programming

Chapman, Barbara, Gabriele Jost, and Ruud Van Der Pas. Using OpenMP: portable shared memory parallel programming. V ol. 10. MIT press, 2008

2008
[11]

Avail- able at: http://hlit.jinr.ru (accessed 24.08.2019)

Multipurpose Information and Computing Complex of the Laboratory of IT of JINR. Avail- able at: http://hlit.jinr.ru (accessed 24.08.2019)

2019

Pith tools

Reviewed August 14, 2026 · model on record in the stance chip above.

[1] [5]

”Clean numerical simulation for some chaotic systems using the parallel multiple-precision Taylor scheme.” Chinese science bulletin 59.33 (2014): 4465-4472

Wang, Pengfei, Yong Liu, and Jianping Li. ”Clean numerical simulation for some chaotic systems using the parallel multiple-precision Taylor scheme.” Chinese science bulletin 59.33 (2014): 4465-4472

work page 2014

[2] [1]

”On the reliability of computed chaotic solutions of non-linear differential equa- tions.” Tellus A: Dynamic Meteorology and Oceanography 61.4 (2008): 550-564

Liao, Shijun. ”On the reliability of computed chaotic solutions of non-linear differential equa- tions.” Tellus A: Dynamic Meteorology and Oceanography 61.4 (2008): 550-564

work page 2008

[3] [2]

”A software package for the numerical integration of ODEs by means of high-order Taylor methods.” Experimental Mathematics 14.1 (2005): 99-117

Jorba, Angel, and Maorong Zou. ”A software package for the numerical integration of ODEs by means of high-order Taylor methods.” Experimental Mathematics 14.1 (2005): 99-117

work page 2005

[4] [3]

”Breaking the limits: the Taylor series method.” Applied mathematics and computation 217.20 (2011): 7940-7954

Barrio, Roberto, et al. ”Breaking the limits: the Taylor series method.” Applied mathematics and computation 217.20 (2011): 7940-7954. 8 S. Dimova, I. Hristov, et al

work page 2011

[5] [4]

Wang, Pengfei, Jianping Li, and Qian Li. ”Computational uncertainty and the application of a high-performance multiple precision scheme to obtaining the correct reference solution of Lorenz equations.” Numerical Algorithms 59.1 (2012): 147-159

work page 2012

[6] [6]

Liao, ShiJun, and PengFei Wang. ”On the mathematically reliable long-term simulation of chaotic solutions of Lorenz equation in the interval [0, 10000].” Science China Physics, Me- chanics and Astronomy 57.2 (2014): 330-335

work page 2014

[7] [7]

”Deterministic nonperiodic ﬂow.” Journal of the atmospheric sciences 20.2 (1963): 130-141

Lorenz, Edward N. ”Deterministic nonperiodic ﬂow.” Journal of the atmospheric sciences 20.2 (1963): 130-141

work page 1963

[8] [8]

Methods and applications of interval analysis

Moore, Ramon E. Methods and applications of interval analysis. Society for Industrial and Applied Mathematics, 1979

work page 1979