March 2017 - mathfun@mailman.xmission.com

[math-fun] More power series Q's
by Henry Baker 28 Mar '17

28 Mar '17

We have sin(3*x)=3*sin(x)-4*sin(x)^3 and sinh(3*x)=3*sinh(x)+4*sinh(x)^3 These give us recursive approximations to sin(x), sinh(x), thanks to RWG/Hakmem. So we have f(3*x)=3*f(x)+-4*f(x)^3 What are the functions f(3*x)=3*f(x)+b*f(x)^3, for some real b ? They seem to live in the land between sin(x) and sinh(x). Actually, I'm more interested in f(x), s.t. f(2*x) = 2*f(x)+b*f(x)^3, for some real (probably integer?) constant b, because I want a representation in which binary shifting is easy(er). If f(x) = x + ..., then we have f(x) = 8 17 7 15 6 13 5 11 4 9 3 7 3035153449 b x 32516717 b x 151 b x 163 b x 37 b x b x --------------------- + ----------------- + ---------- + ---------- + -------- + ----- 137019759699497875200 29401699409082000 3204726525 97389600 771120 945 2 5 3 b x b x + ----- + ---- + x 60 6 How to easily compute these coefficients? If b=4, the coefficient of x^3 becomes 2/3, which matches that in sinh(2x)/2, although f(x) doesn't converge quite as fast. If b=1, f(x) converges somewhat faster, but how fast? Does it converge everywhere? What is its radius of convergence? If f(2*x)=2*x+b*x^3, what is the/a differential equation for f(x)?

2 2

[math-fun] ries strikes again. Mathematica answers!
by Bill Gosper 28 Mar '17

28 Mar '17

Nice! Tnx! --rwg Date: 2017-03-28 10:15 From: Dan Asimov <asimov(a)msri.org> To: math-fun <math-fun(a)mailman.xmission.com> Reply-To: math-fun <math-fun(a)mailman.xmission.com> I took this as a personal challenge. First note that 1/2 ((-I)^(1 + √2) + I^(1 + √2)) = cos(pi(1+sqrt(2))/2) . Let T = pi(1+sqrt(2))/2. Also note that 1/4 ((-I)^(1 + √2) + I^(1 + √2))^2 = (cos T)^2 so that sqrt(-1 + 1/4 ((-I)^(1 + √2) + I^(1 + √2))^2) = i sin T . Then the whole expression is (cos T - i sin T)^(1 + sqrt(2)) + (cos T + i sin T)^(1 + sqrt(2)) or exp( -i T (1 + sqrt(2)) + exp( i T (1 + sqrt(2)) = 2 cos(T(1 + sqrt(2))) = 2 cos(pi(1 + sqrt(2))^2 / 2) = 2 cos((pi/2) (3 + 2 sqrt(2)) = 2 cos(3 pi/2 + pi sqrt(2)) = 2 [cos(3 pi/2) cos(pi sqrt(2)) - sin(3 pi/2) sin(pi sqrt(2))] = 2 sin(pi sqrt(2)) Q.E.D. —Dan Bill Gosper wrote (March 27, 2017 at 5:25 pm): ----- In[69]:= ((1/2 ((-I)^(1 + √2) + I^(1 + √2)) - Sqrt[-1 + 1/4 ((-I)^(1 + √2) + I^(1 + √2))^2])^(-1 + √2) + (1/2 ((-I)^(1 + √2) + I^(1 + √2)) + Sqrt[-1 + 1/4 ((-I)^(1 + √2) + I^(1 + √2))^2])^(-1 + √2)) // FullSimplify Out[69]= 2 Sin[π 2√2] (How??) -----

2 1

[math-fun] ries strikes again. Mathematica answers!
by Bill Gosper 28 Mar '17

28 Mar '17

In[69]:= ((1/2 ((-I)^(1 + √2) + I^(1 + √2)) - Sqrt[-1 + 1/4 ((-I)^(1 + √2) + I^(1 + √2))^2])^(-1 + √2) + (1/2 ((-I)^(1 + √2) + I^(1 + √2)) + Sqrt[-1 + 1/4 ((-I)^(1 + √2) + I^(1 + √2))^2])^(-1 + √2)) // FullSimplify Out[69]= 2 Sin[π 2√2] (How??) --rwg

3 2

Re: [math-fun] posits, IEEE floats, asinh numbers
by W. Edwin Clark 28 Mar '17

28 Mar '17

> > Is there a write-up of posits that actually specifies how the regime bits > work--I couldn't really tell from the Gustafson slides. Is there > aparameter besides the total number of bits that needs to be specified? It should be in this book: The End of Error: Unum Computing (Chapman & Hall/CRC Computational Science) – February 5, 2015 by John L. Gustafson (Author) <https://www.amazon.com/End-Error-Computing-Chapman-Computational/dp/1482239…>

2 2

Re: [math-fun] posits, IEEE floats, asinh numbers
by Leo Broukhis 27 Mar '17

27 Mar '17

I've explained how the regime bits work in one of my recent messages. There are two parameters: total size S and the size of the exponent tail E. The most significant bits of the exponent are expressed in unary ((N+1) 1s, then 0 for non-negative, -N 0s then 1 for negative). The value of the exponent is N*2^E + tail (if abs(N) is too large for the whole tail to fit in the representation, the trailing bits of the tail are assumed to be 0). The mantissa is 1.{remaining bits to the right of the exponent tail} Leo On Mar 27, 2017 12:17, "Mike Speciner" <ms(a)alum.mit.edu> wrote: You're correct if the integer addition produces "mathematically correct" results and we're talking about 32-bit integers and 64-bit IEEE floats. And the point about javascript is right on. I agree that asinh numbers (particularly when scaled so that they handle denormalized floats) would be problematic for integers. I do like the idea of using asinh in place of log for mapping axes on some graphs/plots. Is there a write-up of posits that actually specifies how the regime bits work--I couldn't really tell from the Gustafson slides. Is there a parameter besides the total number of bits that needs to be specified? --ms On 27-Mar-17 15:06, Dave Dyer wrote: > At 11:50 AM 3/27/2017, Mike Speciner wrote: > >> No, e.g., int32 x = -1<<31, y = -1 >> In that case x+y == -1 >> is definitely not equal to (double)x+(double)y >> > That's not well behaved in the integer domain either - you > ought to get a range exception, but of course you never do. > > In cases where the pure integer arithmetic is well behaved, > is there any case where the same calculation in IEEE doubles > gives a different result? > > > > _______________________________________________ > math-fun mailing list > math-fun(a)mailman.xmission.com > https://mailman.xmission.com/cgi-bin/mailman/listinfo/math-fun > > _______________________________________________ math-fun mailing list math-fun(a)mailman.xmission.com https://mailman.xmission.com/cgi-bin/mailman/listinfo/math-fun

1 0

[math-fun] Pafnuti piles on.
by Bill Gosper 27 Mar '17

27 Mar '17

We can play D Wilson's game with Chebyshevs because T_m ○T_n = T_mn . Are these just conjugates of x² - 2? Then why can't I guess the formula for the hyperbolish boundary gosper.org/chebcheb.png ? A possible approach to simplification: Even when it fails at a specific point, e.g., In[432]:= %370 /. x -> -1/Sqrt[5] // FullSimplify Out[432]= 1/10 (6 + (1/ 2 ((-1 - 2 I)^Sqrt[2] + (-1 + 2 I)^Sqrt[2] - . . . (where %370 is T_√2 ○T_√2(x) - 2 x² +1, as in the png), Series can sometimes cough up any number of 0 terms: In[438]:= Series[%370,{x,-1/Sqrt[5],4}]//tim During evaluation of In[438]:= Series::ztest: Unable to decide whether numeric quantities {3/5+1/2 (1/2 5^ . . . (-1+2 I)^Sqrt[2]]} are equal to zero. Assuming they are. During evaluation of In[438]:= 12.774193,6 Out[438]= O[x+1/Sqrt[5]]^5 But if you accept that the expression vanishes at specific points, you might believe that it does for a region peppered with such points. David's remark that there should be branch choices extending the validity of these identities must be true, because of analytic continuation, But it's not clear to me where to switch branches. --rwg gosper.org/cheb1.png shows T_2 ○T_½(x) - x with a point discontinuity at the origin, and T_½ ○T_2(x) - x with an isolated ray discontinuity along the canonical cut. Bruiser: Using the radicals for In[451]:= ChebyshevT[m, x] Out[452]= 1/2 ((x - Sqrt[-1 + x^2])^m + (x + Sqrt[-1 + x^2])^m) show T_mn = T_nm.

1 0

[math-fun] about Viète formula and patterns in the binary expansion of the series
by Simon Plouffe 27 Mar '17

27 Mar '17

Hello everybody, here is something peculiar about Viète formula for Pi. Vieta or Viète formula here : https://proofwiki.org/wiki/Vieta%27s_Formula_for_Pi Viète was the most famous person from Vendée (France). as you may know the formula converges to Pi (or Pi/2 or 2/Pi depending on which side you put the square roots in the expression), Here is the peculiar thing : there is a definite pattern in the successive terms of the binary expansion of those numbers when they converge to Pi. Here is a sample code for Pari-gp for doing this at high precision : \p 4000 p0 = 2000 p1 = sqrt(2^p0) p2 = 2.0 v = p1/p2 for (n=1, 8000 , p1=sqrt(2.0+p1); v=v*p1/2) print(binary(v)) By doing so, it produces a stream of binary digits at each step of the formula. from these lines of bits, I can do 1 image , colourized for a better view of the contrast. And here is the result : http://plouffe.fr/dessins/vi%c3%a8te.png that image is 20000 x 5000 , which is 20000 binary digits of the first 5000 terms of Viète formula. we see the convergence of the formula, if we look from a distance : we see nothing at all BUT if we take a closer look at the patterns in the binary expansion we see this : http://plouffe.fr/dessins/formule%20de%20vi%c3%a8te.PNG Many things can be said about this : - There is a definite pattern in the binary expansion of those numbers that we do not see individually but when they are grouped we see something. - This pattern appears also with the obvious generalization of Viète formula replacing the '2' by 1, 2, 3, 4, 5, ... the pattern is the same. so, numbers like 4*Pi*sqrt(3)/2, Pi/2, 4/5*log(phi)*sqrt(5) , where phi is the golden ratio and sqrt(3)/2*log(2+sqrt(3)). etc. shows similar patterns. - If we do the same method for series like the Newton-Leibniz for Pi/4 : nothing at all appears, the same with sum(1/(n*2^n), n=1..infinity) = log(2), nothing to see either. - Viète process appears so far as being the only one able to produce a pattern, I have tried many , this is the only example so far. - What is surprising is the way the pattern persist, in diagonal as well as vertical. NOTE : to produce your own image : 1- run the small program with Pari-GP. 2- Clean the output to have only ones and zeros for each line 3- Find a way to transform this text file into a pure Unix file with LF at the end of each line OR a pure binary file of let 's say 10000000 million bits. 4- rename the file by adding .raw at the end. 5- Open Photoshop and open the file, the program will detect that it is something in binary and propose to open it in raw format :, just specify the width and depth of the file properly, 6- convert the file to 'indexed colors', 7- change it to a black body color scheme, this will color the ones and the zeros in blue and black or light brown and white depending of your choise, and voilà the image. Any comment on this ? best regards, Simon Plouffe

1 0

Re: [math-fun] about a certain François Mendzina Essomba and surprising pi formulas
by Dan Asimov 27 Mar '17

27 Mar '17

Unfortunately, you have to be logged on to Facebook to see those formulas, but I am not a member (and don't want to be one). Can the formulas be copied and reposted to an open-access web page? —Dan >From: Simon Plouffe <simon.plouffe(a)gmail.com> >Hello funsters and for some Pi fans, > > I was contacted by a certain François Mendzina Essomba from >Cameroun. > > He showed to me some formulas for pi, the first one was kind of >surprising, >then he sent me a couple of others , > > They are surprising, at first sight : it resembles the ones about the >AGM, >but I am not certain, > > you may judge by yourself here, : > > https://www.facebook.com/photo.php?fbid=1261279953954407&set=gm.13739836826… >https://www.facebook.com/photo.php?fbid=1261281097287626&set=pcb.1373984792… > >and here is the facebook page : > >https://www.facebook.com/search/top/?q=Fran%C3%A7ois%20Mendzina%20Essomba > > > Does anyhone has ever see any of this ? not me.

3 3

[math-fun] posits, IEEE floats, asinh numbers
by Henry Baker 27 Mar '17

27 Mar '17

Over the years, I've bored you all to tears talking about asinh numbers. Basically, asinh numbers represent a real number x as the integer: round(alpha*asinh(beta*x)) For example, if alpha=1/asinh(beta), then the asinh # for 1 is also 1. Since asinh(x) is an odd function, we get smooth transitions about 0. This means we can capture a brief segment of the integers within this representation, while also getting a logarithmic representation of numbers with very large absolute values. In particular, if beta ~ 2^(-22) and we have a 32-bit asinh # representation, then all of the signed 16-bit integers are mapped 1-1, but the overall range is ~ +-10^228. If beta ~ 2^(-46) and we have a 64-bit asinh # representation, then all of the signed 32-bit integers are mapped 1-1, but the overall range is ~ +-10^56937. asinh numbers are thus the *smooth* integration of Kahan IEEE "denormalized" floats with normalized floats. asinh numbers don't have a separate exponent and mantissa, and therefore don't require a "regime" number to specify the size of the exponent. Of course, we can adjust alpha & beta so that asinh #'s have the same denorm range as IEEE denorms. I claim that asinh numbers provide a much smoother (and therefore more accurate) set of numbers than Gustafson's "posits". If I have any time, I'll try to re-do some of Gustafson's calculations with asinh numbers.

5 11

Re: [math-fun] Beyond Floats: Next Gen Computer Arithmetic
by Henry Baker 27 Mar '17

27 Mar '17

If I understood Gustafson's talk correctly, IEEE doubles aren't good enough. He showed some examples where IEEE quads (128 bits) weren't good enough. Is there any way you could convert to Maxima bigfloats, do the operation, and convert back? Bigfloats can be set to any precision you want. Alternatively, could you convert to rationals, do the operation, and convert back? (This should be relatively easy -- if slow -- in most Common Lisp implementations.) I'm not interested in speed, but in precision of the *representation*, as divorced from the precision of the operation itself. I'll talk about this again in another post about "asinh numbers", which Yonemoto (Gustafson's programmer) dismissed out of hand to me as requiring table lookup. At 11:03 PM 3/25/2017, you wrote: >For every operation, my implementation converts posits to IEEE doubles, >performs the operation, and converts the result back to a posit by figuring >out the number of fraction bits based on the exponent, and rounding the >scaled fraction to an integer.

3 6