Jonathan Borwein, FRSC cs.dal/~jborwein Canada Research Chair in Collaborative Technology

1. MAA Summer Seminar on Experimental Math in Action (Carleton College July 15-20, 2007) Jonathan Borwein, FRSCwww.cs.dal.ca/~jborwein Canada Research Chair in Collaborative Technology “The question of the ultimate foundations and the ultimate meaning of mathematics remains open: we do not know in what direction it will find its final solution or even whether a final objective answer can be expected at all. 'Mathematizing' may well be a creative activity of man, like language or music, of primary originality, whose historical decisions defy complete objective rationalisation.”Hermann Weyl Revised 18/07/2007 TexPoint fonts used in EMF. Read the TexPoint manual before you delete this box.: AAAAA

2. L2-3 Computer-assisted Discovery & Proof Jonathan Borwein, FRSCwww.cs.dal.ca/~jborwein Canada Research Chairin Collaborative Technology “Elsewhere Kronecker said ``In mathematics, I recognize true scientific value only in concrete mathematical truths, or to put it more pointedly, only in mathematical formulas." ... I would rather say ``computations" than ``formulas", but my view is essentially the same.” Harold Edwards, Essays in Constructive Mathematics, 2004

3. Chapters Two and Three

4. Computer-assisted Discovery and Proof Gfun, Coupon Collecting, Formulae for Pi, and Generating Functions for Riemann’s Zeta Jonathan M. Borwein Dalhousie D-Drive David H Bailey Lawrence Berkeley National Lab Details in Experimental Mathematics in Action, Bailey, Borwein et al, A.K. Peters, 2007. “All truths are easy to understand once they are discovered; the point is to discover them.” – Galileo Galilei

5. Ex 1. Generating Functions in Maple Sums of 2, 3, 4 squares: what we can tell the easy way

6. Ex 2. Generating Functions in ‘gfun’ Identifying and confirming in Maple

7. And what Sloane tells us …

8. And what Sloane tells us … • What a wonderful resource! • The more technical the result, the less we will learn from Sloane and the more from Salvy-Zimmerman

9. Ex 3. Coupon Collecting and Convexity B. The origin of the problem. This arose as the objective function in a 1999 PhD on coupon collection. Ian Affleck wished to show pN was convex on the positive orthant.I hoped not!

10. Coupon Collecting and Convexity B. Doing What is Easy.

11. Coupon Collecting and Convexity B. A Very Convex Integrand. (Is there a direct proof?) A year later, Omar Hijab suggested re-expressing pN as the joint expectation of Poisson distributions. This leads to: Now yi xi yi and standard techniques show 1/pN is concave, as the integrand is. [We can now ignore probability if we wish!] Q.“inclusion-exclusion” convexity: OK for 1/g(x) > 0, g concave.

12. Algorithms Used in Experimental Mathematics • Symbolic computation for algebraic and calculus manipulations (as in Ex. 1, 2 & 3). • Integer-relation methods, especially the “PSLQ” algorithm. • High-precision integer and floating-point arithmetic. • High-precision evaluation of integrals and infinite series summations. • The Wilf-Zeilberger algorithm for proving summation identities. • Iterative approximations to continuous functions. • Identification of functions based on graph characteristics. • Graphics and visualization methods targeted to mathematical objects.

13. Goals for Today • First lecture: fast arithmetic, Newton’s method, PSLQ and numerical quadrature • Second lecture: more quadrature, Wilf-Zeilberger and applications to physics, binomial series and pi. “As you know Mahlburg proved that for every prime p > 3 there are infinitely many pairs (A,B) such that M(r, p, A*n + B) = 0 (mod p) for  r=0,1,...p-1. Here M(r,p,n)= number of partitions of n with crank congruent to r mod p. Actually he proved more. Recently, I have proved the rank-analog. I found a couple of nontrivial examples:(1) N(r, 11, 5^4*11*19^4*n + 4322599) = 0 (mod 11)(2) N(r, 11, 11^2*19^4*n + 172904)  = 0  (mod 11)  (r=0,1,..10)Here N(r,p,n)= number of partitions of n with rank congruent to r mod p. Anyway, I have verified (2) for n=0 using Fortran (Maple has no hope of doing such a computation). I have not even verified (1) for n=0 (I hope it is correct). I hope to find examples for higher primes, next case is p=13. Frank Garvan”

14. Question #1. • Continued fraction is = [1,2,3,4,5,….] • “google” arithmetic continued fraction (or “sloane” ) • http://mathworld.wolfram.com/ContinuedFractionConstant.html tells us that = ratio of Bessel functions indeed: “The general infinite continued fraction with partial quotients that are in arithmetic progression is given by (Schroeppel 1972).” Here I is the modified Bessel function of the first kind

15. Typical Scheme for High-Precision Floating-Point Arithmetic A high-precision number is represented as a string of n + 4 integers (or a string of n + 4 floating-point numbers with integer values): • First word contains sign and n, the number of words. • Second word contains p, the exponent (power of 2b). • Words three through n + 2 contain mantissas m1 through mn. • Words n + 3 and n + 4 are for convenience in arithmetic. • The value is then given by: • For basic arithmetic operations, conventional schemes suffice up to about 1000 digits. Beyond that level, Karatsuba’s algorithm (next slide) or FFTs can be used for significantly faster multiply performance. • Division and square roots can be performed by Newton iterations (next slide). • For transcendental functions, Taylor’s series or (for higher precision) quadratically convergent elliptic algorithms can be used.

16. ….all based on Fast Arithmetic(Complexity Reduction in Action) • Multiplication • Karatsuba multiplication (200 digits +) or Fast Fourier Transform (FFT) … in ranges from 100 to 1,000,000,000,000 digits • The other operations via Newton’s method • Elementary and special functions via Elliptic integrals and Gauss AGM • For example: Karatsuba replaces one ‘times’ by many ‘plus’ FFT multiplication of multi-billion digit numbers reduces centuries to minutes. Trillions must be done with Karatsuba!

17. Newton’s Method for Elementary Operations and Functions Initial guess • 1. Doubles precision at each step Newton is self correcting and quadratically convergent 2. Consequences for work needed: Now multiply by A 3. For the logarithm we approximate by elliptic integrals (AGM) which admit quadratic transformations: near zero 4. We use Newton to obtain the complex exponential So allelementary functions are fast computable Newton’s arcsin

18. DHB’s Arbitrary Precision Computation (ARPREC) Package • Low-level routines written in C++. • C++ and F-90 translation modules permit use with existing programs with only minor code changes. • Double-double (32 digits), quad-double, (64 digits) and arbitrary precision (>64 digits) available. • Special routines for extra-high precision (>1000 dig). • Includes common math functions: sqrt, cos, exp, etc. • PSLQ, root finding, numerical integration. • An interactive “Experimental Mathematician’s Toolkit” employing this software is also available. Available at: http://www.experimentalmath.info Also recommended: GMP/MPFR package, available at http://www.mpfr.org

19. A Matrix Example

20. The key discovery:

21. … the discovery

22. Log Convexity A Full Case Study

23. The PSLQ Integer Relation Algorithm Let (xn) be a vector of real numbers. An integer relation algorithm finds (or excludes) integers (an)such that • At the present time, the PSLQ algorithm of mathematician-sculptor Helaman Fergusonis the best-known integer relation algorithm. • PSLQ was named one of ten “algorithms of the century” by Computing in Science and Engineering. • High precision arithmetic software is required: at least d £ n digits, where d is the size (in digits) of the largest of the integers ak. [APPENDIX II on PSLQ]

24. Ferguson’s Sculpture Time for a movie?

25. Decrease of error = minj |Aj x| in PSLQ Number of iterates

26. Application of PSLQ: Bifurcation Points in Chaos Theory B3 = 3.54409035955… is third bifurcation point of the logistic iteration of chaos theory: In other words, B3 is the smallest r such that the iteration exhibits 8-way periodicity instead of 4-way periodicity. In 1990, a predecessor to PSLQ found that B3 is a root of • Recently B4 was identified as the root of a 240-degree polynomial by a much more challenging computation. These results have subsequently been proven formally (by Groebner basis methods). • An iterative approximation scheme to calculate high-precision values of these constants is described in EMA.

27. Evaluation of Ten Constants from Quantum Field Theory where

28. PSLQ and Sculpture The complement of the figure-eight knot, when viewed in hyperbolic space, has finite volume 2.029883212819307250042… David Broadhurst found, using PSLQ, that this constant is given by the “BBP” formula:

29. Some Supercomputer-Class PSLQ Solutions • Identification of B4, the fourth bifurcation point of the logistic iteration. (earlier) • Integer relation of size 121; 10,000 digit arithmetic. • Identification of Apery sums (later). • 15 integer relation problems, with size up to 118, requiring up to 5,000 digit arithmetic. • Identification of Euler-zeta sums. • Hundreds of integer relation problems, each of size 145 and requiring 5,000 digit arithmetic. • Run on IBM SP parallel system. • Finding relation for root of Lehmer’s polynomial. • Integer relation of size 125; 50,000 digit arithmetic. • Utilizes 3-level, multi-pair parallel PSLQ program. • Run on IBM SP using ARPEC; 16 hours on 64 CPUs.

30. Numerical Integration and PSLQ where is a primitive Dirichlet series modulo three. [Note homogeneity of evaluation.]

31. Numerical Integration: Example 2 This arises in mathematical physics, from analysis of the volumes of ideal tetrahedra in hyperbolic space. This “identity” (one of 998) has now been verified numerically to 20,000 digits, but no proof is known. Note that the integrand function has a nasty singularity.

32. Numerical Integration: Example 3 (Jan 2006) The following integrals arise in Ising theory of mathematical physics: We first showed that this can be transformed to a 1-D integral: where K0 is a modified Bessel function. We then computed 400-digit numerical values, from which we found these results (now proven):

33. Identifying the Limit Using the Inverse Symbolic Calculator We discovered the limit result as follows: We first calculated: We then used the Inverse Symbolic Calculator, anor online numerical constant recognition facility available at: http://oldweb.cecm.sfu.ca/projects/ISC/ISCmain.html or http://ddrive.cs.dal.ca/~isc/ (ISC2.0) Output: Mixed constants, 2 with elementary transforms. 6304735033743867 = sr(2)^2/exp(gamma)^2 In other words, For full details see “An Integral of the Ising Class,” (J. Phys. A, 2007) available at http://crd.lbl.gov/~dhbailey/dhbpapers/IsingBBC.pdf

37. New Ramanujan-Like Identities Guillera has recently found Ramanujan-like identities, including: where Guillera proved the first two of these using the Wilf-Zeilberger algorithm. He ascribed the third to Gourevich, who found it using integer relation methods. Are there any higher-order analogues? Not as far as we can tell

38. Part II JM Borwein and DH Bailey “Anyone who is not shocked by quantum theory has not understood a single word." - Niels Bohr

39. Algorithms Used in Experimental Mathematics (again) • Symbolic computation for algebraic and calculus manipulations. • Integer-relation methods, especially the “PSLQ” algorithm. • High-precision integer and floating-point arithmetic. • High-precision evaluation of integrals and infinite series summations. • The Wilf-Zeilberger algorithm for proving summation identities. • Iterative approximations to continuous functions. • Identification of functions based on graph characteristics. • Graphics and visualization methods targeted to mathematical objects.

40. The Wilf-Zeilberger Algorithmfor Proving Identities • A slick, computer-assisted proof scheme to prove certain types of (holonomic) identities • Provides a nice complement to PSLQ • PSLQ and the like permit one to discover new identitiesbut do not constitute rigorous proof • W-Z methods permit one to prove certain types of identitiesbut do not suggest any means to discover the identity

41. Example Usage of W-Z Recall these experimentally-discovered identities (from Part I): Guillera cunningly started by defining He then used the EKHAD software package to obtain the companion

42. Example Usage of W-Z, II When we define Zeilberger's theorem yields the identity which when written out is A limit argument completes the proof of Guillera’s identities.

43. History of Numerical Quadrature • 1670: Newton devises Newton-Coates integration. • 1740: Thomas Simpson develops Simpson's rule. • 1820: Gauss develops Gaussian quadrature. • 1950-1980: Adaptive quadrature, Romberg integration, Clenshaw-Curtis integration, others. • 1985-1990: Maple and Mathematica feature built-in numerical quadrature facilities. • 2000: Very high-precision (exp-exp) quadrature (1000+ digits). With these extreme-precision values, we can now use PSLQ to obtain analytical evaluations of integrals (not just sums).

44. Goals for Today • First lecture: numerical quadrature more quadrature, Wilf-Zeilberger and applications to physics, binomial series and pi • Second lecture: Continued and exploring strange functions 15million popsicle sticks

45. The Euler-Maclaurin Formula [Here h = (b - a)/n and xj = a + j h. Dm f(x) means m-th derivative of f(x).] Note when f(t) and all of its derivatives are zero at a and b, the error E(h) of a simple block-function approximation to the integral goes to zero more rapidly than any power of h. • See also Borwein-Calkin-Manna, Euler-Boole summation revisited, preprint, July 2007.

46. Block-Function Approximation to the Integral of a Bell-Shaped Function

47. Quadrature and the Euler-Maclaurin Formula Given f(x) defined on (-1,1), employ a function g(t) such that g(t) goes from -1 to 1 over the real line, with g’(t) going to zero for large |t|. Then substituting x = g(t) yields Change of perspective: change the function not the rule [Here xj = g(hj) and wj = g’(hj).] If g’(t) goes to zero rapidly enough for large t, then even if f(x) has an infinite derivative or blow-up singularity at an endpoint, f(g(t)) g’(t) often is a nice bell-shaped function for which the E-M formula applies.

48. Three Suitable ‘g’ Functions • The third & fourth are known as “tanh-sinh” quadrature.

49. Original and Transformed Integrand Function Original function (on [-1,1]): Transformed function using g(t) = erf t:

50. Tanh-Sinh Quadrature Example 1 Let Then PSLQ yields Several general results have now been proven, including