.. _examples: Example programs =============================================================================== .. highlight:: text The *examples* directory (https://github.com/fredrik-johansson/arb/tree/master/examples) contains several complete C programs, which are documented below. Running:: make examples will compile the programs and place the binaries in ``build/examples``. pi.c ------------------------------------------------------------------------------- This program computes `\pi` to an accuracy of roughly *n* decimal digits by calling the :func:`arb_const_pi` function with a working precision of roughly `n \log_2(10)` bits. Sample output, computing `\pi` to one million digits:: > build/examples/pi 1000000 precision = 3321933 bits... cpu/wall(s): 0.243 0.244 virt/peak/res/peak(MB): 24.46 30.44 8.73 14.42 [3.14159265358979323846{...999959 digits...}42209010610577945815 +/- 1.38e-1000000] The program prints an interval guaranteed to contain `\pi`, and where all displayed digits are correct up to an error of plus or minus one unit in the last place (see :func:`arb_printn`). By default, only the first and last few digits are printed. Pass 0 as a second argument to print all digits (or pass *m* to print *m* + 1 leading and *m* trailing digits, as above with the default *m* = 20). The program can optionally compute various other constants, and can use multiple threads:: > build/examples/pi 1000000 -threads 4 precision = 3321933 bits... cpu/wall(s): 0.265 0.147 virt/peak/res/peak(MB): 241.95 422.15 13.33 17.54 [3.14159265358979323846{...999959 digits...}42209010610577945815 +/- 1.38e-1000000] > build/examples/pi 1000000 -constant e precision = 3321933 bits... cpu/wall(s): 0.09 0.09 virt/peak/res/peak(MB): 25.56 29.19 9.58 13.11 [2.71828182845904523536{...999959 digits...}01379817644769422819 +/- 1.39e-1000000] zeta_zeros.c ------------------------------------------------------------------------------- This program computes one or several consecutive zeros of the Riemann zeta function on the critical line:: > build/examples/zeta_zeros -n 1 -count 10 -digits 30 1 14.1347251417346937904572519836 2 21.0220396387715549926284795939 3 25.0108575801456887632137909926 4 30.4248761258595132103118975306 5 32.9350615877391896906623689641 6 37.5861781588256712572177634807 7 40.9187190121474951873981269146 8 43.3270732809149995194961221654 9 48.0051508811671597279424727494 10 49.7738324776723021819167846786 cpu/wall(s): 0.01 0.01 virt/peak/res/peak(MB): 21.28 21.28 7.29 7.29 Five zeros starting with the millionth:: > build/examples/zeta_zeros -n 1000000 -count 5 -digits 20 1000000 600269.67701244495552 1000001 600270.30109071169866 1000002 600270.74787059436613 1000003 600271.48637367364820 1000004 600271.76148042593778 cpu/wall(s): 0.03 0.03 virt/peak/res/peak(MB): 21.41 21.41 7.41 7.41 The program supports the following options:: zeta_zeros [-n n] [-count n] [-prec n] [-digits n] [-threads n] [-platt] [-noplatt] [-v] [-verbose] [-h] [-help] With ``-platt``, Platt's algorithm is used, which may be faster when computing many zeros of large index simultaneously. bernoulli.c ------------------------------------------------------------------------------- This program benchmarks computing the nth Bernoulli number exactly:: > build/examples/bernoulli 1000000 -threads 8 cpu/wall(s): 27.227 5.836 virt/peak/res/peak(MB): 573.47 731.39 73.23 165.13 class_poly.c ------------------------------------------------------------------------------- This program benchmarks computing Hilbert class polynomials:: > build/examples/class_poly -1000004 -threads 8 cpu/wall(s): 6.932 1.478 virt/peak/res/peak(MB): 535.27 653.18 71.02 100.65 degree = 624, bits = -37823 hilbert_matrix.c ------------------------------------------------------------------------------- Given an input integer *n*, this program accurately computes the determinant of the *n* by *n* Hilbert matrix. Hilbert matrices are notoriously ill-conditioned: although the entries are close to unit magnitude, the determinant `h_n` decreases superexponentially (nearly as `1/4^{n^2}`) as a function of *n*. This program automatically doubles the working precision until the ball computed for `h_n` by :func:`arb_mat_det` does not contain zero. Sample output:: $ build/examples/hilbert_matrix 200 prec=20: [+/- 1.32e-335] prec=40: [+/- 1.63e-545] prec=80: [+/- 1.30e-933] prec=160: [+/- 3.62e-1926] prec=320: [+/- 1.81e-4129] prec=640: [+/- 3.84e-8838] prec=1280: [2.955454297e-23924 +/- 8.29e-23935] success! cpu/wall(s): 8.494 8.513 virt/peak/res/peak(MB): 134.98 134.98 111.57 111.57 Called with ``-eig n``, instead of computing the determinant, the program computes the smallest eigenvalue of the Hilbert matrix (in fact, it isolates all eigenvalues and prints the smallest eigenvalue):: $ build/examples/hilbert_matrix -eig 50 prec=20: nan prec=40: nan prec=80: nan prec=160: nan prec=320: nan prec=640: [1.459157797e-74 +/- 2.49e-84] success! cpu/wall(s): 1.84 1.841 virt/peak/res/peak(MB): 33.97 33.97 10.51 10.51 keiper_li.c ------------------------------------------------------------------------------- Given an input integer *n*, this program rigorously computes numerical values of the Keiper-Li coefficients `\lambda_0, \ldots, \lambda_n`. The Keiper-Li coefficients have the property that `\lambda_n > 0` for all `n > 0` if and only if the Riemann hypothesis is true. This program was used for the record computations described in [Joh2013]_ (the paper describes the algorithm in some more detail). The program takes the following parameters:: keiper_li n [-prec prec] [-threads num_threads] [-out out_file] The program prints the first and last few coefficients. It can optionally write all the computed data to a file. The working precision defaults to a value that should give all the coefficients to a few digits of accuracy, but can optionally be set higher (or lower). On a multicore system, using several threads results in faster execution. Sample output:: > build/examples/keiper_li 1000 -threads 2 zeta: cpu/wall(s): 0.4 0.244 virt/peak/res/peak(MB): 167.98 294.69 5.09 7.43 log: cpu/wall(s): 0.03 0.038 gamma: cpu/wall(s): 0.02 0.016 binomial transform: cpu/wall(s): 0.01 0.018 0: -0.69314718055994530941723212145817656807550013436026 +/- 6.5389e-347 1: 0.023095708966121033814310247906495291621932127152051 +/- 2.0924e-345 2: 0.046172867614023335192864243096033943387066108314123 +/- 1.674e-344 3: 0.0692129735181082679304973488726010689942120263932 +/- 5.0219e-344 4: 0.092197619873060409647627872409439018065541673490213 +/- 2.0089e-343 5: 0.11510854289223549048622128109857276671349132303596 +/- 1.0044e-342 6: 0.13792766871372988290416713700341666356138966078654 +/- 6.0264e-342 7: 0.16063715965299421294040287257385366292282442046163 +/- 2.1092e-341 8: 0.18321945964338257908193931774721859848998098273432 +/- 8.4368e-341 9: 0.20565733870917046170289387421343304741236553410044 +/- 7.5931e-340 10: 0.22793393631931577436930340573684453380748385942738 +/- 7.5931e-339 991: 2.3196617961613367928373899656994682562101430813341 +/- 2.461e-11 992: 2.3203766239254884035349896518332550233162909717288 +/- 9.5363e-11 993: 2.321092061239733282811659116333262802034375592414 +/- 1.8495e-10 994: 2.3218073540188462110258826121503870112747188888893 +/- 3.5907e-10 995: 2.3225217392815185726928702951225314023773358152533 +/- 6.978e-10 996: 2.3232344485814623873333223609413703912358283071281 +/- 1.3574e-09 997: 2.3239447114886014522889542667580382034526509232475 +/- 2.6433e-09 998: 2.3246517591032700808344143240352605148856869322209 +/- 5.1524e-09 999: 2.3253548275861382119812576052060526988544993162101 +/- 1.0053e-08 1000: 2.3260531616864664574065046940832238158044982041872 +/- 3.927e-08 virt/peak/res/peak(MB): 170.18 294.69 7.51 7.51 logistic.c ------------------------------------------------------------------------------- This program computes the *n*-th iterate of the logistic map defined by `x_{n+1} = r x_n (1 - x_n)` where `r` and `x_0` are given. It takes the following parameters:: logistic n [x_0] [r] [digits] The inputs `x_0`, *r* and *digits* default to 0.5, 3.75 and 10 respectively. The computation is automatically restarted with doubled precision until the result is accurate to *digits* decimal digits. Sample output:: > build/examples/logistic 10 Trying prec=64 bits...success! cpu/wall(s): 0 0.001 x_10 = [0.6453672908 +/- 3.10e-11] > build/examples/logistic 100 Trying prec=64 bits...ran out of accuracy at step 18 Trying prec=128 bits...ran out of accuracy at step 53 Trying prec=256 bits...success! cpu/wall(s): 0 0 x_100 = [0.8882939923 +/- 1.60e-11] > build/examples/logistic 10000 Trying prec=64 bits...ran out of accuracy at step 18 Trying prec=128 bits...ran out of accuracy at step 53 Trying prec=256 bits...ran out of accuracy at step 121 Trying prec=512 bits...ran out of accuracy at step 256 Trying prec=1024 bits...ran out of accuracy at step 525 Trying prec=2048 bits...ran out of accuracy at step 1063 Trying prec=4096 bits...ran out of accuracy at step 2139 Trying prec=8192 bits...ran out of accuracy at step 4288 Trying prec=16384 bits...ran out of accuracy at step 8584 Trying prec=32768 bits...success! cpu/wall(s): 0.859 0.858 x_10000 = [0.8242048008 +/- 4.35e-11] > build/examples/logistic 1234 0.1 3.99 30 Trying prec=64 bits...ran out of accuracy at step 0 Trying prec=128 bits...ran out of accuracy at step 10 Trying prec=256 bits...ran out of accuracy at step 76 Trying prec=512 bits...ran out of accuracy at step 205 Trying prec=1024 bits...ran out of accuracy at step 461 Trying prec=2048 bits...ran out of accuracy at step 974 Trying prec=4096 bits...success! cpu/wall(s): 0.009 0.009 x_1234 = [0.256445391958651410579677945635 +/- 3.92e-31] real_roots.c ------------------------------------------------------------------------------- This program isolates the roots of a function on the interval `(a,b)` (where *a* and *b* are input as double-precision literals) using the routines in the :ref:`arb_calc ` module. The program takes the following arguments:: real_roots function a b [-refine d] [-verbose] [-maxdepth n] [-maxeval n] [-maxfound n] [-prec n] The following functions (specified by an integer code) are implemented: * 0 - `Z(x)` (Riemann-Siegel Z-function) * 1 - `\sin(x)` * 2 - `\sin(x^2)` * 3 - `\sin(1/x)` * 4 - `\operatorname{Ai}(x)` (Airy function) * 5 - `\operatorname{Ai}'(x)` (Airy function) * 6 - `\operatorname{Bi}(x)` (Airy function) * 7 - `\operatorname{Bi}'(x)` (Airy function) The following options are available: * ``-refine d``: If provided, after isolating the roots, attempt to refine the roots to *d* digits of accuracy using a few bisection steps followed by Newton's method with adaptive precision, and then print them. * ``-verbose``: Print more information. * ``-maxdepth n``: Stop searching after *n* recursive subdivisions. * ``-maxeval n``: Stop searching after approximately *n* function evaluations (the actual number evaluations will be a small multiple of this). * ``-maxfound n``: Stop searching after having found *n* isolated roots. * ``-prec n``: Working precision to use for the root isolation. With *function* 0, the program isolates roots of the Riemann zeta function on the critical line, and guarantees that no roots are missed (see `zeta_zeros.c` for a far more efficient way to do this):: > build/examples/real_roots 0 0.0 50.0 -verbose interval: [0, 50] maxdepth = 30, maxeval = 100000, maxfound = 100000, low_prec = 30 found isolated root in: [14.111328125, 14.16015625] found isolated root in: [20.99609375, 21.044921875] found isolated root in: [25, 25.048828125] found isolated root in: [30.419921875, 30.4443359375] found isolated root in: [32.91015625, 32.958984375] found isolated root in: [37.548828125, 37.59765625] found isolated root in: [40.91796875, 40.966796875] found isolated root in: [43.310546875, 43.3349609375] found isolated root in: [47.998046875, 48.0224609375] found isolated root in: [49.755859375, 49.7802734375] --------------------------------------------------------------- Found roots: 10 Subintervals possibly containing undetected roots: 0 Function evaluations: 3058 cpu/wall(s): 0.202 0.202 virt/peak/res/peak(MB): 26.12 26.14 2.76 2.76 Find just one root and refine it to approximately 75 digits:: > build/examples/real_roots 0 0.0 50.0 -maxfound 1 -refine 75 interval: [0, 50] maxdepth = 30, maxeval = 100000, maxfound = 1, low_prec = 30 refined root (0/8): [14.134725141734693790457251983562470270784257115699243175685567460149963429809 +/- 2.57e-76] --------------------------------------------------------------- Found roots: 1 Subintervals possibly containing undetected roots: 7 Function evaluations: 761 cpu/wall(s): 0.055 0.056 virt/peak/res/peak(MB): 26.12 26.14 2.75 2.75 Find the first few roots of an Airy function and refine them to 50 digits each:: > build/examples/real_roots 4 -10 0 -refine 50 interval: [-10, 0] maxdepth = 30, maxeval = 100000, maxfound = 100000, low_prec = 30 refined root (0/6): [-9.022650853340980380158190839880089256524677535156083 +/- 4.85e-52] refined root (1/6): [-7.944133587120853123138280555798268532140674396972215 +/- 1.92e-52] refined root (2/6): [-6.786708090071758998780246384496176966053882477393494 +/- 3.84e-52] refined root (3/6): [-5.520559828095551059129855512931293573797214280617525 +/- 1.05e-52] refined root (4/6): [-4.087949444130970616636988701457391060224764699108530 +/- 2.46e-52] refined root (5/6): [-2.338107410459767038489197252446735440638540145672388 +/- 1.48e-52] --------------------------------------------------------------- Found roots: 6 Subintervals possibly containing undetected roots: 0 Function evaluations: 200 cpu/wall(s): 0.003 0.003 virt/peak/res/peak(MB): 26.12 26.14 2.24 2.24 Find roots of `\sin(x^2)` on `(0,100)`. The algorithm cannot isolate the root at `x = 0` (it is at the endpoint of the interval, and in any case a root of multiplicity higher than one). The failure is reported:: > build/examples/real_roots 2 0 100 interval: [0, 100] maxdepth = 30, maxeval = 100000, maxfound = 100000, low_prec = 30 --------------------------------------------------------------- Found roots: 3183 Subintervals possibly containing undetected roots: 1 Function evaluations: 34058 cpu/wall(s): 0.032 0.032 virt/peak/res/peak(MB): 26.32 26.37 2.04 2.04 This does not miss any roots:: > build/examples/real_roots 2 1 100 interval: [1, 100] maxdepth = 30, maxeval = 100000, maxfound = 100000, low_prec = 30 --------------------------------------------------------------- Found roots: 3183 Subintervals possibly containing undetected roots: 0 Function evaluations: 34039 cpu/wall(s): 0.023 0.023 virt/peak/res/peak(MB): 26.32 26.37 2.01 2.01 Looking for roots of `\sin(1/x)` on `(0,1)`, the algorithm finds many roots, but will never find all of them since there are infinitely many:: > build/examples/real_roots 3 0.0 1.0 interval: [0, 1] maxdepth = 30, maxeval = 100000, maxfound = 100000, low_prec = 30 --------------------------------------------------------------- Found roots: 10198 Subintervals possibly containing undetected roots: 24695 Function evaluations: 202587 cpu/wall(s): 0.171 0.171 virt/peak/res/peak(MB): 28.39 30.38 4.05 4.05 Remark: the program always computes rigorous containing intervals for the roots, but the accuracy after refinement could be less than *d* digits. poly_roots.c ------------------------------------------------------------------------------- This program finds the complex roots of an integer polynomial by calling :func:`arb_fmpz_poly_complex_roots`, which in turn calls :func:`acb_poly_find_roots` with increasing precision until the roots certainly have been isolated. The program takes the following arguments:: poly_roots [-refine d] [-print d] Isolates all the complex roots of a polynomial with integer coefficients. If -refine d is passed, the roots are refined to a relative tolerance better than 10^(-d). By default, the roots are only computed to sufficient accuracy to isolate them. The refinement is not currently done efficiently. If -print d is passed, the computed roots are printed to d decimals. By default, the roots are not printed. The polynomial can be specified by passing the following as : a Easy polynomial 1 + 2x + ... + (n+1)x^n t Chebyshev polynomial T_n u Chebyshev polynomial U_n p Legendre polynomial P_n c Cyclotomic polynomial Phi_n s Swinnerton-Dyer polynomial S_n b Bernoulli polynomial B_n w Wilkinson polynomial W_n e Taylor series of exp(x) truncated to degree n m The Mignotte-like polynomial x^n + (100x+1)^m, n > m coeffs c0 + c1 x + ... + cn x^n Concatenate to multiply polynomials, e.g.: p 5 t 6 coeffs 1 2 3 for P_5(x)*T_6(x)*(1+2x+3x^2) This finds the roots of the Wilkinson polynomial with roots at the positive integers 1, 2, ..., 100:: > build/examples/poly_roots -print 15 w 100 computing squarefree factorization... cpu/wall(s): 0.001 0.001 roots with multiplicity 1 searching for 100 roots, 100 deflated prec=32: 0 isolated roots | cpu/wall(s): 0.098 0.098 prec=64: 0 isolated roots | cpu/wall(s): 0.247 0.247 prec=128: 0 isolated roots | cpu/wall(s): 0.498 0.497 prec=256: 0 isolated roots | cpu/wall(s): 0.713 0.713 prec=512: 100 isolated roots | cpu/wall(s): 0.104 0.105 done! [1.00000000000000 +/- 3e-20] [2.00000000000000 +/- 3e-19] [3.00000000000000 +/- 1e-19] [4.00000000000000 +/- 1e-19] [5.00000000000000 +/- 1e-19] ... [96.0000000000000 +/- 1e-17] [97.0000000000000 +/- 1e-17] [98.0000000000000 +/- 3e-17] [99.0000000000000 +/- 3e-17] [100.000000000000 +/- 3e-17] cpu/wall(s): 1.664 1.664 This finds the roots of a Bernoulli polynomial which has both real and complex roots:: > build/examples/poly_roots -refine 100 -print 20 b 16 computing squarefree factorization... cpu/wall(s): 0.001 0 roots with multiplicity 1 searching for 16 roots, 16 deflated prec=32: 16 isolated roots | cpu/wall(s): 0.006 0.006 prec=64: 16 isolated roots | cpu/wall(s): 0.001 0.001 prec=128: 16 isolated roots | cpu/wall(s): 0.001 0.001 prec=256: 16 isolated roots | cpu/wall(s): 0.001 0.002 prec=512: 16 isolated roots | cpu/wall(s): 0.002 0.001 done! [-0.94308706466055783383 +/- 2.02e-21] [-0.75534059252067985752 +/- 2.70e-21] [-0.24999757119077421009 +/- 4.27e-21] [0.24999757152512726002 +/- 4.43e-21] [0.75000242847487273998 +/- 4.43e-21] [1.2499975711907742101 +/- 1.43e-20] [1.7553405925206798575 +/- 1.74e-20] [1.9430870646605578338 +/- 3.21e-20] [-0.99509334829256233279 +/- 9.42e-22] + [0.44547958157103608805 +/- 3.59e-21]*I [-0.99509334829256233279 +/- 9.42e-22] + [-0.44547958157103608805 +/- 3.59e-21]*I [1.9950933482925623328 +/- 1.10e-20] + [0.44547958157103608805 +/- 3.59e-21]*I [1.9950933482925623328 +/- 1.10e-20] + [-0.44547958157103608805 +/- 3.59e-21]*I [-0.92177327714429290564 +/- 4.68e-21] + [-1.0954360955079385542 +/- 1.71e-21]*I [-0.92177327714429290564 +/- 4.68e-21] + [1.0954360955079385542 +/- 1.71e-21]*I [1.9217732771442929056 +/- 3.54e-20] + [1.0954360955079385542 +/- 1.71e-21]*I [1.9217732771442929056 +/- 3.54e-20] + [-1.0954360955079385542 +/- 1.71e-21]*I cpu/wall(s): 0.011 0.012 Roots are automatically separated by multiplicity by performing an initial squarefree factorization:: > build/examples/poly_roots -print 5 p 5 p 5 t 7 coeffs 1 5 10 10 5 1 computing squarefree factorization... cpu/wall(s): 0 0 roots with multiplicity 1 searching for 6 roots, 3 deflated prec=32: 3 isolated roots | cpu/wall(s): 0 0.001 done! [-0.97493 +/- 2.10e-6] [-0.78183 +/- 1.49e-6] [-0.43388 +/- 3.75e-6] [0.43388 +/- 3.75e-6] [0.78183 +/- 1.49e-6] [0.97493 +/- 2.10e-6] roots with multiplicity 2 searching for 4 roots, 2 deflated prec=32: 2 isolated roots | cpu/wall(s): 0 0 done! [-0.90618 +/- 1.56e-7] [-0.53847 +/- 6.91e-7] [0.53847 +/- 6.91e-7] [0.90618 +/- 1.56e-7] roots with multiplicity 3 searching for 1 roots, 0 deflated prec=32: 0 isolated roots | cpu/wall(s): 0 0 done! 0 roots with multiplicity 5 searching for 1 roots, 1 deflated prec=32: 1 isolated roots | cpu/wall(s): 0 0 done! -1.0000 cpu/wall(s): 0 0.001 zeta_zeros.c ------------------------------------------------------------------------------- This program finds the imaginary parts of consecutive nontrivial zeros of the Riemann zeta function by calling either :func:`acb_dirichlet_hardy_z_zeros` or :func:`acb_dirichlet_platt_local_hardy_z_zeros` depending on the height of the zeros and the number of zeros requested. The program takes the following arguments:: zeta_zeros [-n n] [-count n] [-prec n] [-threads n] [-platt] [-noplatt] [-v] [-verbose] [-h] [-help] > build/examples/zeta_zeros -n 1048449114 -count 2 1048449114 [388858886.0022851217767970582 +/- 7.46e-20] 1048449115 [388858886.0023936897027167201 +/- 7.59e-20] cpu/wall(s): 0.255 0.255 virt/peak/res/peak(MB): 26.77 26.77 7.88 7.88 complex_plot.c ------------------------------------------------------------------------------- This program plots one of the predefined functions over a complex interval `[x_a, x_b] + [y_a, y_b]i` using domain coloring, at a resolution of *xn* times *yn* pixels. The program takes the parameters:: complex_plot [-range xa xb ya yb] [-size xn yn] [-color n] [-threads n] Defaults parameters are `[-10,10] + [-10,10]i` and *xn* = *yn* = 512. A color function can be selected with -color. Valid options are 0 (phase=hue, magnitude=brightness) and 1 (phase only, white-gold-black-blue-white counterclockwise). The output is written to ``arbplot.ppm``. If you have ImageMagick, run ``convert arbplot.ppm arbplot.png`` to get a PNG. Function codes ```` are: * ``gamma`` - Gamma function * ``digamma`` - Digamma function * ``lgamma`` - Logarithmic gamma function * ``zeta`` - Riemann zeta function * ``erf`` - Error function * ``ai`` - Airy function Ai * ``bi`` - Airy function Bi * ``besselj`` - Bessel function `J_0` * ``bessely`` - Bessel function `Y_0` * ``besseli`` - Bessel function `I_0` * ``besselk`` - Bessel function `K_0` * ``modj`` - Modular j-function * ``modeta`` - Dedekind eta function * ``barnesg`` - Barnes G-function * ``agm`` - Arithmetic geometric mean The function is just sampled at point values; no attempt is made to resolve small features by adaptive subsampling. For example, the following plots the Riemann zeta function around a portion of the critical strip with imaginary part between 100 and 140:: > build/examples/complex_plot zeta -range -10 10 100 140 -size 256 512 For parallel computation on a multicore system, use ``-threads n``. lvalue.c ------------------------------------------------------------------------------- This program evaluates Dirichlet L-functions. It takes the following input:: > build/examples/lvalue lvalue [-character q n] [-re a] [-im b] [-prec p] [-z] [-deflate] [-len l] Print value of Dirichlet L-function at s = a+bi. Default a = 0.5, b = 0, p = 53, (q, n) = (1, 0) (Riemann zeta) [-z] - compute Z(s) instead of L(s) [-deflate] - remove singular term at s = 1 [-len l] - compute l terms in Taylor series at s Evaluating the Riemann zeta function and the Dirichlet beta function at `s = 2`:: > build/examples/lvalue -re 2 -prec 128 L(s) = [1.64493406684822643647241516664602518922 +/- 4.37e-39] cpu/wall(s): 0.001 0.001 virt/peak/res/peak(MB): 26.86 26.88 2.05 2.05 > build/examples/lvalue -character 4 3 -re 2 -prec 128 L(s) = [0.91596559417721901505460351493238411077 +/- 7.86e-39] cpu/wall(s): 0.002 0.003 virt/peak/res/peak(MB): 26.86 26.88 2.31 2.31 Evaluating the L-function for character number 101 modulo 1009 at `s = 1/2` and `s = 1`:: > build/examples/lvalue -character 1009 101 L(s) = [-0.459256562383872 +/- 5.24e-16] + [1.346937111206009 +/- 3.03e-16]*I cpu/wall(s): 0.012 0.012 virt/peak/res/peak(MB): 26.86 26.88 2.30 2.30 > build/examples/lvalue -character 1009 101 -re 1 L(s) = [0.657952586112728 +/- 6.02e-16] + [1.004145273214022 +/- 3.10e-16]*I cpu/wall(s): 0.017 0.018 virt/peak/res/peak(MB): 26.86 26.88 2.30 2.30 Computing the first few coefficients in the Laurent series of the Riemann zeta function at `s = 1`:: > build/examples/lvalue -re 1 -deflate -len 8 L(s) = [0.577215664901532861 +/- 5.29e-19] L'(s) = [0.072815845483676725 +/- 2.68e-19] [x^2] L(s+x) = [-0.004845181596436159 +/- 3.87e-19] [x^3] L(s+x) = [-0.000342305736717224 +/- 4.20e-19] [x^4] L(s+x) = [9.6890419394471e-5 +/- 2.40e-19] [x^5] L(s+x) = [-6.6110318108422e-6 +/- 4.51e-20] [x^6] L(s+x) = [-3.316240908753e-7 +/- 3.85e-20] [x^7] L(s+x) = [1.0462094584479e-7 +/- 7.78e-21] cpu/wall(s): 0.003 0.004 virt/peak/res/peak(MB): 26.86 26.88 2.30 2.30 Evaluating the Riemann zeta function near the first nontrivial root:: > build/examples/lvalue -re 0.5 -im 14.134725 L(s) = [1.76743e-8 +/- 1.93e-14] + [-1.110203e-7 +/- 2.84e-14]*I cpu/wall(s): 0.001 0.001 virt/peak/res/peak(MB): 26.86 26.88 2.31 2.31 > build/examples/lvalue -z -re 14.134725 -prec 200 Z(s) = [-1.12418349839417533300111494358128257497862927935658e-7 +/- 4.62e-58] cpu/wall(s): 0.001 0.001 virt/peak/res/peak(MB): 26.86 26.88 2.57 2.57 > build/examples/lvalue -z -re 14.134725 -len 4 Z(s) = [-1.124184e-7 +/- 7.00e-14] Z'(s) = [0.793160414884 +/- 4.09e-13] [x^2] Z(s+x) = [0.065164586492 +/- 5.39e-13] [x^3] Z(s+x) = [-0.020707762705 +/- 5.37e-13] cpu/wall(s): 0.002 0.003 virt/peak/res/peak(MB): 26.86 26.88 2.57 2.57 lcentral.c ------------------------------------------------------------------------------- This program computes the central value `L(1/2)` for each Dirichlet L-function character modulo *q* for each *q* in the range *qmin* to *qmax*. Usage:: > build/examples/lcentral Computes central values (s = 0.5) of Dirichlet L-functions. usage: build/examples/lcentral [--quiet] [--check] [--prec ] qmin qmax The first few values:: > build/examples/lcentral 1 8 3,2: [0.48086755769682862618122006324 +/- 7.35e-30] 4,3: [0.66769145718960917665869092930 +/- 1.62e-30] 5,2: [0.76374788011728687822451215264 +/- 2.32e-30] + [0.21696476751886069363858659310 +/- 3.06e-30]*I 5,4: [0.23175094750401575588338366176 +/- 2.21e-30] 5,3: [0.76374788011728687822451215264 +/- 2.32e-30] + [-0.21696476751886069363858659310 +/- 3.06e-30]*I 7,3: [0.71394334376831949285993820742 +/- 1.21e-30] + [0.47490218277139938263745243935 +/- 4.52e-30]*I 7,2: [0.31008936259836766059195052534 +/- 5.29e-30] + [-0.07264193137017790524562171245 +/- 5.48e-30]*I 7,6: [1.14658566690370833367712697646 +/- 1.95e-30] 7,4: [0.31008936259836766059195052534 +/- 5.29e-30] + [0.07264193137017790524562171245 +/- 5.48e-30]*I 7,5: [0.71394334376831949285993820742 +/- 1.21e-30] + [-0.47490218277139938263745243935 +/- 4.52e-30]*I 8,5: [0.37369171291254730738158695002 +/- 4.01e-30] 8,3: [1.10042140952554837756713576997 +/- 3.37e-30] cpu/wall(s): 0.002 0.003 virt/peak/res/peak(MB): 26.32 26.34 2.35 2.35 Testing a large *q*:: > build/examples/lcentral --quiet --check --prec 256 100000 100000 cpu/wall(s): 1.668 1.667 virt/peak/res/peak(MB): 35.67 46.66 11.67 22.61 It is conjectured that the central value never vanishes. Running with ``--check`` verifies that the interval certainly is nonzero. This can fail with insufficient precision:: > build/examples/lcentral --check --prec 15 100000 100000 100000,71877: [0.1 +/- 0.0772] + [+/- 0.136]*I 100000,90629: [2e+0 +/- 0.106] + [+/- 0.920]*I 100000,28133: [+/- 0.811] + [-2e+0 +/- 0.501]*I 100000,3141: [0.8 +/- 0.0407] + [-0.1 +/- 0.0243]*I 100000,53189: [4.0 +/- 0.0826] + [+/- 0.107]*I 100000,53253: [1.9 +/- 0.0855] + [-3.9 +/- 0.0681]*I Value could be zero! 100000,53381: [+/- 0.0329] + [+/- 0.0413]*I Aborted integrals.c ------------------------------------------------------------------------------- This program computes integrals using :func:`acb_calc_integrate`. Invoking the program without parameters shows usage:: > build/examples/integrals Compute integrals using acb_calc_integrate. Usage: integrals -i n [-prec p] [-tol eps] [-twice] [...] -i n - compute integral n (0 <= n <= 23), or "-i all" -prec p - precision in bits (default p = 64) -goal p - approximate relative accuracy goal (default p) -tol eps - approximate absolute error goal (default 2^-p) -twice - run twice (to see overhead of computing nodes) -heap - use heap for subinterval queue -verbose - show information -verbose2 - show more information -deg n - use quadrature degree up to n -eval n - limit number of function evaluations to n -depth n - limit subinterval queue size to n -threads n - use parallel computation with n threads Implemented integrals: I0 = int_0^100 sin(x) dx I1 = 4 int_0^1 1/(1+x^2) dx I2 = 2 int_0^{inf} 1/(1+x^2) dx (using domain truncation) I3 = 4 int_0^1 sqrt(1-x^2) dx I4 = int_0^8 sin(x+exp(x)) dx I5 = int_1^101 floor(x) dx I6 = int_0^1 |x^4+10x^3+19x^2-6x-6| exp(x) dx I7 = 1/(2 pi i) int zeta(s) ds (closed path around s = 1) I8 = int_0^1 sin(1/x) dx (slow convergence, use -heap and/or -tol) I9 = int_0^1 x sin(1/x) dx (slow convergence, use -heap and/or -tol) I10 = int_0^10000 x^1000 exp(-x) dx I11 = int_1^{1+1000i} gamma(x) dx I12 = int_{-10}^{10} sin(x) + exp(-200-x^2) dx I13 = int_{-1020}^{-1010} exp(x) dx (use -tol 0 for relative error) I14 = int_0^{inf} exp(-x^2) dx (using domain truncation) I15 = int_0^1 sech(10(x-0.2))^2 + sech(100(x-0.4))^4 + sech(1000(x-0.6))^6 dx I16 = int_0^8 (exp(x)-floor(exp(x))) sin(x+exp(x)) dx (use higher -eval) I17 = int_0^{inf} sech(x) dx (using domain truncation) I18 = int_0^{inf} sech^3(x) dx (using domain truncation) I19 = int_0^1 -log(x)/(1+x) dx (using domain truncation) I20 = int_0^{inf} x exp(-x)/(1+exp(-x)) dx (using domain truncation) I21 = int_C wp(x)/x^(11) dx (contour for 10th Laurent coefficient of Weierstrass p-function) I22 = N(1000) = count zeros with 0 < t <= 1000 of zeta(s) using argument principle I23 = int_0^{1000} W_0(x) dx I24 = int_0^pi max(sin(x), cos(x)) dx I25 = int_{-1}^1 erf(x/sqrt(0.0002)*0.5+1.5)*exp(-x) dx I26 = int_{-10}^10 Ai(x) dx I27 = int_0^10 (x-floor(x)-1/2) max(sin(x),cos(x)) dx I28 = int_{-1-i}^{-1+i} sqrt(x) dx I29 = int_0^{inf} exp(-x^2+ix) dx (using domain truncation) I30 = int_0^{inf} exp(-x) Ai(-x) dx (using domain truncation) I31 = int_0^pi x sin(x) / (1 + cos(x)^2) dx A few examples:: build/examples/integrals -i 4 I4 = int_0^8 sin(x+exp(x)) dx ... cpu/wall(s): 0.02 0.02 I4 = [0.34740017265725 +/- 3.95e-15] > build/examples/integrals -i 3 -prec 333 -tol 1e-80 I3 = 4 int_0^1 sqrt(1-x^2) dx ... cpu/wall(s): 0.024 0.024 I3 = [3.141592653589793238462643383279502884197169399375105820974944592307816406286209 +/- 4.24e-79] > build/examples/integrals -i 9 -heap I9 = int_0^1 x sin(1/x) dx (slow convergence, use -heap and/or -tol) ... cpu/wall(s): 0.019 0.018 I9 = [0.3785300 +/- 3.17e-8] fpwrap.c ------------------------------------------------------------------------------- This program demonstrates calling the floating-point wrapper:: > build/examples/fpwrap zeta(2) = 1.644934066848226 zeta(0.5 + 123i) = 0.006252861175594465 + 0.08206030514520983i functions_benchmark.c ------------------------------------------------------------------------------- This program benchmarks performance of some standard functions. .. highlight:: c