March 2017 - mathfun@mailman.xmission.com

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

27 Mar '17

I think the early largest Mersenne primes were found using a bignum package built on floating point. (E.g. Brillhart's 11213(?)). -------- Original Message -------- Date: 2017-03-26 18:07 From: Henry Baker <hbaker1(a)pipeline.com> To: math-fun <math-fun(a)mailman.xmission.com> Reply-To: math-fun <math-fun(a)mailman.xmission.com> At 01:02 PM 3/26/2017, Gareth McCaughan wrote: This has the drawback that essentially no integers other than >0 and +-1 are *exactly* represented, unless I'm missing something. >I'm not sure how much that matters, but it would be quite a >departure from conventional floating-point representations. I'm not so sure that anyone cares about this anymore. I may have been one of the few people trying to do integer calculations using a floating point unit: January 1992 http://home.pipeline.com/~hbaker1/AB-mod-N.html http://home.pipeline.com/~hbaker1/AB-mod-N.pdf "We show how to compute A*B (mod N) efficiently, for single-precision A,B, and N, on a modern RISC architecture (Intel 80860) in ANSI C. On this architecture, our method computes A*B (mod N) faster than ANSI C computes A%N, for unsigned longs A and N."

2 1

[math-fun] Merely algebraic FullSimplify challenge
by Bill Gosper 27 Mar '17

27 Mar '17

In[50]:= Assuming[Re@z > Sqrt[1/4 + (Im@z)^2/3], FullSimplify[1/2 ((-3 z + 4 z^3 - Sqrt[-1 + (3 z - 4 z^3)^2])^( 1/3) + (-3 z + 4 z^3 + Sqrt[-1 + (3 z - 4 z^3)^2])^(1/3))]] Out[50]= 1/2 ((-3 z + 4 z^3 - Sqrt[-1 + (3 z - 4 z^3)^2])^( 1/3) + (-3 z + 4 z^3 + Sqrt[-1 + (3 z - 4 z^3)^2])^(1/3)) "Should" be z. Some evidence: In[51]:= FullSimplify[%/.z->1/2] Out[51]= (-1)^(1/3) (*Isolated discontinuity*) In[52]:= FullSimplify[%50/.z->11/21] Out[52]= 11/21 In[53]:= FullSimplify[%50/.z->10001/20001] Out[53]= 10001/20001 In[71]:= Series[%50,{z,1/2,9}]//FullSimplify Out[71]= 1/2+(z-1/2)+O[z-1/2]^10 Not f(1/2) + ... ! In[74]:= Series[%50,{z,1,9}]//FullSimplify Out[74]= 1+(z-1)+O[z-1]^(19/2) In[75]:= Series[%50,{z,3/2,9}]//FullSimplify Out[75]= 3/2+(z-3/2)+O[z-3/2]^10 In[76]:= Series[%50,{z,1+I,9}]//FullSimplify Out[76]= (1+I)+(z-(1+I))+O[z-(1+I)]^10. --Bill Gosper

1 0

Re: [math-fun] Brutal FullSimplify challenge II
by Bill Gosper 27 Mar '17

27 Mar '17

Yes, the Computer Algebra viewpoint is that multivalued functions are a thought crime, and we must all agree where the branch cuts are in all the named functions. This just gave me an isolated point singularity in a merely algebraic function that I almost reported as a bug: Evidently, -2 z + (-3 z + 4 z^3 + (1 - 4 z^2) Sqrt[-1 + z^2])^( 1/3) + (-3 z + 4 z^3 + Sqrt[-1 + z^2] (-1 + 4 z^2))^(1/3) is 0 within the semihyperbola Re@z > Sqrt[1/4 + (Im@z)^2/3], discontinuous on the boundary, and smooth elsewhere. But boundary pt z=1/2 is discontinuous with both regions. Almost impossible to find with Plot... --rwg -------- Original Message -------- Date: 2017-03-26 08:43 From: "David Wilson" <davidwwilson(a)comcast.net> To: "'math-fun'" <math-fun(a)mailman.xmission.com> Reply-To: math-fun <math-fun(a)mailman.xmission.com> Apologies if I am seem pedantic I programmed this in C++, and got essentially the same values for c(c(x)) plotted by rwg. For any x < 2, the evaluation of c(c(x)) ultimately ends up computing a^√2 for complex a. Unfortunate, for complex a, a^√2 = exp(√2 (log(a) + 2πi k)) has an infinitude of values (for integer k). For C++, and presumably Mma, log(a) returns the principle value with -π < Im(log(a)) <= π. This happily coincides with real log(a) for positive real a. Regarding the c(c(x)) calculation: If x >= 2, we can restrict calculations to the real realm, take the principle (real) value of a^√2, and end up with c(c(x)) = x^2 -2. If x < 2, however, a^√2 must be evaluated for non-positive-real x, and it is not clear that the principle value is appropriate in any particular calculation. I submit that, when c(c(x)) differs from x^2 - 2, we might restore c(c(x)) = x^2 -2 by using the appropriate (non-principle) value of a^√2 in the appropriate calculations, though it is unclear which value of a^√2 (if any) works for any particular evaluation of c(x). > -----Original Message----- > From: math-fun [mailto:math-fun-bounces@mailman.xmission.com] On > Behalf Of Bill Gosper > Sent: Saturday, March 25, 2017 11:35 PM > To: Wolfram Technical Support > Cc: math-fun(a)mailman.xmission.com > Subject: [math-fun] Brutal FullSimplify challenge II > > In[540]:= $Version > > Out[540]= "11.0.1 for Mac OS X x86 (64-bit) (September 21, 2016)" > > cc[x_] := Nest[((# + √(#^2 - 4))/2)^√2 + ((# - √(#^2 - 4))/2)^√2 &, x, 2] > > In[653]:=FullSimplify@cc@-Sqrt[2 (1 - Sin[√2 π])] == 2 Sin[√2 π] > > Out[653]= > 2^(-1 - Sqrt[ > 2]) (((-Sqrt[1 - Sin[Sqrt[2] \[Pi]]] - > I Sqrt[1 + Sin[Sqrt[2] \[Pi]]])^Sqrt[ > 2] + (-Sqrt[1 - Sin[Sqrt[2] \[Pi]]] + > I Sqrt[1 + Sin[Sqrt[2] \[Pi]]])^Sqrt[ > 2] - \[Sqrt](-2^( > 1 + Sqrt[2]) + (-Sqrt[1 - Sin[Sqrt[2] \[Pi]]] - > I Sqrt[1 + Sin[Sqrt[2] \[Pi]]])^( > 2 Sqrt[2]) + (-Sqrt[1 - Sin[Sqrt[2] \[Pi]]] + > I Sqrt[1 + Sin[Sqrt[2] \[Pi]]])^(2 Sqrt[2])))^Sqrt[ > 2] + ((-Sqrt[1 - Sin[Sqrt[2] \[Pi]]] - > I Sqrt[1 + Sin[Sqrt[2] \[Pi]]])^Sqrt[ > 2] + (-Sqrt[1 - Sin[Sqrt[2] \[Pi]]] + > I Sqrt[1 + Sin[Sqrt[2] \[Pi]]])^Sqrt[ > 2] + \[Sqrt](-2^( > 1 + Sqrt[2]) + (-Sqrt[1 - Sin[Sqrt[2] \[Pi]]] - > I Sqrt[1 + Sin[Sqrt[2] \[Pi]]])^( > 2 Sqrt[2]) + (-Sqrt[1 - Sin[Sqrt[2] \[Pi]]] + > I Sqrt[1 + Sin[Sqrt[2] \[Pi]]])^(2 Sqrt[2])))^Sqrt[2]) == > 2 Sin[Sqrt[2] \[Pi]] > > In[654]:= Chop@N@# & /@ % > > Out[654]= True > --Bill Gosper

1 0

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

26 Mar '17

asinh(x) ~ x for |x|<C, C some a priori constant. You can't get much nicer than linear (in this case, the identity function). %i25) taylor(asinh(x),x,0,9); 3 5 7 9 x 3 x 5 x 35 x (%o25)/T/ x - -- + ---- - ---- + ----- + . . . 6 40 112 1152 Thus, for alpha=1/asinh(beta) and beta ~ 2^(-46), round(alpha*asinh(beta*x)) = x for integer |x|<2^32. At 11:45 AM 3/26/2017, Warren D Smith wrote: >Another problem: >While asinh(x) is nice for large x, it is not so nice for very tiny x. >I think that should cause a very bad precision loss in many real computations.

1 0

[math-fun] about a certain François Mendzina Essomba and surprising pi formulas
by Simon Plouffe 26 Mar '17

26 Mar '17

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. Best regards, Simon Plouffe

2 1

[math-fun] Brutal FullSimplify challenge II
by Bill Gosper 26 Mar '17

26 Mar '17

In[540]:= $Version Out[540]= "11.0.1 for Mac OS X x86 (64-bit) (September 21, 2016)" cc[x_] := Nest[((# + √(#^2 - 4))/2)^√2 + ((# - √(#^2 - 4))/2)^√2 &, x, 2] In[653]:=FullSimplify@cc@-Sqrt[2 (1 - Sin[√2 π])] == 2 Sin[√2 π] Out[653]= 2^(-1 - Sqrt[ 2]) (((-Sqrt[1 - Sin[Sqrt[2] \[Pi]]] - I Sqrt[1 + Sin[Sqrt[2] \[Pi]]])^Sqrt[ 2] + (-Sqrt[1 - Sin[Sqrt[2] \[Pi]]] + I Sqrt[1 + Sin[Sqrt[2] \[Pi]]])^Sqrt[ 2] - \[Sqrt](-2^( 1 + Sqrt[2]) + (-Sqrt[1 - Sin[Sqrt[2] \[Pi]]] - I Sqrt[1 + Sin[Sqrt[2] \[Pi]]])^( 2 Sqrt[2]) + (-Sqrt[1 - Sin[Sqrt[2] \[Pi]]] + I Sqrt[1 + Sin[Sqrt[2] \[Pi]]])^(2 Sqrt[2])))^Sqrt[ 2] + ((-Sqrt[1 - Sin[Sqrt[2] \[Pi]]] - I Sqrt[1 + Sin[Sqrt[2] \[Pi]]])^Sqrt[ 2] + (-Sqrt[1 - Sin[Sqrt[2] \[Pi]]] + I Sqrt[1 + Sin[Sqrt[2] \[Pi]]])^Sqrt[ 2] + \[Sqrt](-2^( 1 + Sqrt[2]) + (-Sqrt[1 - Sin[Sqrt[2] \[Pi]]] - I Sqrt[1 + Sin[Sqrt[2] \[Pi]]])^( 2 Sqrt[2]) + (-Sqrt[1 - Sin[Sqrt[2] \[Pi]]] + I Sqrt[1 + Sin[Sqrt[2] \[Pi]]])^(2 Sqrt[2])))^Sqrt[2]) == 2 Sin[Sqrt[2] \[Pi]] In[654]:= Chop@N@# & /@ % Out[654]= True --Bill Gosper

2 2

Re: [math-fun] Beyond Floats: Next Gen Computer Arithmetic
by Leo Broukhis 26 Mar '17

26 Mar '17

On Sat, 25 Mar 2017 09:50:46 -0400, "J.P. Grossman" <jpg(a)alum.mit.edu> wrote: > Could you expand on what "flaws, defects, and failures" you found? Also, > did you use es = 1 or es = 3? > > I can vary size and es arbitrarily. 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. For very large numbers, when the number of fraction bits is 0, and the exponent may have to be rounded as well, the results may be wrong. I'll have to look further for that; but I doubt if it would have made any difference for this test. (In the following, I'm quoting tersely, as I'm retyping from my work laptop which is a walled garden, but I've verified the signs of all numbers and exponents). For Posit<32,3>: 27 significant digits found Addition/Subtraction neither rounds nor chops. Sticky bit used incorrectly or not at all. FLAW: lack(s) of guard digits or failure(s) to correctly round or chop. Square root is neither chopped nor correctly rounded. Observed errors run from -5.0e-01 to 8.5e+00 ulps. SERIOUS DEFECT: sqrt gets too many last digits wrong. FAILURE: multiplication gets too many last digits wrong. FLAW: Underflow can stick at an allegedly positive value PseudZero that prints out as 2.37389e-66 [according to Gistafson, the above is not a flaw, we'll discount that - LB] Since comparison denies Z = 0, evaluating (Z+Z)/Z should be safe. What the machine gets for (Z+Z)/Z is 1.0e+00 This is a DEFECT. FAILURE: Comparisons involving +--5.26561e+64, +-1.05312e+65 and +-9.49557e-66 are confused by Overflow. SERIOUS DEFECT: Comparison alleges that Z = -5.265615e+64 is too far from sqrt(Z) ^ 2 (2.49999985098838806e-01). [the above I don't understand; the sign seems to be wrong] DEFECT: Badly unbalanced range; UfThold * V = -5.0e-01 is too far from 1. Using Posit<32,1> adds one defect and one serious defect: SERIOUS DEFECT: Range is too narrow, U1^4 Underflows. [U1 is min x : 1-x not equal to 1] DEFECT: computing 5.0e+00 ^ 1.1e+01 yielded 4.8827904e+07 which compared unequal to correct 4.8828416e+07; they differ by -5.12e+02 Indeed, the precision loss with increasing exponent is too sharp in Posit<32,1>. Leo

1 0

Re: [math-fun] Beyond Floats: Next Gen Computer Arithmetic
by Axel Vogt 26 Mar '17

26 Mar '17

Kahan on that approach: https://people.eecs.berkeley.edu/~wkahan/ (scroll down to the bottom) > Subject: [math-fun] Beyond Floats: Next Gen Computer Arithmetic > FYI -- > > https://www.youtube.com/watch?v=aP0Y1uAA-2Y > > Stanford Seminar: Beyond Floating Point: Next Generation Computer Arithmetic > > John L. Gustafson, Natl Univ of Singapore > > https://web.stanford.edu/class/ee380/Abstracts/170201.html > > A new data type called a "posit" is designed for direct drop-in replacement for IEEE Standard 754 floats. Unlike unum arithmetic, posits do not require interval-type mathematics or variable size operands, and they round if an answer is inexact, much the way floats do. However, they provide compelling advantages over floats, including simpler hardware implementation that scales from as few as two-bit operands to thousands of bits. For any bit width, they have a larger dynamic range, higher accuracy, better closure under arithmetic operations, and simpler exception-handling. For example, posits never overflow to infinity or underflow to zero, and there is no "Not-a-Number" (NaN) value. Posits should take up less space to implement in silicon than an IEEE float of the same size. With fewer gate delays per operation as well as lower silicon footprint, the posit operations per second (POPS) supported by a chip can be significantly higher than the FLOPs using similar hardware reso > urces. GPU accelerators, in particular, could do more arithmetic per watt and per dollar yet deliver superior answer quality. > > A series of comprehensive benchmarks compares how many decimals of accuracy can be produced for a set number of bits-per-value, using various number formats. Low-precision posits provide a better solution than "approximate computing" methods that try to tolerate decreases in answer quality. High-precision posits provide better answers (more correct decimals) than floats of the same size, suggesting that in some cases, a 32-bit posit may do a better job than a 64-bit float. In other words, posits beat floats at their own game. > > http://ee380.stanford.edu/Abstracts/170201-slides.pdf

3 2

Re: [math-fun] [seqfan] Re: Is 4 a semi-Fibonacci number?
by David Wilson 26 Mar '17

26 Mar '17

Apologies, then, for being a day late. At any rate, I was curious as to the growth rate of (odd indexed elements) of sF, so I computed a few values of sF(2^k - 1). Looks to grow pretty fast, perhaps like n^(c log(n)). n sF(n) 1 1 3 2 7 5 15 16 31 69 63 430 127 4137 255 64436 511 1676353 1023 74555322 2047 5777029421 4095 792086153688 8191 194591768192733 16383 86534148901444102 32767 70244955881077121873 65535 104827174339054175240700 131071 289320796542222620694103961 262143 1484554051525798547483095566354 524287 14226732533670039175251226319995157 1048575 255658362329675830634139832019198358176 2097151 8646156383344527622677186248446917786683381 4194303 552073085216495053730117030388698424955408298334 8388607 66749153256722612438012197449644764555811889274235129 16777215 15322089889433495509821892166057427813245442064758958688836 33554431 6693576684186674844436583641024041055004958124882988723321728369 67108863 5577272394534357747969836555204334929788023962243060637883527988942186 134217727 8881483115234363964483793627146185513254764624755802869033890445863234590077 268435455 27080472748266980783930080876503920295698702610405257020769093391323735487516698472 536870911 158372376484792407948495752863706348514435815566310070102912509445802547425058287520342413 From: maximilian(a)hasler.fr [mailto:maximilian@hasler.fr] On Behalf Of M. F. Hasler Sent: Saturday, March 25, 2017 5:44 PM To: David Wilson Subject: Re: [seqfan] Re: Is 4 a semi-Fibonacci number? On Sat, Mar 25, 2017 at 12:09 AM, David Wilson <davidwwilson(a)comcast.net> wrote: Also, sF(n) is strictly increasing on odd n, so {sF(2n-1): n >= 1} is the ordered sequence of all range values of sF. 4 is not in that sequence. That's what I said a day ago... <http://list.seqfan.eu/pipermail/seqfan/2017-March/017390.html> ;-) M. > -----Original Message----- > From: SeqFan [mailto:seqfan-bounces@list.seqfan.eu] On Behalf Of > israel(a)math.ubc.ca > Sent: Thursday, March 23, 2017 11:24 AM > To: Sequence Fanatics Discussion list > Subject: [seqfan] Re: Is 4 a semi-Fibonacci number? > > The least n for which any value occurs must be odd, since sF(n) = sF(n/2) > for even n. So if you have ruled out sF(n) = 4 for odd n, it follows that 4 > can never occur. > > Cheers, > Robert > > On Mar 23 2017, Alonso Del Arte wrote:t > > >As you know, 4 is not a Fibonacci number. The Fibonacci function can be > >extended to all real numbers, but to get Fibonacci(x) = 4 requires x be > >what looks like a transcendental number. > > > >But could 4 be a semi-Fibonacci number? (A030067) The definition is sF(1) = > >1, sF(n) = sF(n/2) if n is even, sF(n) = sF(n - 1) + sF(n - 2) if n is odd. > > > >Couple of years ago, Roberg G Wilson v determined that 4 does not occur > >among the first million terms. It's not a rigorous proof, of course, but it > >does suggest that 4 never occurs. > > > >It's fairly easy to prove that sF(n) = 4 is impossible if n is odd. But I > >haven't been able to rule out sF(n) = 4 for n even. That would mean sF(n) = > >8, but the parity of n seems like it could be anything. I can visualize a > >whole tree but many of the branches of that tree might not even exist. > > > >Can anyone make a determination on this question, or is this another one > of > >those plausible but unproven conjectures? > > > >Al

1 0

Re: [math-fun] Beyond Floats: Next Gen Computer Arithmetic
by Dave Dyer 25 Mar '17

25 Mar '17

> >> https://www.youtube.com/watch?v=aP0Y1uAA-2Y >> >> Stanford Seminar: Beyond Floating Point: Next Generation Computer >> Arithmetic Sounds pretty solid to me. The unintended consequence if this were widely adopted as 'plug compatible' for IEEE floats, would be that data dumps which contained binary floats would be aliased. IEEE floats are so ubiquitous, everyone treats them as if they were as god-given as integers. It's dangerous to change them.

1 0