From 30c32c39ab45d2102369f530e303ae7711407687 Mon Sep 17 00:00:00 2001 From: Evan Miller Date: Thu, 17 Sep 2026 06:50:53 -0400 Subject: [PATCH 1/3] Form the Jacobi theta exponents without amplified rounding Every term of the theta series is exp(-E) where E is a product or quotient of tau, pi and a small integer, or (for tau < 1, after the modular transformation) the square of z plus a multiple of pi/2, divided by pi*tau. Rounding E to working precision costs |E| ulps in the term, and the exponents run into the hundreds before the terms underflow, so the functions lost tens to hundreds of ulps wherever the exponents were large: near tau = 1 in the direct series, and for small tau throughout, where the argument reduction by a rounded pi added an absolute error that the Gaussians amplified by 2z/(pi*tau). The exponents are now formed by a small object that holds pi and a = pi*tau (or -ln q) as unevaluated sums hi + lo and tracks the rounding of every product, sum and quotient with Dekker's and Knuth's error-free transformations, which work for any binary floating-point type (an fma is used where the standard library advertises a fast one). The error is applied to each term as a factor (1 - dE). The argument reduction uses the same split pi and passes the remainder's error into the Gaussians, and the z = 0 shortcuts form the nome of 1/tau without rounding the reciprocal. The q parameterization builds a = -ln q directly instead of going through tau = -ln(q)/pi and back, halving its remaining error. Against a 50-digit reference (double, z in [0, 5], tau in [0.005, 31]): tau parameterization, tau < 1: worst 177-360 ulp -> 3-6, mean 6-12 -> 0.5 tau parameterization, tau >= 1: worst 20-76 ulp -> 2.5, mean 2-8 -> 0.5 q parameterization, q > e^-pi: worst 380-790 ulp -> 64-81, mean 12-14 -> 3 What remains in the q parameterization is the rounding of ln q, which is inherent to it. The "minus one" functions below tau = 1 still cancel against 1, as documented. The Gaussian-sum paths cost about twice what they did; the direct series about a quarter more. The periodicity test shifts z by the rounded two_pi, which for large q now changes the function by more than the evaluation error, so its tolerance accounts for that; the Mellin transform integrals get a couple more ulps of slack. The q sweep of the ULP plot generator started at q = 0, which has been a domain error since the last round of fixes. The ULP graphs are regenerated. Co-Authored-By: Claude Fable 5.1 --- doc/graphs/jacobi_theta1_float.svg | 5037 +++++++++-------- doc/graphs/jacobi_theta1q_float.svg | 5003 ++++++++-------- doc/graphs/jacobi_theta2_float.svg | 5012 ++++++++-------- doc/graphs/jacobi_theta2q_float.svg | 5025 ++++++++-------- doc/graphs/jacobi_theta3_float.svg | 4969 ++++++++-------- doc/graphs/jacobi_theta3q_float.svg | 5025 ++++++++-------- doc/graphs/jacobi_theta4_float.svg | 4978 ++++++++-------- doc/graphs/jacobi_theta4q_float.svg | 5003 ++++++++-------- doc/sf/jacobi_theta.qbk | 12 +- .../math/special_functions/jacobi_theta.hpp | 531 +- reporting/accuracy/plot_jacobi_theta_q.cpp | 8 +- test/test_jacobi_theta.cpp | 11 +- 12 files changed, 20499 insertions(+), 20115 deletions(-) diff --git a/doc/graphs/jacobi_theta1_float.svg b/doc/graphs/jacobi_theta1_float.svg index e1b87a2c1f..b08f0a03ff 100644 --- a/doc/graphs/jacobi_theta1_float.svg +++ b/doc/graphs/jacobi_theta1_float.svg @@ -6,23 +6,27 @@ svg { background-color:black; } jacobi_theta1(x, 0.5) ULP plot at float precision - - --75 - --50 - --25 - -0 - -25 - -50 - -75 - -100 + + +-2.5 + +-2 + +-1.5 + +-1 + +-0.5 + +0.5 + +1 + +1.5 + +2 + +2.5 0.6283 @@ -43,2498 +47,2511 @@ svg { background-color:black; } 5.655 6.283 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/doc/graphs/jacobi_theta1q_float.svg b/doc/graphs/jacobi_theta1q_float.svg index 63dcb428df..41402c7561 100644 --- a/doc/graphs/jacobi_theta1q_float.svg +++ b/doc/graphs/jacobi_theta1q_float.svg @@ -6,23 +6,23 @@ svg { background-color:black; } jacobi_theta1(5.0, q) ULP plot at float precision - + --87.28 +-14.75 --74.57 +-9.682 --61.85 +-4.611 --49.13 +0.4611 --36.42 +5.533 --23.7 +10.6 --10.98 +15.68 -1.733 +20.75 0.1 @@ -43,2486 +43,2509 @@ svg { background-color:black; } 0.9 1 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/doc/graphs/jacobi_theta2_float.svg b/doc/graphs/jacobi_theta2_float.svg index 0d6f4b748f..ff27cf0d52 100644 --- a/doc/graphs/jacobi_theta2_float.svg +++ b/doc/graphs/jacobi_theta2_float.svg @@ -6,21 +6,21 @@ svg { background-color:black; } jacobi_theta2(x, 0.5) ULP plot at float precision - + --75 +-12.81 --50 +3.307 --25 +19.42 -0 +35.54 -25 +51.65 -50 +67.77 -75 +83.88 100 @@ -43,2495 +43,2511 @@ svg { background-color:black; } 5.655 6.283 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/doc/graphs/jacobi_theta2q_float.svg b/doc/graphs/jacobi_theta2q_float.svg index 367d85e5b2..4d576aec55 100644 --- a/doc/graphs/jacobi_theta2q_float.svg +++ b/doc/graphs/jacobi_theta2q_float.svg @@ -6,23 +6,23 @@ svg { background-color:black; } jacobi_theta2(0.4, q) ULP plot at float precision - + --85.25 +-47.88 --70.49 +-40.64 --55.74 +-33.4 --40.99 +-26.16 --26.24 +-18.92 --11.48 +-11.68 -3.27 +-4.439 -18.02 +2.802 0.1 @@ -43,2508 +43,2509 @@ svg { background-color:black; } 0.9 1 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/doc/graphs/jacobi_theta3_float.svg b/doc/graphs/jacobi_theta3_float.svg index 0490844afb..4587497615 100644 --- a/doc/graphs/jacobi_theta3_float.svg +++ b/doc/graphs/jacobi_theta3_float.svg @@ -43,2474 +43,2515 @@ svg { background-color:black; } 5.655 6.283 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + - - - - - - - - - - - - - - - - - - - - - - - - - - + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/doc/graphs/jacobi_theta3q_float.svg b/doc/graphs/jacobi_theta3q_float.svg index a9753c0ff2..cc1e145bf8 100644 --- a/doc/graphs/jacobi_theta3q_float.svg +++ b/doc/graphs/jacobi_theta3q_float.svg @@ -6,23 +6,23 @@ svg { background-color:black; } jacobi_theta3(0.4, q) ULP plot at float precision - + --87.17 +-28.75 --74.35 +-23.4 --61.52 +-18.06 --48.69 +-12.71 --35.87 +-7.364 --23.04 +-2.018 --10.21 +3.328 -2.613 +8.674 0.1 @@ -43,2508 +43,2509 @@ svg { background-color:black; } 0.9 1 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/doc/graphs/jacobi_theta4_float.svg b/doc/graphs/jacobi_theta4_float.svg index 4b25f89eed..fb170075b0 100644 --- a/doc/graphs/jacobi_theta4_float.svg +++ b/doc/graphs/jacobi_theta4_float.svg @@ -43,2481 +43,2515 @@ svg { background-color:black; } 5.655 6.283 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + - - - - - - - - - - - - - - - - - - - - + + + + + + + + + + + + + + + + + + + + + + + diff --git a/doc/graphs/jacobi_theta4q_float.svg b/doc/graphs/jacobi_theta4q_float.svg index 2910ee71e3..3a2e759d1b 100644 --- a/doc/graphs/jacobi_theta4q_float.svg +++ b/doc/graphs/jacobi_theta4q_float.svg @@ -6,23 +6,23 @@ svg { background-color:black; } jacobi_theta4(5.0, q) ULP plot at float precision - + -10.99 +-7.544 -23.71 +-3.926 -36.42 +-0.307 -49.14 +3.312 -61.85 +6.93 -74.57 +10.55 -87.28 +14.17 -100 +17.79 0.1 @@ -43,2486 +43,2509 @@ svg { background-color:black; } 0.9 1 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/doc/sf/jacobi_theta.qbk b/doc/sf/jacobi_theta.qbk index d5973d8d5d..ae88040f12 100644 --- a/doc/sf/jacobi_theta.qbk +++ b/doc/sf/jacobi_theta.qbk @@ -55,6 +55,10 @@ A more accurate computation will take advantage of [tau]: Internally, when /q/ is larger than exp(-[pi]) (that is, when [tau] is less than 1), Boost implements the /q/ parameterization by taking the logarithm of /q/ and passing it to the [tau] parameterization; as such, using the [tau] parameterization directly will generally yield greater precision in that regime. When /q/ is smaller than exp(-[pi]), the Fourier series is summed directly in terms of /q/, and the two parameterizations are equally accurate. + +Every term of the series is an exponential whose argument is a product or quotient of [tau], [pi] and a small integer, or, in the transformed series used when [tau] < 1, the square of /z/ plus a multiple of [pi]/2. +Rounding such an argument to working precision would cost as many ulps in the term as the argument is large, so the rounding errors of [pi], of the products, of the division and of the reduction of /z/ are tracked exactly (using error-free transformations that work for any binary floating-point type) and applied as a correction factor to each term. +The [tau] parameterization is therefore accurate to a few ulps throughout, apart from the conditioning of the function itself. In the /q/ parameterization the rounding of log(/q/) remains, and its effect grows as /q/ approaches 1. As another example, if the complement of /q/ is known with great accuracy, then instead of: jacobi_theta1(x, 1-q_complement); @@ -125,7 +129,7 @@ Fixing /x/=5 and varying /q/, the ULPs plot looks like: [graph jacobi_theta1q_float] -Accuracy tends to degenerate near /q/=1 (small [tau]). +Accuracy degrades gently as /q/ approaches 1 (small [tau]), where the rounding of log(/q/) is amplified by the size of the exponents in the transformed series; the [tau] parameterization does not have this problem. [heading Implementation] @@ -194,7 +198,7 @@ Fixing /x/=0.4 and varying /q/, the ULPs plot looks like: [graph jacobi_theta2q_float] -Accuracy tends to degenerate near /q/=1 (small [tau]). +Accuracy degrades gently as /q/ approaches 1 (small [tau]), where the rounding of log(/q/) is amplified by the size of the exponents in the transformed series; the [tau] parameterization does not have this problem. [heading Implementation] @@ -276,7 +280,7 @@ Fixing /x/=0.4 and varying /q/, the ULPs plot looks like: [graph jacobi_theta3q_float] -Accuracy tends to degenerate near /q/=1 (small [tau]). +Accuracy degrades gently as /q/ approaches 1 (small [tau]), where the rounding of log(/q/) is amplified by the size of the exponents in the transformed series; the [tau] parameterization does not have this problem. [heading Implementation] @@ -358,7 +362,7 @@ Fixing /x/=5 and varying /q/, the ULPs plot looks like: [graph jacobi_theta4q_float] -Accuracy tends to degenerate near /q/=1 (small [tau]). +Accuracy degrades gently as /q/ approaches 1 (small [tau]), where the rounding of log(/q/) is amplified by the size of the exponents in the transformed series; the [tau] parameterization does not have this problem. [heading Implementation] diff --git a/include/boost/math/special_functions/jacobi_theta.hpp b/include/boost/math/special_functions/jacobi_theta.hpp index ab8aa8db97..4512002cf3 100644 --- a/include/boost/math/special_functions/jacobi_theta.hpp +++ b/include/boost/math/special_functions/jacobi_theta.hpp @@ -100,6 +100,7 @@ #ifndef BOOST_MATH_JACOBI_THETA_HPP #define BOOST_MATH_JACOBI_THETA_HPP +#include #include #include #include @@ -228,12 +229,231 @@ struct _jacobi_theta_q_power { } }; +// When the caller supplies tau, every term is an exponential exp(-E) whose +// argument E is a product or quotient of tau, pi and a small integer (or, for +// tau < 1, the square of z plus a multiple of pi/2). Rounding E to working +// precision costs |E| ulps in the exponential, so instead the rounding +// errors of pi, of the products and of the division are tracked exactly with +// Dekker's product and Knuth's sum (which work for any binary floating-point +// type) and applied as a correction factor exp(-dE) = 1 - dE. What remains +// is the rounding of tau itself, which is the caller's. +// (Cached per type: for a type whose precision can change at run time the +// cached value may split at the wrong place, in which case the corrections +// below are merely inexact, i.e. no worse than not applying them.) +template +inline RealType _jacobi_theta_splitter() { + BOOST_MATH_STD_USING + static const RealType splitter = ldexp(RealType(1), (tools::digits() + 1) / 2) + 1; + return splitter; +} + +// Veltkamp's splitting: t = hi + lo where hi holds only the leading half of +// the significand, so that products of hi parts are exact. +template +inline void _jacobi_theta_split(RealType t, RealType splitter, RealType& hi, RealType& lo) { + RealType g = splitter * t; + hi = g - (g - t); + lo = t - hi; +} + +// The exact product a*b is ab + (return value), where ab = fl(a*b). Where +// the standard library advertises a fast fused multiply-add it is used +// instead, since fma(a, b, -ab) is exactly this quantity. +template +inline RealType _jacobi_theta_product_error(RealType a, RealType b, RealType ab, RealType splitter) { + RealType ah, al, bh, bl; + _jacobi_theta_split(a, splitter, ah, al); + _jacobi_theta_split(b, splitter, bh, bl); + return (((ah * bh - ab) + ah * bl) + al * bh) + al * bl; +} +#ifdef FP_FAST_FMAF +inline float _jacobi_theta_product_error(float a, float b, float ab, float) { + return std::fma(a, b, -ab); +} +#endif +#ifdef FP_FAST_FMA +inline double _jacobi_theta_product_error(double a, double b, double ab, double) { + return std::fma(a, b, -ab); +} +#endif +#ifdef FP_FAST_FMAL +inline long double _jacobi_theta_product_error(long double a, long double b, long double ab, long double) { + return std::fma(a, b, -ab); +} +#endif + +// The exact sum a+b is s + (return value), where s = fl(a+b). +template +inline RealType _jacobi_theta_sum_error(RealType a, RealType b, RealType s) { + RealType bb = s - a; + return (a - (s - bb)) + (b - bb); +} + +// Pi as hi + lo: hi is exactly representable and lo completes it to about +// 120 bits, as five pieces of 24 significant bits each. Types with more +// precision than that use their own rounding of pi with lo = 0. +template +inline RealType _jacobi_theta_pi(RealType& lo) { + BOOST_MATH_STD_USING + if (tools::digits() > 116) { + lo = 0; + return constants::pi(); + } + static const RealType pi_lo = ldexp(RealType(10625384), -46) + ldexp(RealType(12727492), -70) + ldexp(RealType(13001355), -94) + ldexp(RealType(8444956), -118); + static const RealType pi_hi = ldexp(RealType(13176794), -22); + lo = pi_lo; + return pi_hi; +} + +// The quantity everything is expressed in is a = pi*tau = -ln(q), held as +// a_hi + a_lo. Built from tau it is exact (to second order); built from q it +// carries the rounding of the logarithm, which is inherent to that +// parameterization, but nothing else. +template +struct _jacobi_theta_exponents { + RealType tau; // only used for the scale factor 1/sqrt(tau) + RealType splitter; + RealType pi_hi, pi_lo; // pi = pi_hi + pi_lo + RealType a_hi, a_lo; // pi * tau = -ln(q) = a_hi + a_lo + RealType inv_a_hi; // 1 / a_hi, only needed when tau < 1 + + static _jacobi_theta_exponents from_tau(RealType tau) { + _jacobi_theta_exponents x; + x.tau = tau; + x.a_hi = x.pi_times(tau, x.a_lo); + x.inv_a_hi = (tau < 1) ? 1 / x.a_hi : RealType(0); + return x; + } + + static _jacobi_theta_exponents from_nome(RealType q) { + BOOST_MATH_STD_USING + _jacobi_theta_exponents x; + x.a_hi = -log(q); + x.a_lo = 0; + x.tau = x.a_hi / constants::pi(); + x.inv_a_hi = (x.tau < 1) ? 1 / x.a_hi : RealType(0); + return x; + } + + // pi * e as hi + lo for a small exact e + RealType pi_times(RealType e, RealType& lo) const { + RealType p = pi_hi * e; + RealType q = pi_lo * e; + RealType hi = p + q; + // |q| << |p|, so the sum's error is simply q - (hi - p) + lo = (q - (hi - p)) + _jacobi_theta_product_error(pi_hi, e, p, splitter); + return hi; + } + + // (N + dN) / (a_hi + a_lo) as E + dE, to second order + RealType divide_by_a(RealType N, RealType dN, RealType& dE) const { + RealType E = N * inv_a_hi; + RealType Ea = E * a_hi; + RealType r = (N - Ea) - _jacobi_theta_product_error(E, a_hi, Ea, splitter); + dE = (r + dN - E * a_lo) * inv_a_hi; + return E; + } + + // exp(-a * e) = q^e for a small exact e such as n^2 or (n+1/2)^2 + RealType direct(RealType e) const { + BOOST_MATH_STD_USING + RealType E = a_hi * e; + RealType result = exp(-E); + if (result == 0) + return result; + RealType dE = _jacobi_theta_product_error(a_hi, e, E, splitter) + a_lo * e; + return result * (1 - dE); + } + + // exp(-pi * e / tau) = exp(-pi^2 * e / a), the nome of 1/tau raised to e + RealType inverse(RealType e) const { + BOOST_MATH_STD_USING + RealType dN, dP; + RealType N = pi_times(e, dN); + RealType P = pi_times(N, dP); + dP += pi_hi * dN; + RealType dE; + RealType E = divide_by_a(P, dP, dE); + RealType result = exp(-E); + if (result == 0) + return result; + return result * (1 - dE); + } + + // Reduces z by the nearest multiple of the period (half_pis * pi/2), + // returning the remainder along with its error dz (the remainder is only + // exact with respect to the rounded pi, and the Gaussians below amplify + // that error by 2z/(pi*tau)), and the multiple k for the caller's sign. + RealType reduce(RealType z, int half_pis, RealType& dz, RealType& k) const { + BOOST_MATH_STD_USING + RealType period = pi_hi * RealType(half_pis) / 2; + k = floor(z / period + RealType(0.5)); + RealType dc; + RealType c = pi_times(k * RealType(half_pis) / 2, dc); + RealType r = z - c; + dz = _jacobi_theta_sum_error(z, -c, r) - dc; + return r; + } + + // exp(-(z + m*pi/2)^2 / (pi * tau)) for an integer m, where dz is the + // error in z from the argument reduction + RealType gaussian(RealType z, RealType dz, int m) const { + BOOST_MATH_STD_USING + RealType z_n = z; + if (m != 0) { + RealType dc; + RealType c = pi_times(RealType(m), dc) / 2; + z_n = z + c; + dz += _jacobi_theta_sum_error(z, c, z_n) + dc / 2; + } + RealType s = z_n * z_n; + RealType ds = _jacobi_theta_product_error(z_n, z_n, s, splitter) + 2 * z_n * dz; + RealType dE; + RealType E = divide_by_a(s, ds, dE); + RealType result = exp(-E); + if (result == 0) + return result; + return result * (1 - dE); + } + + // expm1(-2 * z * k / tau) = expm1(-2 * pi * z * k / a) for a small exact + // k, where dz is the error in z + RealType expm1_scaled(RealType z, RealType dz, RealType k) const { + BOOST_MATH_STD_USING + RealType P = z * k; + RealType dP = _jacobi_theta_product_error(z, k, P, splitter) + dz * k; + RealType dN; + RealType N = pi_times(P, dN); + dN += pi_hi * dP; + RealType dQ; + RealType Q = divide_by_a(N, dN, dQ); + RealType result = boost::math::expm1(-2 * Q); + if (result == -1) + return result; + return result - 2 * dQ * (1 + result); + } + +private: + _jacobi_theta_exponents() : splitter(_jacobi_theta_splitter()) { + pi_hi = _jacobi_theta_pi(pi_lo); + } +}; + template struct _jacobi_theta_tau_power { - RealType tau; + const _jacobi_theta_exponents& x; RealType operator()(RealType exponent) const { - BOOST_MATH_STD_USING - return exp(-tau * constants::pi() * exponent); + return x.direct(exponent); + } +}; + +// exp(-pi * e / tau): the nome of 1/tau, used by the modular transformation +// at z = 0 without forming the rounded reciprocal. +template +struct _jacobi_theta_inverse_tau_power { + const _jacobi_theta_exponents& x; + RealType operator()(RealType exponent) const { + return x.inverse(exponent); } }; @@ -327,18 +547,19 @@ _jacobi_theta4m1_series(RealType z, const NomePower& q_pow, const Policy&) { return result; } +// SUM exp(-(z + m*pi/2)^2 / (pi*tau)) over m = m0, m0 + m_step, ... until the +// terms are negligible. template inline RealType -_jacobi_theta_sum(RealType tau, RealType z_n, RealType z_increment, RealType eps) { - BOOST_MATH_STD_USING +_jacobi_theta_sum(const _jacobi_theta_exponents& x, RealType z, RealType dz, int m, int m_step, RealType eps) { RealType delta = 0, partial_result = 0; RealType last_delta = 0; do { last_delta = delta; - delta = exp(-tau*z_n*z_n/constants::pi()); + delta = x.gaussian(z, dz, m); partial_result += delta; - z_n += z_increment; + m += m_step; } while (!_jacobi_theta_converged(last_delta, delta, eps)); return partial_result; @@ -347,10 +568,10 @@ _jacobi_theta_sum(RealType tau, RealType z_n, RealType z_increment, RealType eps // The following _IMAGINARY theta functions assume imaginary z and are for // internal use only. They are designed to increase accuracy and reduce the // number of iterations required for convergence for large |q|. The z argument -// is scaled by tau, and the summations are rewritten to be double-sided +// is scaled by 1/tau, and the summations are rewritten to be double-sided // following DLMF 20.13.4 and 20.13.5. Each term is a Gaussian -// exp(-tau*(z - c)^2/Pi) centered at a multiple of Pi or Pi/2, and the -// results are scaled by sqrt(tau). +// exp(-(z - c)^2/(Pi*tau)) centered at a multiple of Pi or Pi/2, and the +// results are scaled by 1/sqrt(tau). // // These functions are triggered when tau < 1, i.e. |q| > exp(-Pi) = 0.043 // @@ -358,27 +579,27 @@ _jacobi_theta_sum(RealType tau, RealType z_n, RealType z_increment, RealType eps // vice-versa). jacobi_theta1 and jacobi_theta3 use the imaginary versions of // themselves, following DLMF 20.7.30 - 20.7.33. -// theta1(z|i/tau) = sqrt(tau) * SUM_{n>=0} (-1)^n [G(z - c_n) - G(z + c_n)] -// with c_n = Pi*(n+1/2) and G(x) = exp(-tau*x^2/Pi). +// theta1(z|tau) = 1/sqrt(tau) * SUM_{n>=0} (-1)^n [G(z - c_n) - G(z + c_n)] +// with c_n = Pi*(n+1/2) and G(x) = exp(-x^2/(Pi*tau)). // // Each bracket is a difference of two Gaussians which nearly cancel when z is // small, so it is evaluated instead through the exact identity -// G(z - c) - G(z + c) = -G(z - c) * expm1(-2*tau*z*(2n+1)), +// G(z - c) - G(z + c) = -G(z - c) * expm1(-2*z*(2n+1)/tau), // which keeps full relative precision all the way down to z -> 0. -// Requires 0 <= z <= Pi/2; the caller reduces z into this range. +// Requires 0 <= z <= Pi/2; the caller reduces z into this range, and +// passes the error dz of the reduced z. template inline RealType -_IMAGINARY_jacobi_theta1tau(RealType z, RealType tau, const Policy& pol) { +_IMAGINARY_jacobi_theta1tau(RealType z, RealType dz, const _jacobi_theta_exponents& x, const Policy&) { BOOST_MATH_STD_USING RealType eps = policies::get_epsilon(); - RealType result = 0, g = 0, last_g, c, pair; + RealType result = 0, g = 0, last_g, pair; unsigned n = 0; do { last_g = g; - c = constants::pi() * RealType(n + 0.5); - g = exp(-tau * (z - c) * (z - c) / constants::pi()); - pair = -g * boost::math::expm1(RealType(-2 * tau * z * RealType(2*n + 1)), pol); + g = x.gaussian(z, dz, -static_cast(2*n + 1)); + pair = -g * x.expm1_scaled(z, dz, RealType(2*n + 1)); if (n%2) pair = -pair; @@ -386,65 +607,163 @@ _IMAGINARY_jacobi_theta1tau(RealType z, RealType tau, const Policy& pol) { n++; } while (!_jacobi_theta_converged(last_g, g, eps)); - return result * sqrt(tau); + return result / sqrt(x.tau); } template inline RealType -_IMAGINARY_jacobi_theta2tau(RealType z, RealType tau, const Policy&) { +_IMAGINARY_jacobi_theta2tau(RealType z, RealType dz, const _jacobi_theta_exponents& x, const Policy&) { BOOST_MATH_STD_USING RealType eps = policies::get_epsilon(); RealType result = RealType(0); - // n>=0 - result += _jacobi_theta_sum(tau, RealType(z + constants::half_pi()), constants::pi(), eps); + // n>=0: centers at z + Pi/2 + n*Pi + result += _jacobi_theta_sum(x, z, dz, 1, 2, eps); // n<0 - result += _jacobi_theta_sum(tau, RealType(z - constants::half_pi()), RealType (-constants::pi()), eps); + result += _jacobi_theta_sum(x, z, dz, -1, -2, eps); - return result * sqrt(tau); + return result / sqrt(x.tau); } template inline RealType -_IMAGINARY_jacobi_theta3tau(RealType z, RealType tau, const Policy&) { +_IMAGINARY_jacobi_theta3tau(RealType z, RealType dz, const _jacobi_theta_exponents& x, const Policy&) { BOOST_MATH_STD_USING RealType eps = policies::get_epsilon(); RealType result = 0; // n=0 - result += exp(-z*z*tau/constants::pi()); - // n>0 - result += _jacobi_theta_sum(tau, RealType(z + constants::pi()), constants::pi(), eps); + result += x.gaussian(z, dz, 0); + // n>0: centers at z + n*Pi + result += _jacobi_theta_sum(x, z, dz, 2, 2, eps); // n<0 - result += _jacobi_theta_sum(tau, RealType(z - constants::pi()), RealType(-constants::pi()), eps); + result += _jacobi_theta_sum(x, z, dz, -2, -2, eps); - return result * sqrt(tau); + return result / sqrt(x.tau); } template inline RealType -_IMAGINARY_jacobi_theta4tau(RealType z, RealType tau, const Policy&) { +_IMAGINARY_jacobi_theta4tau(RealType z, RealType dz, const _jacobi_theta_exponents& x, const Policy&) { BOOST_MATH_STD_USING RealType eps = policies::get_epsilon(); RealType result = 0; // n = 0 - result += exp(-z*z*tau/constants::pi()); + result += x.gaussian(z, dz, 0); - // n > 0 odd - result -= _jacobi_theta_sum(tau, RealType(z + constants::pi()), constants::two_pi(), eps); + // n > 0 odd: centers at z + Pi + 2n*Pi + result -= _jacobi_theta_sum(x, z, dz, 2, 4, eps); // n < 0 odd - result -= _jacobi_theta_sum(tau, RealType(z - constants::pi()), RealType (-constants::two_pi()), eps); - // n > 0 even - result += _jacobi_theta_sum(tau, RealType(z + constants::two_pi()), constants::two_pi(), eps); + result -= _jacobi_theta_sum(x, z, dz, -2, -4, eps); + // n > 0 even: centers at z + 2*Pi + 2n*Pi + result += _jacobi_theta_sum(x, z, dz, 4, 4, eps); // n < 0 even - result += _jacobi_theta_sum(tau, RealType(z - constants::two_pi()), RealType (-constants::two_pi()), eps); + result += _jacobi_theta_sum(x, z, dz, -4, -4, eps); - return result * sqrt(tau); + return result / sqrt(x.tau); } -// First Jacobi theta function (Parameterized by tau - assumed imaginary) +// Dispatch on the size of tau (i.e. of the nome): the direct Fourier series +// for tau >= 1, otherwise the modular transformation to 1/tau. At z = 0 the +// transformed series is evaluated directly (single-sided, with the nome of +// 1/tau formed without rounding the reciprocal); otherwise as double-sided +// Gaussian sums. + // = 2 * SUM (-1)^n * exp(i*Pi*Tau*(n+1/2)^2) * sin((2n+1)z) +template +inline RealType +_jacobi_theta1_dispatch(RealType z, const _jacobi_theta_exponents& x, const Policy& pol) { + BOOST_MATH_STD_USING + if (x.tau < 1.0) { + // Reduce to -Pi/2 <= z <= Pi/2 using theta1(z + Pi) = -theta1(z)... + RealType dz, k; + z = x.reduce(z, 2, dz, k); + RealType sign = (fmod(k, RealType(2)) == 0) ? 1 : -1; + // ...and then to 0 <= z <= Pi/2 since theta1 is odd. + if (z < 0) { + z = -z; + dz = -dz; + sign = -sign; + } + return sign * _IMAGINARY_jacobi_theta1tau(z, dz, x, pol); + } + return _jacobi_theta1_series(z, _jacobi_theta_tau_power{x}, pol); +} + +// = 2 * SUM exp(i*Pi*Tau*(n+1/2)^2) * cos((2n+1)z) +template +inline RealType +_jacobi_theta2_dispatch(RealType z, const _jacobi_theta_exponents& x, const Policy& pol) { + BOOST_MATH_STD_USING + if (x.tau < 1.0 && abs(z) == 0.0) { // theta4(0|1/tau)/sqrt(tau) + return (RealType(1) + _jacobi_theta4m1_series(z, _jacobi_theta_inverse_tau_power{x}, pol)) / sqrt(x.tau); + } else if (x.tau < 1.0) { // DLMF 20.7.31 + // Reduce to -Pi <= z <= Pi (theta2 has period 2*Pi) + RealType dz, k; + z = x.reduce(z, 4, dz, k); + return _IMAGINARY_jacobi_theta4tau(z, dz, x, pol); + } + return _jacobi_theta2_series(z, _jacobi_theta_tau_power{x}, pol); +} + +// = 1 + 2 * SUM exp(i*Pi*Tau*(n)^2) * cos(2nz) +template +inline RealType +_jacobi_theta3_dispatch(RealType z, const _jacobi_theta_exponents& x, const Policy& pol) { + BOOST_MATH_STD_USING + if (x.tau < 1.0 && abs(z) == 0.0) { // theta3(0|1/tau)/sqrt(tau) + return (RealType(1) + _jacobi_theta3m1_series(z, _jacobi_theta_inverse_tau_power{x}, pol)) / sqrt(x.tau); + } else if (x.tau < 1.0) { // DLMF 20.7.32 + // Reduce to -Pi/2 <= z <= Pi/2 (theta3 has period Pi) + RealType dz, k; + z = x.reduce(z, 2, dz, k); + return _IMAGINARY_jacobi_theta3tau(z, dz, x, pol); + } + return RealType(1) + _jacobi_theta3m1_series(z, _jacobi_theta_tau_power{x}, pol); +} + +// = 2 * SUM exp(i*Pi*Tau*(n)^2) * cos(2nz), n >= 1 (theta3 minus one) +// This preserves accuracy for small values of q (i.e. tau > 1). For larger +// values of q, the minus one version usually won't help. +template +inline RealType +_jacobi_theta3m1_dispatch(RealType z, const _jacobi_theta_exponents& x, const Policy& pol) { + if (x.tau < 1.0) + return _jacobi_theta3_dispatch(z, x, pol) - RealType(1); + return _jacobi_theta3m1_series(z, _jacobi_theta_tau_power{x}, pol); +} + +// = 1 + 2 * SUM (-1)^n exp(i*Pi*Tau*(n)^2) * cos(2nz) +template +inline RealType +_jacobi_theta4_dispatch(RealType z, const _jacobi_theta_exponents& x, const Policy& pol) { + BOOST_MATH_STD_USING + if (x.tau < 1.0 && abs(z) == 0.0) { // theta2(0|1/tau)/sqrt(tau) + return _jacobi_theta2_series(z, _jacobi_theta_inverse_tau_power{x}, pol) / sqrt(x.tau); + } else if (x.tau < 1.0) { // DLMF 20.7.33 + // Reduce to -Pi/2 <= z <= Pi/2 (theta4 has period Pi) + RealType dz, k; + z = x.reduce(z, 2, dz, k); + return _IMAGINARY_jacobi_theta2tau(z, dz, x, pol); + } + return RealType(1) + _jacobi_theta4m1_series(z, _jacobi_theta_tau_power{x}, pol); +} + +// = 2 * SUM (-1)^n exp(i*Pi*Tau*(n)^2) * cos(2nz), n >= 1 (theta4 minus one) +// This preserves accuracy for small values of q (i.e. tau > 1). +template +inline RealType +_jacobi_theta4m1_dispatch(RealType z, const _jacobi_theta_exponents& x, const Policy& pol) { + if (x.tau < 1.0) + return _jacobi_theta4_dispatch(z, x, pol) - RealType(1); + return _jacobi_theta4m1_series(z, _jacobi_theta_tau_power{x}, pol); +} + +// The twelve _imp functions below validate their arguments and then hand +// over to the dispatchers above. The q versions use the direct series with +// pow() when q < exp(-Pi), and otherwise go through a = -ln(q). + template inline RealType jacobi_theta1tau_imp(RealType z, RealType tau, const Policy& pol, const char *function) @@ -456,42 +775,12 @@ jacobi_theta1tau_imp(RealType z, RealType tau, const Policy& pol, const char *fu return result; if (!_jacobi_theta_check_z(z, pol, function, &result)) return result; - if (abs(z) == 0.0) return result; - if (tau < 1.0) { - // Reduce to -Pi <= z <= Pi (theta1 has period 2*Pi)... - z = fmod(z, constants::two_pi()); - while (z > constants::pi()) { - z -= constants::two_pi(); - } - while (z < -constants::pi()) { - z += constants::two_pi(); - } - // ...then to -Pi/2 <= z <= Pi/2 using theta1(z + Pi) = -theta1(z)... - RealType sign = 1; - if (z > constants::half_pi()) { - z -= constants::pi(); - sign = -sign; - } else if (z < -constants::half_pi()) { - z += constants::pi(); - sign = -sign; - } - // ...and finally to 0 <= z <= Pi/2 since theta1 is odd. - if (z < 0) { - z = -z; - sign = -sign; - } - - return sign * _IMAGINARY_jacobi_theta1tau(z, RealType(1/tau), pol); - } - - return _jacobi_theta1_series(z, _jacobi_theta_tau_power{tau}, pol); + return _jacobi_theta1_dispatch(z, _jacobi_theta_exponents::from_tau(tau), pol); } -// First Jacobi theta function (Parameterized by q) -// = 2 * SUM (-1)^n * q^(n+1/2)^2 * sin((2n+1)z) template inline RealType jacobi_theta1_imp(RealType z, RealType q, const Policy& pol, const char *function) { @@ -502,20 +791,19 @@ jacobi_theta1_imp(RealType z, RealType q, const Policy& pol, const char *functio return result; if (!_jacobi_theta_check_z(z, pol, function, &result)) return result; + if (abs(z) == 0.0) + return result; if (q < exp(-constants::pi())) return _jacobi_theta1_series(z, _jacobi_theta_q_power{q}, pol); - return jacobi_theta1tau_imp(z, RealType (-log(q)/constants::pi()), pol, function); + return _jacobi_theta1_dispatch(z, _jacobi_theta_exponents::from_nome(q), pol); } -// Second Jacobi theta function (Parameterized by tau - assumed imaginary) -// = 2 * SUM exp(i*Pi*Tau*(n+1/2)^2) * cos((2n+1)z) template inline RealType jacobi_theta2tau_imp(RealType z, RealType tau, const Policy& pol, const char *function) { - BOOST_MATH_STD_USING RealType result = 0; if (!_jacobi_theta_check_tau(tau, pol, function, &result)) @@ -523,25 +811,9 @@ jacobi_theta2tau_imp(RealType z, RealType tau, const Policy& pol, const char *fu if (!_jacobi_theta_check_z(z, pol, function, &result)) return result; - if (tau < 1.0 && abs(z) == 0.0) { - return jacobi_theta4tau(z, RealType(1/tau), pol) / sqrt(tau); - } else if (tau < 1.0) { // DLMF 20.7.31 - z = fmod(z, constants::two_pi()); - while (z > constants::pi()) { - z -= constants::two_pi(); - } - while (z < -constants::pi()) { - z += constants::two_pi(); - } - - return _IMAGINARY_jacobi_theta4tau(z, RealType(1/tau), pol); - } - - return _jacobi_theta2_series(z, _jacobi_theta_tau_power{tau}, pol); + return _jacobi_theta2_dispatch(z, _jacobi_theta_exponents::from_tau(tau), pol); } -// Second Jacobi theta function, parameterized by q -// = 2 * SUM q^(n+1/2)^2 * cos((2n+1)z) template inline RealType jacobi_theta2_imp(RealType z, RealType q, const Policy& pol, const char *function) { @@ -556,16 +828,13 @@ jacobi_theta2_imp(RealType z, RealType q, const Policy& pol, const char *functio if (q < exp(-constants::pi())) return _jacobi_theta2_series(z, _jacobi_theta_q_power{q}, pol); - return jacobi_theta2tau_imp(z, RealType (-log(q)/constants::pi()), pol, function); + return _jacobi_theta2_dispatch(z, _jacobi_theta_exponents::from_nome(q), pol); } -// Third Jacobi theta function, parameterized by tau -// = 1 + 2 * SUM exp(i*Pi*Tau*(n)^2) * cos(2nz) template inline RealType jacobi_theta3tau_imp(RealType z, RealType tau, const Policy& pol, const char *function) { - BOOST_MATH_STD_USING RealType result = 0; if (!_jacobi_theta_check_tau(tau, pol, function, &result)) @@ -573,30 +842,13 @@ jacobi_theta3tau_imp(RealType z, RealType tau, const Policy& pol, const char *fu if (!_jacobi_theta_check_z(z, pol, function, &result)) return result; - if (tau < 1.0 && abs(z) == 0.0) { - return jacobi_theta3tau(z, RealType(1/tau), pol) / sqrt(tau); - } else if (tau < 1.0) { // DLMF 20.7.32 - z = fmod(z, constants::pi()); - while (z > constants::half_pi()) { - z -= constants::pi(); - } - while (z < -constants::half_pi()) { - z += constants::pi(); - } - return _IMAGINARY_jacobi_theta3tau(z, RealType(1/tau), pol); - } - return RealType(1) + _jacobi_theta3m1_series(z, _jacobi_theta_tau_power{tau}, pol); + return _jacobi_theta3_dispatch(z, _jacobi_theta_exponents::from_tau(tau), pol); } -// Third Jacobi theta function, minus one (Parameterized by tau - assumed imaginary) -// This function preserves accuracy for small values of q (i.e. |q| < exp(-Pi) = 0.043) -// For larger values of q, the minus one version usually won't help. -// = 2 * SUM exp(i*Pi*Tau*(n)^2) * cos(2nz) template inline RealType jacobi_theta3m1tau_imp(RealType z, RealType tau, const Policy& pol, const char *function) { - BOOST_MATH_STD_USING RealType result = 0; if (!_jacobi_theta_check_tau(tau, pol, function, &result)) @@ -604,14 +856,9 @@ jacobi_theta3m1tau_imp(RealType z, RealType tau, const Policy& pol, const char * if (!_jacobi_theta_check_z(z, pol, function, &result)) return result; - if (tau < 1.0) - return jacobi_theta3tau_imp(z, tau, pol, function) - RealType(1); - - return _jacobi_theta3m1_series(z, _jacobi_theta_tau_power{tau}, pol); + return _jacobi_theta3m1_dispatch(z, _jacobi_theta_exponents::from_tau(tau), pol); } -// Third Jacobi theta function, minus one (parameterized by q) -// = 2 * SUM q^n^2 * cos(2nz) template inline RealType jacobi_theta3m1_imp(RealType z, RealType q, const Policy& pol, const char *function) { @@ -626,11 +873,9 @@ jacobi_theta3m1_imp(RealType z, RealType q, const Policy& pol, const char *funct if (q < exp(-constants::pi())) return _jacobi_theta3m1_series(z, _jacobi_theta_q_power{q}, pol); - return jacobi_theta3m1tau_imp(z, RealType (-log(q)/constants::pi()), pol, function); + return _jacobi_theta3m1_dispatch(z, _jacobi_theta_exponents::from_nome(q), pol); } -// Third Jacobi theta function (parameterized by q) -// = 1 + 2 * SUM q^n^2 * cos(2nz) template inline RealType jacobi_theta3_imp(RealType z, RealType q, const Policy& pol, const char *function) { @@ -645,16 +890,13 @@ jacobi_theta3_imp(RealType z, RealType q, const Policy& pol, const char *functio if (q < exp(-constants::pi())) return RealType(1) + _jacobi_theta3m1_series(z, _jacobi_theta_q_power{q}, pol); - return jacobi_theta3tau_imp(z, RealType (-log(q)/constants::pi()), pol, function); + return _jacobi_theta3_dispatch(z, _jacobi_theta_exponents::from_nome(q), pol); } -// Fourth Jacobi theta function (Parameterized by tau) -// = 1 + 2 * SUM (-1)^n exp(i*Pi*Tau*(n)^2) * cos(2nz) template inline RealType jacobi_theta4tau_imp(RealType z, RealType tau, const Policy& pol, const char *function) { - BOOST_MATH_STD_USING RealType result = 0; if (!_jacobi_theta_check_tau(tau, pol, function, &result)) @@ -662,30 +904,13 @@ jacobi_theta4tau_imp(RealType z, RealType tau, const Policy& pol, const char *fu if (!_jacobi_theta_check_z(z, pol, function, &result)) return result; - if (tau < 1.0 && abs(z) == 0.0) { - return jacobi_theta2tau(z, RealType(1/tau), pol) / sqrt(tau); - } else if (tau < 1.0) { // DLMF 20.7.33 - z = fmod(z, constants::pi()); - while (z > constants::half_pi()) { - z -= constants::pi(); - } - while (z < -constants::half_pi()) { - z += constants::pi(); - } - return _IMAGINARY_jacobi_theta2tau(z, RealType(1/tau), pol); - } - - return RealType(1) + _jacobi_theta4m1_series(z, _jacobi_theta_tau_power{tau}, pol); + return _jacobi_theta4_dispatch(z, _jacobi_theta_exponents::from_tau(tau), pol); } -// Fourth Jacobi theta function, minus one (Parameterized by tau) -// This function preserves accuracy for small values of q (i.e. tau > 1) -// = 2 * SUM (-1)^n exp(i*Pi*Tau*(n)^2) * cos(2nz) template inline RealType jacobi_theta4m1tau_imp(RealType z, RealType tau, const Policy& pol, const char *function) { - BOOST_MATH_STD_USING RealType result = 0; if (!_jacobi_theta_check_tau(tau, pol, function, &result)) @@ -693,15 +918,9 @@ jacobi_theta4m1tau_imp(RealType z, RealType tau, const Policy& pol, const char * if (!_jacobi_theta_check_z(z, pol, function, &result)) return result; - if (tau < 1.0) - return jacobi_theta4tau_imp(z, tau, pol, function) - RealType(1); - - return _jacobi_theta4m1_series(z, _jacobi_theta_tau_power{tau}, pol); + return _jacobi_theta4m1_dispatch(z, _jacobi_theta_exponents::from_tau(tau), pol); } -// Fourth Jacobi theta function, minus one (Parameterized by q) -// This function preserves accuracy for small values of q -// = 2 * SUM (-1)^n q^n^2 * cos(2nz) template inline RealType jacobi_theta4m1_imp(RealType z, RealType q, const Policy& pol, const char *function) { @@ -716,11 +935,9 @@ jacobi_theta4m1_imp(RealType z, RealType q, const Policy& pol, const char *funct if (q < exp(-constants::pi())) return _jacobi_theta4m1_series(z, _jacobi_theta_q_power{q}, pol); - return jacobi_theta4m1tau_imp(z, RealType (-log(q)/constants::pi()), pol, function); + return _jacobi_theta4m1_dispatch(z, _jacobi_theta_exponents::from_nome(q), pol); } -// Fourth Jacobi theta function, parameterized by q -// = 1 + 2 * SUM (-1)^n q^n^2 * cos(2nz) template inline RealType jacobi_theta4_imp(RealType z, RealType q, const Policy& pol, const char *function) { @@ -735,7 +952,7 @@ jacobi_theta4_imp(RealType z, RealType q, const Policy& pol, const char *functio if (q < exp(-constants::pi())) return RealType(1) + _jacobi_theta4m1_series(z, _jacobi_theta_q_power{q}, pol); - return jacobi_theta4tau_imp(z, RealType(-log(q)/constants::pi()), pol, function); + return _jacobi_theta4_dispatch(z, _jacobi_theta_exponents::from_nome(q), pol); } // Begin public API diff --git a/reporting/accuracy/plot_jacobi_theta_q.cpp b/reporting/accuracy/plot_jacobi_theta_q.cpp index 437fc3b593..be9d0264ba 100644 --- a/reporting/accuracy/plot_jacobi_theta_q.cpp +++ b/reporting/accuracy/plot_jacobi_theta_q.cpp @@ -47,7 +47,7 @@ int main() { PreciseReal clip = 100; std::string filename1 = "jacobi_theta1q_" + boost::core::demangle(typeid(CoarseReal).name()) + ".svg"; - auto plot1 = ulps_plot(jacobi_theta1_precise, CoarseReal(0), CoarseReal(0.999999), samples); + auto plot1 = ulps_plot(jacobi_theta1_precise, CoarseReal(0.000001), CoarseReal(0.999999), samples); plot1.clip(clip).width(width); std::string title1 = "jacobi_theta1(5.0, q) ULP plot at " + boost::core::demangle(typeid(CoarseReal).name()) + " precision"; plot1.title(title1); @@ -56,7 +56,7 @@ int main() { plot1.write(filename1); std::string filename2 = "jacobi_theta2q_" + boost::core::demangle(typeid(CoarseReal).name()) + ".svg"; - auto plot2 = ulps_plot(jacobi_theta2_precise, CoarseReal(0), CoarseReal(0.999999), samples); + auto plot2 = ulps_plot(jacobi_theta2_precise, CoarseReal(0.000001), CoarseReal(0.999999), samples); plot2.clip(clip).width(width); std::string title2 = "jacobi_theta2(0.4, q) ULP plot at " + boost::core::demangle(typeid(CoarseReal).name()) + " precision"; plot2.title(title2); @@ -65,7 +65,7 @@ int main() { plot2.write(filename2); std::string filename3 = "jacobi_theta3q_" + boost::core::demangle(typeid(CoarseReal).name()) + ".svg"; - auto plot3 = ulps_plot(jacobi_theta3_precise, CoarseReal(0), CoarseReal(0.999999), samples); + auto plot3 = ulps_plot(jacobi_theta3_precise, CoarseReal(0.000001), CoarseReal(0.999999), samples); plot3.clip(clip).width(width); std::string title3 = "jacobi_theta3(0.4, q) ULP plot at " + boost::core::demangle(typeid(CoarseReal).name()) + " precision"; plot3.title(title3); @@ -74,7 +74,7 @@ int main() { plot3.write(filename3); std::string filename4 = "jacobi_theta4q_" + boost::core::demangle(typeid(CoarseReal).name()) + ".svg"; - auto plot4 = ulps_plot(jacobi_theta4_precise, CoarseReal(0), CoarseReal(0.999999), samples); + auto plot4 = ulps_plot(jacobi_theta4_precise, CoarseReal(0.000001), CoarseReal(0.999999), samples); plot4.clip(clip).width(width); std::string title4 = "jacobi_theta4(5.0, q) ULP plot at " + boost::core::demangle(typeid(CoarseReal).name()) + " precision"; plot4.title(title4); diff --git a/test/test_jacobi_theta.cpp b/test/test_jacobi_theta.cpp index 3a803d3aea..b2430e9d20 100644 --- a/test/test_jacobi_theta.cpp +++ b/test/test_jacobi_theta.cpp @@ -167,8 +167,14 @@ BOOST_AUTO_TEST_CASE( test_main ) } for (double q=0.0078125; q<1.0; q += 0.0078125) { // = 1/128 + // The periodicity test shifts z by the rounded constant two_pi, which + // differs from the true period by about eps. For large q the theta + // functions are steep enough (their logarithmic derivative is of + // order 1/tau = -pi/ln q) that this shift changes them by more than + // the rounding of the evaluation itself, so allow for it. + double periodicity_tol = 100 * eps + 4 * constants::pi() * constants::pi() * eps / -log(q); for (double z=-8.0; z<=8.0; z += 0.125) { - test_periodicity(z, q, 100 * eps); + test_periodicity(z, q, periodicity_tol); test_argument_translation(z, q, 100 * eps); test_sums_of_squares(z, q, 100 * eps); // The addition formula is complicated, cut it some extra slack @@ -188,7 +194,8 @@ BOOST_AUTO_TEST_CASE( test_main ) test_special_values(eps); for (double s=0.125; s<3.0; s+=0.125) { - test_mellin_transforms(2.0 + s, eps, 3 * eps); + // The integrals sum thousands of theta values, so allow a few ulps + test_mellin_transforms(2.0 + s, eps, 6 * eps); test_laplace_transforms(s, eps, 4 * eps); } From 1d1ea52141274949d2c38ba87ad10811aa299253 Mon Sep 17 00:00:00 2001 From: Evan Miller Date: Thu, 17 Sep 2026 10:25:03 -0400 Subject: [PATCH 2/3] Benchmark the theta evaluation regimes; fix expression-template types Add cases to the Jacobi theta benchmark that exercise each evaluation path (the Gaussian sums for tau < 1, the z = 0 shortcut, the direct series) at representative points, for the tau and q parameterizations, in float, double and cpp_bin_float_50. The existing theta1 cases draw q from (0, 0.01), which only exercises the direct series in q. Guard the float128 cases so the benchmark builds where quadmath is unavailable. Building it with mpfr types showed that two of the new exponent helpers were called with expression templates rather than numbers, which fails to deduce; convert explicitly. Co-Authored-By: Claude Fable 5.1 --- .../math/special_functions/jacobi_theta.hpp | 4 +- .../performance/jacobi_theta_performance.cpp | 54 +++++++++++++++++++ 2 files changed, 56 insertions(+), 2 deletions(-) diff --git a/include/boost/math/special_functions/jacobi_theta.hpp b/include/boost/math/special_functions/jacobi_theta.hpp index 4512002cf3..52a4fedfcb 100644 --- a/include/boost/math/special_functions/jacobi_theta.hpp +++ b/include/boost/math/special_functions/jacobi_theta.hpp @@ -391,7 +391,7 @@ struct _jacobi_theta_exponents { RealType dc; RealType c = pi_times(k * RealType(half_pis) / 2, dc); RealType r = z - c; - dz = _jacobi_theta_sum_error(z, -c, r) - dc; + dz = _jacobi_theta_sum_error(z, RealType(-c), r) - dc; return r; } @@ -427,7 +427,7 @@ struct _jacobi_theta_exponents { dN += pi_hi * dP; RealType dQ; RealType Q = divide_by_a(N, dN, dQ); - RealType result = boost::math::expm1(-2 * Q); + RealType result = boost::math::expm1(RealType(-2 * Q)); if (result == -1) return result; return result - 2 * dQ * (1 + result); diff --git a/reporting/performance/jacobi_theta_performance.cpp b/reporting/performance/jacobi_theta_performance.cpp index a1eea6ad0f..71aecd8e3f 100644 --- a/reporting/performance/jacobi_theta_performance.cpp +++ b/reporting/performance/jacobi_theta_performance.cpp @@ -6,13 +6,17 @@ #include #include #include +#ifdef BOOST_HAS_FLOAT128 #include +#endif #include #include using boost::multiprecision::number; using boost::multiprecision::mpfr_float_backend; +#ifdef BOOST_HAS_FLOAT128 using boost::multiprecision::float128; +#endif using boost::multiprecision::cpp_bin_float_50; using boost::multiprecision::cpp_bin_float_100; using boost::math::jacobi_theta1; @@ -37,7 +41,9 @@ void JacobiTheta1(benchmark::State& state) BENCHMARK_TEMPLATE(JacobiTheta1, float); BENCHMARK_TEMPLATE(JacobiTheta1, double); BENCHMARK_TEMPLATE(JacobiTheta1, long double); +#ifdef BOOST_HAS_FLOAT128 BENCHMARK_TEMPLATE(JacobiTheta1, float128); +#endif BENCHMARK_TEMPLATE(JacobiTheta1, number>); BENCHMARK_TEMPLATE(JacobiTheta1, number>); BENCHMARK_TEMPLATE(JacobiTheta1, number>); @@ -65,7 +71,9 @@ void JacobiTheta1Tau(benchmark::State& state) BENCHMARK_TEMPLATE(JacobiTheta1Tau, float); BENCHMARK_TEMPLATE(JacobiTheta1Tau, double); BENCHMARK_TEMPLATE(JacobiTheta1Tau, long double); +#ifdef BOOST_HAS_FLOAT128 BENCHMARK_TEMPLATE(JacobiTheta1Tau, float128); +#endif BENCHMARK_TEMPLATE(JacobiTheta1Tau, number>); BENCHMARK_TEMPLATE(JacobiTheta1Tau, number>); BENCHMARK_TEMPLATE(JacobiTheta1Tau, number>); @@ -74,4 +82,50 @@ BENCHMARK_TEMPLATE(JacobiTheta1Tau, number>); BENCHMARK_TEMPLATE(JacobiTheta1Tau, cpp_bin_float_50); BENCHMARK_TEMPLATE(JacobiTheta1Tau, cpp_bin_float_100); +// The evaluation regimes, as (z, tau) pairs: the double-sided Gaussian sums +// used for tau < 1 and z != 0, the single-sided transformed series used for +// tau < 1 and z == 0, and the direct Fourier series used for tau >= 1. The +// argument index selects the pair. +static const double jacobi_theta_regimes[][2] = { + { 0.5, 0.3 }, // Gaussian sums + { 0.3, 0.05 }, // Gaussian sums, small tau + { 0.0, 0.3 }, // z == 0 shortcut + { 0.5, 3.0 }, // direct series + { 0.0, 5.0 }, // direct series, one or two terms +}; + +struct theta1tau { template Real operator()(Real z, Real tau) const { return boost::math::jacobi_theta1tau(z, tau); } }; +struct theta2tau { template Real operator()(Real z, Real tau) const { return boost::math::jacobi_theta2tau(z, tau); } }; +struct theta3tau { template Real operator()(Real z, Real tau) const { return boost::math::jacobi_theta3tau(z, tau); } }; +struct theta4tau { template Real operator()(Real z, Real tau) const { return boost::math::jacobi_theta4tau(z, tau); } }; +struct theta4m1tau { template Real operator()(Real z, Real tau) const { return boost::math::jacobi_theta4m1tau(z, tau); } }; +// theta4(z, q) with q = exp(-pi tau), i.e. the q parameterization of the same point +struct theta4q { template Real operator()(Real z, Real tau) const { return boost::math::jacobi_theta4(z, exp(-boost::math::constants::pi() * tau)); } }; + +template +void JacobiThetaRegime(benchmark::State& state) +{ + const double* regime = jacobi_theta_regimes[state.range(0)]; + Real z = static_cast(regime[0]); + Real tau = static_cast(regime[1]); + F f; + for (auto _ : state) + { + benchmark::DoNotOptimize(f(z, tau)); + tau += std::numeric_limits::epsilon(); + } +} + +#define JACOBI_THETA_REGIMES(F, Real) BENCHMARK_TEMPLATE(JacobiThetaRegime, Real, F)->DenseRange(0, 4) + +JACOBI_THETA_REGIMES(theta1tau, float); +JACOBI_THETA_REGIMES(theta1tau, double); +JACOBI_THETA_REGIMES(theta2tau, double); +JACOBI_THETA_REGIMES(theta3tau, double); +JACOBI_THETA_REGIMES(theta4tau, float); +JACOBI_THETA_REGIMES(theta4tau, double); +JACOBI_THETA_REGIMES(theta4m1tau, double); +JACOBI_THETA_REGIMES(theta4q, double); +JACOBI_THETA_REGIMES(theta4tau, cpp_bin_float_50); + BENCHMARK_MAIN(); From b8a6a4261f5559624a1feb458011b58d886d43a0 Mon Sep 17 00:00:00 2001 From: Evan Miller Date: Thu, 17 Sep 2026 13:37:43 -0400 Subject: [PATCH 3/3] Move the theta helpers into detail and always use std::fma for built-in types Review feedback: the underscore-prefixed helper functions and structs now live in boost::math::detail without the prefix, matching the rest of the library, and the float/double/long double overloads of the product-error helper call std::fma unconditionally instead of only when FP_FAST_FMA* is defined. std::fma is correctly rounded either way; the Dekker splitting remains as the generic template for types without an fma. Co-Authored-By: Claude Fable 5.1 --- .../math/special_functions/jacobi_theta.hpp | 256 +++++++++--------- 1 file changed, 127 insertions(+), 129 deletions(-) diff --git a/include/boost/math/special_functions/jacobi_theta.hpp b/include/boost/math/special_functions/jacobi_theta.hpp index 52a4fedfcb..356cf85c8e 100644 --- a/include/boost/math/special_functions/jacobi_theta.hpp +++ b/include/boost/math/special_functions/jacobi_theta.hpp @@ -173,19 +173,21 @@ inline typename tools::promote_args::type jacobi_theta3m1tau(T z, U tau, c BOOST_MATH_EXPORT template inline typename tools::promote_args::type jacobi_theta4m1tau(T z, U tau, const Policy& pol); +namespace detail { + // Compare the non-oscillating component of the delta to the previous delta. // Both are assumed to be non-negative. Written so that a NaN delta counts as // converged: otherwise a NaN would never satisfy the test and the summation // loops below would never terminate. template inline bool -_jacobi_theta_converged(RealType last_delta, RealType delta, RealType eps) { +jacobi_theta_converged(RealType last_delta, RealType delta, RealType eps) { return !(delta > eps*last_delta); } template inline bool -_jacobi_theta_check_z(RealType z, const Policy& pol, const char* function, RealType* result) { +jacobi_theta_check_z(RealType z, const Policy& pol, const char* function, RealType* result) { if (!(boost::math::isfinite)(z)) { *result = policies::raise_domain_error(function, "z must be finite but got %1%.", z, pol); return false; @@ -195,7 +197,7 @@ _jacobi_theta_check_z(RealType z, const Policy& pol, const char* function, RealT template inline bool -_jacobi_theta_check_tau(RealType tau, const Policy& pol, const char* function, RealType* result) { +jacobi_theta_check_tau(RealType tau, const Policy& pol, const char* function, RealType* result) { // The negated comparison also rejects NaN. if (!(tau > 0)) { *result = policies::raise_domain_error(function, "tau must be greater than 0 but got %1%.", tau, pol); @@ -206,7 +208,7 @@ _jacobi_theta_check_tau(RealType tau, const Policy& pol, const char* function, R template inline bool -_jacobi_theta_check_q(RealType q, const Policy& pol, const char* function, RealType* result) { +jacobi_theta_check_q(RealType q, const Policy& pol, const char* function, RealType* result) { // The negated comparison also rejects NaN. if (!(q > 0 && q < 1)) { *result = policies::raise_domain_error(function, "q must be greater than 0 and less than 1 but got %1%.", q, pol); @@ -221,7 +223,7 @@ _jacobi_theta_check_q(RealType q, const Policy& pol, const char* function, RealT // exp() would multiply the rounding error by |log q|, which is exactly the // small-q regime the "m1" functions exist to serve. template -struct _jacobi_theta_q_power { +struct jacobi_theta_q_power { RealType q; RealType operator()(RealType exponent) const { BOOST_MATH_STD_USING @@ -241,7 +243,7 @@ struct _jacobi_theta_q_power { // cached value may split at the wrong place, in which case the corrections // below are merely inexact, i.e. no worse than not applying them.) template -inline RealType _jacobi_theta_splitter() { +inline RealType jacobi_theta_splitter() { BOOST_MATH_STD_USING static const RealType splitter = ldexp(RealType(1), (tools::digits() + 1) / 2) + 1; return splitter; @@ -250,41 +252,35 @@ inline RealType _jacobi_theta_splitter() { // Veltkamp's splitting: t = hi + lo where hi holds only the leading half of // the significand, so that products of hi parts are exact. template -inline void _jacobi_theta_split(RealType t, RealType splitter, RealType& hi, RealType& lo) { +inline void jacobi_theta_split(RealType t, RealType splitter, RealType& hi, RealType& lo) { RealType g = splitter * t; hi = g - (g - t); lo = t - hi; } -// The exact product a*b is ab + (return value), where ab = fl(a*b). Where -// the standard library advertises a fast fused multiply-add it is used -// instead, since fma(a, b, -ab) is exactly this quantity. +// The exact product a*b is ab + (return value), where ab = fl(a*b). For the +// built-in types std::fma(a, b, -ab) is exactly this quantity, so it is used +// instead; the splitting below is for types without an fma. template -inline RealType _jacobi_theta_product_error(RealType a, RealType b, RealType ab, RealType splitter) { +inline RealType jacobi_theta_product_error(RealType a, RealType b, RealType ab, RealType splitter) { RealType ah, al, bh, bl; - _jacobi_theta_split(a, splitter, ah, al); - _jacobi_theta_split(b, splitter, bh, bl); + jacobi_theta_split(a, splitter, ah, al); + jacobi_theta_split(b, splitter, bh, bl); return (((ah * bh - ab) + ah * bl) + al * bh) + al * bl; } -#ifdef FP_FAST_FMAF -inline float _jacobi_theta_product_error(float a, float b, float ab, float) { +inline float jacobi_theta_product_error(float a, float b, float ab, float) { return std::fma(a, b, -ab); } -#endif -#ifdef FP_FAST_FMA -inline double _jacobi_theta_product_error(double a, double b, double ab, double) { +inline double jacobi_theta_product_error(double a, double b, double ab, double) { return std::fma(a, b, -ab); } -#endif -#ifdef FP_FAST_FMAL -inline long double _jacobi_theta_product_error(long double a, long double b, long double ab, long double) { +inline long double jacobi_theta_product_error(long double a, long double b, long double ab, long double) { return std::fma(a, b, -ab); } -#endif // The exact sum a+b is s + (return value), where s = fl(a+b). template -inline RealType _jacobi_theta_sum_error(RealType a, RealType b, RealType s) { +inline RealType jacobi_theta_sum_error(RealType a, RealType b, RealType s) { RealType bb = s - a; return (a - (s - bb)) + (b - bb); } @@ -293,7 +289,7 @@ inline RealType _jacobi_theta_sum_error(RealType a, RealType b, RealType s) { // 120 bits, as five pieces of 24 significant bits each. Types with more // precision than that use their own rounding of pi with lo = 0. template -inline RealType _jacobi_theta_pi(RealType& lo) { +inline RealType jacobi_theta_pi(RealType& lo) { BOOST_MATH_STD_USING if (tools::digits() > 116) { lo = 0; @@ -310,24 +306,24 @@ inline RealType _jacobi_theta_pi(RealType& lo) { // carries the rounding of the logarithm, which is inherent to that // parameterization, but nothing else. template -struct _jacobi_theta_exponents { +struct jacobi_theta_exponents { RealType tau; // only used for the scale factor 1/sqrt(tau) RealType splitter; RealType pi_hi, pi_lo; // pi = pi_hi + pi_lo RealType a_hi, a_lo; // pi * tau = -ln(q) = a_hi + a_lo RealType inv_a_hi; // 1 / a_hi, only needed when tau < 1 - static _jacobi_theta_exponents from_tau(RealType tau) { - _jacobi_theta_exponents x; + static jacobi_theta_exponents from_tau(RealType tau) { + jacobi_theta_exponents x; x.tau = tau; x.a_hi = x.pi_times(tau, x.a_lo); x.inv_a_hi = (tau < 1) ? 1 / x.a_hi : RealType(0); return x; } - static _jacobi_theta_exponents from_nome(RealType q) { + static jacobi_theta_exponents from_nome(RealType q) { BOOST_MATH_STD_USING - _jacobi_theta_exponents x; + jacobi_theta_exponents x; x.a_hi = -log(q); x.a_lo = 0; x.tau = x.a_hi / constants::pi(); @@ -341,7 +337,7 @@ struct _jacobi_theta_exponents { RealType q = pi_lo * e; RealType hi = p + q; // |q| << |p|, so the sum's error is simply q - (hi - p) - lo = (q - (hi - p)) + _jacobi_theta_product_error(pi_hi, e, p, splitter); + lo = (q - (hi - p)) + jacobi_theta_product_error(pi_hi, e, p, splitter); return hi; } @@ -349,7 +345,7 @@ struct _jacobi_theta_exponents { RealType divide_by_a(RealType N, RealType dN, RealType& dE) const { RealType E = N * inv_a_hi; RealType Ea = E * a_hi; - RealType r = (N - Ea) - _jacobi_theta_product_error(E, a_hi, Ea, splitter); + RealType r = (N - Ea) - jacobi_theta_product_error(E, a_hi, Ea, splitter); dE = (r + dN - E * a_lo) * inv_a_hi; return E; } @@ -361,7 +357,7 @@ struct _jacobi_theta_exponents { RealType result = exp(-E); if (result == 0) return result; - RealType dE = _jacobi_theta_product_error(a_hi, e, E, splitter) + a_lo * e; + RealType dE = jacobi_theta_product_error(a_hi, e, E, splitter) + a_lo * e; return result * (1 - dE); } @@ -391,7 +387,7 @@ struct _jacobi_theta_exponents { RealType dc; RealType c = pi_times(k * RealType(half_pis) / 2, dc); RealType r = z - c; - dz = _jacobi_theta_sum_error(z, RealType(-c), r) - dc; + dz = jacobi_theta_sum_error(z, RealType(-c), r) - dc; return r; } @@ -404,10 +400,10 @@ struct _jacobi_theta_exponents { RealType dc; RealType c = pi_times(RealType(m), dc) / 2; z_n = z + c; - dz += _jacobi_theta_sum_error(z, c, z_n) + dc / 2; + dz += jacobi_theta_sum_error(z, c, z_n) + dc / 2; } RealType s = z_n * z_n; - RealType ds = _jacobi_theta_product_error(z_n, z_n, s, splitter) + 2 * z_n * dz; + RealType ds = jacobi_theta_product_error(z_n, z_n, s, splitter) + 2 * z_n * dz; RealType dE; RealType E = divide_by_a(s, ds, dE); RealType result = exp(-E); @@ -421,7 +417,7 @@ struct _jacobi_theta_exponents { RealType expm1_scaled(RealType z, RealType dz, RealType k) const { BOOST_MATH_STD_USING RealType P = z * k; - RealType dP = _jacobi_theta_product_error(z, k, P, splitter) + dz * k; + RealType dP = jacobi_theta_product_error(z, k, P, splitter) + dz * k; RealType dN; RealType N = pi_times(P, dN); dN += pi_hi * dP; @@ -434,14 +430,14 @@ struct _jacobi_theta_exponents { } private: - _jacobi_theta_exponents() : splitter(_jacobi_theta_splitter()) { - pi_hi = _jacobi_theta_pi(pi_lo); + jacobi_theta_exponents() : splitter(jacobi_theta_splitter()) { + pi_hi = jacobi_theta_pi(pi_lo); } }; template -struct _jacobi_theta_tau_power { - const _jacobi_theta_exponents& x; +struct jacobi_theta_tau_power { + const jacobi_theta_exponents& x; RealType operator()(RealType exponent) const { return x.direct(exponent); } @@ -450,8 +446,8 @@ struct _jacobi_theta_tau_power { // exp(-pi * e / tau): the nome of 1/tau, used by the modular transformation // at z = 0 without forming the rounded reciprocal. template -struct _jacobi_theta_inverse_tau_power { - const _jacobi_theta_exponents& x; +struct jacobi_theta_inverse_tau_power { + const jacobi_theta_exponents& x; RealType operator()(RealType exponent) const { return x.inverse(exponent); } @@ -464,7 +460,7 @@ struct _jacobi_theta_inverse_tau_power { // = 2 * SUM (-1)^n * q^(n+1/2)^2 * sin((2n+1)z) template inline RealType -_jacobi_theta1_series(RealType z, const NomePower& q_pow, const Policy&) { +jacobi_theta1_series(RealType z, const NomePower& q_pow, const Policy&) { BOOST_MATH_STD_USING unsigned n = 0; RealType eps = policies::get_epsilon(); @@ -479,7 +475,7 @@ _jacobi_theta1_series(RealType z, const NomePower& q_pow, const Policy&) { result += delta + delta; n++; - } while (!_jacobi_theta_converged(last_q_n, q_n, eps)); + } while (!jacobi_theta_converged(last_q_n, q_n, eps)); return result; } @@ -487,7 +483,7 @@ _jacobi_theta1_series(RealType z, const NomePower& q_pow, const Policy&) { // = 2 * SUM q^(n+1/2)^2 * cos((2n+1)z) template inline RealType -_jacobi_theta2_series(RealType z, const NomePower& q_pow, const Policy&) { +jacobi_theta2_series(RealType z, const NomePower& q_pow, const Policy&) { BOOST_MATH_STD_USING unsigned n = 0; RealType eps = policies::get_epsilon(); @@ -499,7 +495,7 @@ _jacobi_theta2_series(RealType z, const NomePower& q_pow, const Policy&) { delta = q_n * cos(RealType(2*n+1)*z); result += delta + delta; n++; - } while (!_jacobi_theta_converged(last_q_n, q_n, eps)); + } while (!jacobi_theta_converged(last_q_n, q_n, eps)); return result; } @@ -507,7 +503,7 @@ _jacobi_theta2_series(RealType z, const NomePower& q_pow, const Policy&) { // = 2 * SUM q^n^2 * cos(2nz), n >= 1 (i.e. theta3 minus one) template inline RealType -_jacobi_theta3m1_series(RealType z, const NomePower& q_pow, const Policy&) { +jacobi_theta3m1_series(RealType z, const NomePower& q_pow, const Policy&) { BOOST_MATH_STD_USING unsigned n = 1; RealType eps = policies::get_epsilon(); @@ -519,7 +515,7 @@ _jacobi_theta3m1_series(RealType z, const NomePower& q_pow, const Policy&) { delta = q_n * cos(RealType(2*n)*z); result += delta + delta; n++; - } while (!_jacobi_theta_converged(last_q_n, q_n, eps)); + } while (!jacobi_theta_converged(last_q_n, q_n, eps)); return result; } @@ -527,7 +523,7 @@ _jacobi_theta3m1_series(RealType z, const NomePower& q_pow, const Policy&) { // = 2 * SUM (-1)^n q^n^2 * cos(2nz), n >= 1 (i.e. theta4 minus one) template inline RealType -_jacobi_theta4m1_series(RealType z, const NomePower& q_pow, const Policy&) { +jacobi_theta4m1_series(RealType z, const NomePower& q_pow, const Policy&) { BOOST_MATH_STD_USING unsigned n = 1; RealType eps = policies::get_epsilon(); @@ -542,7 +538,7 @@ _jacobi_theta4m1_series(RealType z, const NomePower& q_pow, const Policy&) { result += delta + delta; n++; - } while (!_jacobi_theta_converged(last_q_n, q_n, eps)); + } while (!jacobi_theta_converged(last_q_n, q_n, eps)); return result; } @@ -551,7 +547,7 @@ _jacobi_theta4m1_series(RealType z, const NomePower& q_pow, const Policy&) { // terms are negligible. template inline RealType -_jacobi_theta_sum(const _jacobi_theta_exponents& x, RealType z, RealType dz, int m, int m_step, RealType eps) { +jacobi_theta_sum(const jacobi_theta_exponents& x, RealType z, RealType dz, int m, int m_step, RealType eps) { RealType delta = 0, partial_result = 0; RealType last_delta = 0; @@ -560,12 +556,12 @@ _jacobi_theta_sum(const _jacobi_theta_exponents& x, RealType z, RealTy delta = x.gaussian(z, dz, m); partial_result += delta; m += m_step; - } while (!_jacobi_theta_converged(last_delta, delta, eps)); + } while (!jacobi_theta_converged(last_delta, delta, eps)); return partial_result; } -// The following _IMAGINARY theta functions assume imaginary z and are for +// The following imaginary theta functions assume imaginary z and are for // internal use only. They are designed to increase accuracy and reduce the // number of iterations required for convergence for large |q|. The z argument // is scaled by 1/tau, and the summations are rewritten to be double-sided @@ -590,7 +586,7 @@ _jacobi_theta_sum(const _jacobi_theta_exponents& x, RealType z, RealTy // passes the error dz of the reduced z. template inline RealType -_IMAGINARY_jacobi_theta1tau(RealType z, RealType dz, const _jacobi_theta_exponents& x, const Policy&) { +imaginary_jacobi_theta1tau(RealType z, RealType dz, const jacobi_theta_exponents& x, const Policy&) { BOOST_MATH_STD_USING RealType eps = policies::get_epsilon(); RealType result = 0, g = 0, last_g, pair; @@ -605,29 +601,29 @@ _IMAGINARY_jacobi_theta1tau(RealType z, RealType dz, const _jacobi_theta_exponen result += pair; n++; - } while (!_jacobi_theta_converged(last_g, g, eps)); + } while (!jacobi_theta_converged(last_g, g, eps)); return result / sqrt(x.tau); } template inline RealType -_IMAGINARY_jacobi_theta2tau(RealType z, RealType dz, const _jacobi_theta_exponents& x, const Policy&) { +imaginary_jacobi_theta2tau(RealType z, RealType dz, const jacobi_theta_exponents& x, const Policy&) { BOOST_MATH_STD_USING RealType eps = policies::get_epsilon(); RealType result = RealType(0); // n>=0: centers at z + Pi/2 + n*Pi - result += _jacobi_theta_sum(x, z, dz, 1, 2, eps); + result += jacobi_theta_sum(x, z, dz, 1, 2, eps); // n<0 - result += _jacobi_theta_sum(x, z, dz, -1, -2, eps); + result += jacobi_theta_sum(x, z, dz, -1, -2, eps); return result / sqrt(x.tau); } template inline RealType -_IMAGINARY_jacobi_theta3tau(RealType z, RealType dz, const _jacobi_theta_exponents& x, const Policy&) { +imaginary_jacobi_theta3tau(RealType z, RealType dz, const jacobi_theta_exponents& x, const Policy&) { BOOST_MATH_STD_USING RealType eps = policies::get_epsilon(); RealType result = 0; @@ -635,16 +631,16 @@ _IMAGINARY_jacobi_theta3tau(RealType z, RealType dz, const _jacobi_theta_exponen // n=0 result += x.gaussian(z, dz, 0); // n>0: centers at z + n*Pi - result += _jacobi_theta_sum(x, z, dz, 2, 2, eps); + result += jacobi_theta_sum(x, z, dz, 2, 2, eps); // n<0 - result += _jacobi_theta_sum(x, z, dz, -2, -2, eps); + result += jacobi_theta_sum(x, z, dz, -2, -2, eps); return result / sqrt(x.tau); } template inline RealType -_IMAGINARY_jacobi_theta4tau(RealType z, RealType dz, const _jacobi_theta_exponents& x, const Policy&) { +imaginary_jacobi_theta4tau(RealType z, RealType dz, const jacobi_theta_exponents& x, const Policy&) { BOOST_MATH_STD_USING RealType eps = policies::get_epsilon(); RealType result = 0; @@ -653,13 +649,13 @@ _IMAGINARY_jacobi_theta4tau(RealType z, RealType dz, const _jacobi_theta_exponen result += x.gaussian(z, dz, 0); // n > 0 odd: centers at z + Pi + 2n*Pi - result -= _jacobi_theta_sum(x, z, dz, 2, 4, eps); + result -= jacobi_theta_sum(x, z, dz, 2, 4, eps); // n < 0 odd - result -= _jacobi_theta_sum(x, z, dz, -2, -4, eps); + result -= jacobi_theta_sum(x, z, dz, -2, -4, eps); // n > 0 even: centers at z + 2*Pi + 2n*Pi - result += _jacobi_theta_sum(x, z, dz, 4, 4, eps); + result += jacobi_theta_sum(x, z, dz, 4, 4, eps); // n < 0 even - result += _jacobi_theta_sum(x, z, dz, -4, -4, eps); + result += jacobi_theta_sum(x, z, dz, -4, -4, eps); return result / sqrt(x.tau); } @@ -673,7 +669,7 @@ _IMAGINARY_jacobi_theta4tau(RealType z, RealType dz, const _jacobi_theta_exponen // = 2 * SUM (-1)^n * exp(i*Pi*Tau*(n+1/2)^2) * sin((2n+1)z) template inline RealType -_jacobi_theta1_dispatch(RealType z, const _jacobi_theta_exponents& x, const Policy& pol) { +jacobi_theta1_dispatch(RealType z, const jacobi_theta_exponents& x, const Policy& pol) { BOOST_MATH_STD_USING if (x.tau < 1.0) { // Reduce to -Pi/2 <= z <= Pi/2 using theta1(z + Pi) = -theta1(z)... @@ -686,41 +682,41 @@ _jacobi_theta1_dispatch(RealType z, const _jacobi_theta_exponents& x, dz = -dz; sign = -sign; } - return sign * _IMAGINARY_jacobi_theta1tau(z, dz, x, pol); + return sign * imaginary_jacobi_theta1tau(z, dz, x, pol); } - return _jacobi_theta1_series(z, _jacobi_theta_tau_power{x}, pol); + return jacobi_theta1_series(z, jacobi_theta_tau_power{x}, pol); } // = 2 * SUM exp(i*Pi*Tau*(n+1/2)^2) * cos((2n+1)z) template inline RealType -_jacobi_theta2_dispatch(RealType z, const _jacobi_theta_exponents& x, const Policy& pol) { +jacobi_theta2_dispatch(RealType z, const jacobi_theta_exponents& x, const Policy& pol) { BOOST_MATH_STD_USING if (x.tau < 1.0 && abs(z) == 0.0) { // theta4(0|1/tau)/sqrt(tau) - return (RealType(1) + _jacobi_theta4m1_series(z, _jacobi_theta_inverse_tau_power{x}, pol)) / sqrt(x.tau); + return (RealType(1) + jacobi_theta4m1_series(z, jacobi_theta_inverse_tau_power{x}, pol)) / sqrt(x.tau); } else if (x.tau < 1.0) { // DLMF 20.7.31 // Reduce to -Pi <= z <= Pi (theta2 has period 2*Pi) RealType dz, k; z = x.reduce(z, 4, dz, k); - return _IMAGINARY_jacobi_theta4tau(z, dz, x, pol); + return imaginary_jacobi_theta4tau(z, dz, x, pol); } - return _jacobi_theta2_series(z, _jacobi_theta_tau_power{x}, pol); + return jacobi_theta2_series(z, jacobi_theta_tau_power{x}, pol); } // = 1 + 2 * SUM exp(i*Pi*Tau*(n)^2) * cos(2nz) template inline RealType -_jacobi_theta3_dispatch(RealType z, const _jacobi_theta_exponents& x, const Policy& pol) { +jacobi_theta3_dispatch(RealType z, const jacobi_theta_exponents& x, const Policy& pol) { BOOST_MATH_STD_USING if (x.tau < 1.0 && abs(z) == 0.0) { // theta3(0|1/tau)/sqrt(tau) - return (RealType(1) + _jacobi_theta3m1_series(z, _jacobi_theta_inverse_tau_power{x}, pol)) / sqrt(x.tau); + return (RealType(1) + jacobi_theta3m1_series(z, jacobi_theta_inverse_tau_power{x}, pol)) / sqrt(x.tau); } else if (x.tau < 1.0) { // DLMF 20.7.32 // Reduce to -Pi/2 <= z <= Pi/2 (theta3 has period Pi) RealType dz, k; z = x.reduce(z, 2, dz, k); - return _IMAGINARY_jacobi_theta3tau(z, dz, x, pol); + return imaginary_jacobi_theta3tau(z, dz, x, pol); } - return RealType(1) + _jacobi_theta3m1_series(z, _jacobi_theta_tau_power{x}, pol); + return RealType(1) + jacobi_theta3m1_series(z, jacobi_theta_tau_power{x}, pol); } // = 2 * SUM exp(i*Pi*Tau*(n)^2) * cos(2nz), n >= 1 (theta3 minus one) @@ -728,38 +724,40 @@ _jacobi_theta3_dispatch(RealType z, const _jacobi_theta_exponents& x, // values of q, the minus one version usually won't help. template inline RealType -_jacobi_theta3m1_dispatch(RealType z, const _jacobi_theta_exponents& x, const Policy& pol) { +jacobi_theta3m1_dispatch(RealType z, const jacobi_theta_exponents& x, const Policy& pol) { if (x.tau < 1.0) - return _jacobi_theta3_dispatch(z, x, pol) - RealType(1); - return _jacobi_theta3m1_series(z, _jacobi_theta_tau_power{x}, pol); + return jacobi_theta3_dispatch(z, x, pol) - RealType(1); + return jacobi_theta3m1_series(z, jacobi_theta_tau_power{x}, pol); } // = 1 + 2 * SUM (-1)^n exp(i*Pi*Tau*(n)^2) * cos(2nz) template inline RealType -_jacobi_theta4_dispatch(RealType z, const _jacobi_theta_exponents& x, const Policy& pol) { +jacobi_theta4_dispatch(RealType z, const jacobi_theta_exponents& x, const Policy& pol) { BOOST_MATH_STD_USING if (x.tau < 1.0 && abs(z) == 0.0) { // theta2(0|1/tau)/sqrt(tau) - return _jacobi_theta2_series(z, _jacobi_theta_inverse_tau_power{x}, pol) / sqrt(x.tau); + return jacobi_theta2_series(z, jacobi_theta_inverse_tau_power{x}, pol) / sqrt(x.tau); } else if (x.tau < 1.0) { // DLMF 20.7.33 // Reduce to -Pi/2 <= z <= Pi/2 (theta4 has period Pi) RealType dz, k; z = x.reduce(z, 2, dz, k); - return _IMAGINARY_jacobi_theta2tau(z, dz, x, pol); + return imaginary_jacobi_theta2tau(z, dz, x, pol); } - return RealType(1) + _jacobi_theta4m1_series(z, _jacobi_theta_tau_power{x}, pol); + return RealType(1) + jacobi_theta4m1_series(z, jacobi_theta_tau_power{x}, pol); } // = 2 * SUM (-1)^n exp(i*Pi*Tau*(n)^2) * cos(2nz), n >= 1 (theta4 minus one) // This preserves accuracy for small values of q (i.e. tau > 1). template inline RealType -_jacobi_theta4m1_dispatch(RealType z, const _jacobi_theta_exponents& x, const Policy& pol) { +jacobi_theta4m1_dispatch(RealType z, const jacobi_theta_exponents& x, const Policy& pol) { if (x.tau < 1.0) - return _jacobi_theta4_dispatch(z, x, pol) - RealType(1); - return _jacobi_theta4m1_series(z, _jacobi_theta_tau_power{x}, pol); + return jacobi_theta4_dispatch(z, x, pol) - RealType(1); + return jacobi_theta4m1_series(z, jacobi_theta_tau_power{x}, pol); } +} // namespace detail + // The twelve _imp functions below validate their arguments and then hand // over to the dispatchers above. The q versions use the direct series with // pow() when q < exp(-Pi), and otherwise go through a = -ln(q). @@ -771,14 +769,14 @@ jacobi_theta1tau_imp(RealType z, RealType tau, const Policy& pol, const char *fu BOOST_MATH_STD_USING RealType result = 0; - if (!_jacobi_theta_check_tau(tau, pol, function, &result)) + if (!detail::jacobi_theta_check_tau(tau, pol, function, &result)) return result; - if (!_jacobi_theta_check_z(z, pol, function, &result)) + if (!detail::jacobi_theta_check_z(z, pol, function, &result)) return result; if (abs(z) == 0.0) return result; - return _jacobi_theta1_dispatch(z, _jacobi_theta_exponents::from_tau(tau), pol); + return detail::jacobi_theta1_dispatch(z, detail::jacobi_theta_exponents::from_tau(tau), pol); } template @@ -787,17 +785,17 @@ jacobi_theta1_imp(RealType z, RealType q, const Policy& pol, const char *functio BOOST_MATH_STD_USING RealType result = 0; - if (!_jacobi_theta_check_q(q, pol, function, &result)) + if (!detail::jacobi_theta_check_q(q, pol, function, &result)) return result; - if (!_jacobi_theta_check_z(z, pol, function, &result)) + if (!detail::jacobi_theta_check_z(z, pol, function, &result)) return result; if (abs(z) == 0.0) return result; if (q < exp(-constants::pi())) - return _jacobi_theta1_series(z, _jacobi_theta_q_power{q}, pol); + return detail::jacobi_theta1_series(z, detail::jacobi_theta_q_power{q}, pol); - return _jacobi_theta1_dispatch(z, _jacobi_theta_exponents::from_nome(q), pol); + return detail::jacobi_theta1_dispatch(z, detail::jacobi_theta_exponents::from_nome(q), pol); } template @@ -806,12 +804,12 @@ jacobi_theta2tau_imp(RealType z, RealType tau, const Policy& pol, const char *fu { RealType result = 0; - if (!_jacobi_theta_check_tau(tau, pol, function, &result)) + if (!detail::jacobi_theta_check_tau(tau, pol, function, &result)) return result; - if (!_jacobi_theta_check_z(z, pol, function, &result)) + if (!detail::jacobi_theta_check_z(z, pol, function, &result)) return result; - return _jacobi_theta2_dispatch(z, _jacobi_theta_exponents::from_tau(tau), pol); + return detail::jacobi_theta2_dispatch(z, detail::jacobi_theta_exponents::from_tau(tau), pol); } template @@ -820,15 +818,15 @@ jacobi_theta2_imp(RealType z, RealType q, const Policy& pol, const char *functio BOOST_MATH_STD_USING RealType result = 0; - if (!_jacobi_theta_check_q(q, pol, function, &result)) + if (!detail::jacobi_theta_check_q(q, pol, function, &result)) return result; - if (!_jacobi_theta_check_z(z, pol, function, &result)) + if (!detail::jacobi_theta_check_z(z, pol, function, &result)) return result; if (q < exp(-constants::pi())) - return _jacobi_theta2_series(z, _jacobi_theta_q_power{q}, pol); + return detail::jacobi_theta2_series(z, detail::jacobi_theta_q_power{q}, pol); - return _jacobi_theta2_dispatch(z, _jacobi_theta_exponents::from_nome(q), pol); + return detail::jacobi_theta2_dispatch(z, detail::jacobi_theta_exponents::from_nome(q), pol); } template @@ -837,12 +835,12 @@ jacobi_theta3tau_imp(RealType z, RealType tau, const Policy& pol, const char *fu { RealType result = 0; - if (!_jacobi_theta_check_tau(tau, pol, function, &result)) + if (!detail::jacobi_theta_check_tau(tau, pol, function, &result)) return result; - if (!_jacobi_theta_check_z(z, pol, function, &result)) + if (!detail::jacobi_theta_check_z(z, pol, function, &result)) return result; - return _jacobi_theta3_dispatch(z, _jacobi_theta_exponents::from_tau(tau), pol); + return detail::jacobi_theta3_dispatch(z, detail::jacobi_theta_exponents::from_tau(tau), pol); } template @@ -851,12 +849,12 @@ jacobi_theta3m1tau_imp(RealType z, RealType tau, const Policy& pol, const char * { RealType result = 0; - if (!_jacobi_theta_check_tau(tau, pol, function, &result)) + if (!detail::jacobi_theta_check_tau(tau, pol, function, &result)) return result; - if (!_jacobi_theta_check_z(z, pol, function, &result)) + if (!detail::jacobi_theta_check_z(z, pol, function, &result)) return result; - return _jacobi_theta3m1_dispatch(z, _jacobi_theta_exponents::from_tau(tau), pol); + return detail::jacobi_theta3m1_dispatch(z, detail::jacobi_theta_exponents::from_tau(tau), pol); } template @@ -865,15 +863,15 @@ jacobi_theta3m1_imp(RealType z, RealType q, const Policy& pol, const char *funct BOOST_MATH_STD_USING RealType result = 0; - if (!_jacobi_theta_check_q(q, pol, function, &result)) + if (!detail::jacobi_theta_check_q(q, pol, function, &result)) return result; - if (!_jacobi_theta_check_z(z, pol, function, &result)) + if (!detail::jacobi_theta_check_z(z, pol, function, &result)) return result; if (q < exp(-constants::pi())) - return _jacobi_theta3m1_series(z, _jacobi_theta_q_power{q}, pol); + return detail::jacobi_theta3m1_series(z, detail::jacobi_theta_q_power{q}, pol); - return _jacobi_theta3m1_dispatch(z, _jacobi_theta_exponents::from_nome(q), pol); + return detail::jacobi_theta3m1_dispatch(z, detail::jacobi_theta_exponents::from_nome(q), pol); } template @@ -882,15 +880,15 @@ jacobi_theta3_imp(RealType z, RealType q, const Policy& pol, const char *functio BOOST_MATH_STD_USING RealType result = 0; - if (!_jacobi_theta_check_q(q, pol, function, &result)) + if (!detail::jacobi_theta_check_q(q, pol, function, &result)) return result; - if (!_jacobi_theta_check_z(z, pol, function, &result)) + if (!detail::jacobi_theta_check_z(z, pol, function, &result)) return result; if (q < exp(-constants::pi())) - return RealType(1) + _jacobi_theta3m1_series(z, _jacobi_theta_q_power{q}, pol); + return RealType(1) + detail::jacobi_theta3m1_series(z, detail::jacobi_theta_q_power{q}, pol); - return _jacobi_theta3_dispatch(z, _jacobi_theta_exponents::from_nome(q), pol); + return detail::jacobi_theta3_dispatch(z, detail::jacobi_theta_exponents::from_nome(q), pol); } template @@ -899,12 +897,12 @@ jacobi_theta4tau_imp(RealType z, RealType tau, const Policy& pol, const char *fu { RealType result = 0; - if (!_jacobi_theta_check_tau(tau, pol, function, &result)) + if (!detail::jacobi_theta_check_tau(tau, pol, function, &result)) return result; - if (!_jacobi_theta_check_z(z, pol, function, &result)) + if (!detail::jacobi_theta_check_z(z, pol, function, &result)) return result; - return _jacobi_theta4_dispatch(z, _jacobi_theta_exponents::from_tau(tau), pol); + return detail::jacobi_theta4_dispatch(z, detail::jacobi_theta_exponents::from_tau(tau), pol); } template @@ -913,12 +911,12 @@ jacobi_theta4m1tau_imp(RealType z, RealType tau, const Policy& pol, const char * { RealType result = 0; - if (!_jacobi_theta_check_tau(tau, pol, function, &result)) + if (!detail::jacobi_theta_check_tau(tau, pol, function, &result)) return result; - if (!_jacobi_theta_check_z(z, pol, function, &result)) + if (!detail::jacobi_theta_check_z(z, pol, function, &result)) return result; - return _jacobi_theta4m1_dispatch(z, _jacobi_theta_exponents::from_tau(tau), pol); + return detail::jacobi_theta4m1_dispatch(z, detail::jacobi_theta_exponents::from_tau(tau), pol); } template @@ -927,15 +925,15 @@ jacobi_theta4m1_imp(RealType z, RealType q, const Policy& pol, const char *funct BOOST_MATH_STD_USING RealType result = 0; - if (!_jacobi_theta_check_q(q, pol, function, &result)) + if (!detail::jacobi_theta_check_q(q, pol, function, &result)) return result; - if (!_jacobi_theta_check_z(z, pol, function, &result)) + if (!detail::jacobi_theta_check_z(z, pol, function, &result)) return result; if (q < exp(-constants::pi())) - return _jacobi_theta4m1_series(z, _jacobi_theta_q_power{q}, pol); + return detail::jacobi_theta4m1_series(z, detail::jacobi_theta_q_power{q}, pol); - return _jacobi_theta4m1_dispatch(z, _jacobi_theta_exponents::from_nome(q), pol); + return detail::jacobi_theta4m1_dispatch(z, detail::jacobi_theta_exponents::from_nome(q), pol); } template @@ -944,15 +942,15 @@ jacobi_theta4_imp(RealType z, RealType q, const Policy& pol, const char *functio BOOST_MATH_STD_USING RealType result = 0; - if (!_jacobi_theta_check_q(q, pol, function, &result)) + if (!detail::jacobi_theta_check_q(q, pol, function, &result)) return result; - if (!_jacobi_theta_check_z(z, pol, function, &result)) + if (!detail::jacobi_theta_check_z(z, pol, function, &result)) return result; if (q < exp(-constants::pi())) - return RealType(1) + _jacobi_theta4m1_series(z, _jacobi_theta_q_power{q}, pol); + return RealType(1) + detail::jacobi_theta4m1_series(z, detail::jacobi_theta_q_power{q}, pol); - return _jacobi_theta4_dispatch(z, _jacobi_theta_exponents::from_nome(q), pol); + return detail::jacobi_theta4_dispatch(z, detail::jacobi_theta_exponents::from_nome(q), pol); } // Begin public API