ASKSAGE: Sage Q&A Forum - RSS feedhttps://ask.sagemath.org/questions/Q&A Forum for SageenCopyright Sage, 2010. Some rights reserved under creative commons license.Sat, 14 Sep 2013 12:01:39 +0200contour_plot issue with large numbershttps://ask.sagemath.org/question/10536/contour_plot-issue-with-large-numbers/I was trying to whip up a variant on Tupper's "self-referential" formula graph (nice exposition [here](http://shreevatsa.wordpress.com/2011/04/12/how-does-tuppers-self-referential-formula-work/) if you're not familiar with it) in Sage, but I ran into a graphing issue with largish numbers. In an attempt to nail down the problem, I've tried to reduce this to a minimal case:
x,y = var('x,y')
def mod2(x): return floor(x)-2*floor(x/2)
def mod4(x): return floor(x)-4*floor(x/4)
N=2^53+1
contour_plot(mod2((N+floor(y/4)) / 2^(4*floor(x)+mod4(y))),(x,0,14),(y,0,4))
(I wrote the mod2/4 functions because I kept on running into issues with the % operator and type coercion.)
If you understand the idea behind the Tupper formula, then you know that the binary representation of N is basically a readout of the "pixels" of the eventual graph generated. Now, 2^53+1 is 1000...001 in binary, so there should be
TWO "pixels" lit.
But, when I actually execute the above code on my home computer or in the Sage Cloud, something funny happens. I get a graph where only ONE "pixel" is lit -- the one corresponding to the least significant bit of 2^53+1, which should be in the bottom left of the graph (i.e., near the origin) is mysteriously absent.
However, if you check for...
N=2^52+1
...the graph correctly shows two "pixels". So I first presumed this was some sort of floating point or rounding error. But if you actually check the values returned by the function being graphed, the returned values are always correct. For example, setting (x,y) to a value squarely within the "missing pixel" and manually evaluating...
x,y=0.5,0.5
mod2((N+floor(y/4)) / 2^(4*floor(x)+mod4(y)))
...yields the correct value of 1, even for N=2^53+1! Similarly, setting (x,y) to (0.5,1.5), i.e., outside said "pixel", appropriately yields 0. So I don't get at all why the graphing isn't working.
Any insight would be appreciated.Fri, 13 Sep 2013 15:25:34 +0200https://ask.sagemath.org/question/10536/contour_plot-issue-with-large-numbers/Answer by sappidus for <p>I was trying to whip up a variant on Tupper's "self-referential" formula graph (nice exposition <a href="http://shreevatsa.wordpress.com/2011/04/12/how-does-tuppers-self-referential-formula-work/">here</a> if you're not familiar with it) in Sage, but I ran into a graphing issue with largish numbers. In an attempt to nail down the problem, I've tried to reduce this to a minimal case:</p>
<pre><code>x,y = var('x,y')
def mod2(x): return floor(x)-2*floor(x/2)
def mod4(x): return floor(x)-4*floor(x/4)
N=2^53+1
contour_plot(mod2((N+floor(y/4)) / 2^(4*floor(x)+mod4(y))),(x,0,14),(y,0,4))
</code></pre>
<p>(I wrote the mod2/4 functions because I kept on running into issues with the % operator and type coercion.)</p>
<p>If you understand the idea behind the Tupper formula, then you know that the binary representation of N is basically a readout of the "pixels" of the eventual graph generated. Now, 2^53+1 is 1000...001 in binary, so there should be
TWO "pixels" lit.</p>
<p>But, when I actually execute the above code on my home computer or in the Sage Cloud, something funny happens. I get a graph where only ONE "pixel" is lit -- the one corresponding to the least significant bit of 2^53+1, which should be in the bottom left of the graph (i.e., near the origin) is mysteriously absent.</p>
<p>However, if you check for...</p>
<pre><code>N=2^52+1
</code></pre>
<p>...the graph correctly shows two "pixels". So I first presumed this was some sort of floating point or rounding error. But if you actually check the values returned by the function being graphed, the returned values are always correct. For example, setting (x,y) to a value squarely within the "missing pixel" and manually evaluating...</p>
<pre><code>x,y=0.5,0.5
mod2((N+floor(y/4)) / 2^(4*floor(x)+mod4(y)))
</code></pre>
<p>...yields the correct value of 1, even for N=2^53+1! Similarly, setting (x,y) to (0.5,1.5), i.e., outside said "pixel", appropriately yields 0. So I don't get at all why the graphing isn't working.</p>
<p>Any insight would be appreciated.</p>
https://ask.sagemath.org/question/10536/contour_plot-issue-with-large-numbers/?answer=15450#post-id-15450Ah, I did not realize plotting is by default done with floats. By setting:
f(x,y)=mod2((N+floor(y/4)) / 2^(4*floor(x)+mod4(y)))
g=fast_callable(f, vars=[x,y], domain=RealField(100))
and contour_plot'ing g instead, which forces the use of higher-precision arithmetic, I get what I want.Fri, 13 Sep 2013 21:09:24 +0200https://ask.sagemath.org/question/10536/contour_plot-issue-with-large-numbers/?answer=15450#post-id-15450Comment by ppurka for <p>Ah, I did not realize plotting is by default done with floats. By setting:</p>
<pre><code>f(x,y)=mod2((N+floor(y/4)) / 2^(4*floor(x)+mod4(y)))
g=fast_callable(f, vars=[x,y], domain=RealField(100))
</code></pre>
<p>and contour_plot'ing g instead, which forces the use of higher-precision arithmetic, I get what I want.</p>
https://ask.sagemath.org/question/10536/contour_plot-issue-with-large-numbers/?comment=17014#post-id-17014That's a very good observation. Your solution is also a potential solution to ticket 13559. There is also ticket 15030 that aims to switch the plot functions to use `fast_callable` instead of `fast_float`.Sat, 14 Sep 2013 12:01:39 +0200https://ask.sagemath.org/question/10536/contour_plot-issue-with-large-numbers/?comment=17014#post-id-17014