plotting cube root function

[Alex Ghitza]

-----BEGIN PGP SIGNED MESSAGE-----
Hash: SHA1

This is really upsetting:

sage: plot(lambda x: x^(1/3), -5, 5)

-
---------------------------------------------------------------------------
<type 'exceptions.TypeError'> Traceback (most recent call last)

/home/ghitza/colby/ma311/<ipython console> in <module>()

/opt/sage/local/lib/python2.5/site-packages/sage/plot/plot.py in
__call__(self, funcs, *args, **kwds)
~ 3349 xmax = args[1]
~ 3350 args = args[2:]
- -> 3351 G = self._call(funcs, (xmin, xmax), *args, **kwds)
~ 3352 else:
~ 3353 print "there were %s extra arguments (besides
%s)" % (n, funcs)

/opt/sage/local/lib/python2.5/site-packages/sage/plot/plot.py in
_call(self, funcs, xrange, parametric, polar, label, **kwds)
~ 3419 del options['plot_division']
~ 3420 while i < len(data) - 1:
- -> 3421 if abs(data[i+1][1] - data[i][1]) > max_bend:
~ 3422 x = (data[i+1][0] + data[i][0])/2
~ 3423 try:

<type 'exceptions.TypeError'>: 'float' object is unsubscriptable
- -------------------------------------------------------------------------

Yes, it's related to tickets #2038 and #2045, but much worse, because in
this case not only is the error message "unfriendly", but it's also
unwarranted: the cube root is a well-defined function on all of RR, so
it should be possible to get a plot of it on the interval [-5,5].

This probably got broken when (-1)^(1/3) stopped being -1.

Alex

- --
Alexandru Ghitza
Assistant Professor
Department of Mathematics
Colby College
Waterville, ME 04901
http://bayes.colby.edu/~ghitza/
-----BEGIN PGP SIGNATURE-----
Version: GnuPG v2.0.7 (GNU/Linux)
Comment: Using GnuPG with Mozilla - http://enigmail.mozdev.org

iD8DBQFHulB+dZTaNFFPILgRApr2AJ4p1Vc9PimTKNsFUS/Lw32tSbAvDgCfYjsF
ol3u1oW3kpRM9dK/o08+UAw=
=ZkrJ
-----END PGP SIGNATURE-----

[William Stein]

1. I agree that the error message sucks.

2. You might want to do:

sage: plot(x^(1/3),(0,5)) + plot(-(-x)^(1/3),(-5,0))

3. Even with Python floats, (-1)^(0.3333) would never work in Python,
so just directly plotting x^(1/3) for negative x isn't going to work:

In [1]: a = float(-1); b = float(0.333333333333333333333333333333333)
In [2]: a**b
---------------------------------------------------------------------------
<type 'exceptions.ValueError'> Traceback (most recent call last)
/Users/was/build/sage-2.10.2.alpha0/devel/sage-review/sage/modular/<ipython
console> in <module>()
<type 'exceptions.ValueError'>: negative number cannot be raised to a
fractional power

[Bobby Moretti]

Sage reals are slightly better than python floats:

sage: RR(-1)^(RR(1/3))
0.500000000000000 + 0.866025403784439*I

This computes the principal value cube root, which by convention gives
us the root in the upper half-plane.

But when you are plotting the cube root as a function like that, you
really do want the real root if at all possible. One idea would be to
have some sort of "real only" option which would try to keep all
results real. I'm not sure how this would be implemented...

-Bobby

--
Bobby Moretti
mor...@u.washington.edu

[Carl Witty]

On Feb 18, 7:43 pm, Alex Ghitza <aghi...@gmail.com> wrote:
> This is really upsetting:
>
> sage: plot(lambda x: x^(1/3), -5, 5)
>
> -
> ---------------------------------------------------------------------------
> <type 'exceptions.TypeError'> Traceback (most recent call last)

This works (slowly):

sage: plot(lambda x: RR(x).nth_root(3), -5, 5)

Carl

[William Stein]

Nice!

This is faster:

sage: plot(lambda x: RR(x).nth_root(3), -5, 5, plot_points=20)

[kcrisman]

>
> > On Feb 18, 7:43 pm, Alex Ghitza <aghi...@gmail.com> wrote:
> > > This is really upsetting:
>
> > > sage: plot(lambda x: x^(1/3), -5, 5)
>

If anyone is interested, here is a potential inconsistency:

The first one works and yields something useful and manifestly real;
the second one we already know yields Alex's slew of error messages:

sage: show(plot(x^(1/2),-5,5))
sage: show(plot(x^(1/3),-5,5))

Does the preparser use math.sqrt() for the square root, or maybe does
this get translated straight to maxima which can handle it? Why
should this not otherwise yield the "primitive" square root of a
positive number, which would be the negative root, and hence not work?

I realize that arguments of a mathematical nature have convinced the
team to let (-1)^(1/3) be complex. But in plot (and perhaps in other
functions?), presumably there is implicit recognition that one desires
real answers by the user giving the range (e.g., -5 to 5). Is there
any way to preparse this sort of thing so that it in this case coerces
to Carl's suggestion of RR(x).nth_root(3), or whatever? And the same
thing would go for any other fractional power, of course, if there is
something similar to coerce it to.

Sage will not catch on in the secondary and low-level postsecondary
markets if this doesn't work, so I hope that something like this might
be possible. I agree with those who've said in the past that having a
"real" option would be cumbersome, though.

- kcrisman

[Dan Grayson]

What were the reasons it was decided to make (-1)^(1/3) complex?
Could that decision be reversed? It seems reasonable that x^(p/q)
should be real-valued for negative x as long as q is odd, unless you
think of a real number as simply a complex number whose imaginary part
happens to be zero.

Even if that decision were reversed at a high level, there may be no
way to inform other components, such as plotters, of the decision, I
realize.

Mathematica underwent a similar change in about version 1.2 or 1.3, by
the way, and since then has made (-1)^(1/3) complex.

[Jason Grout]

For lots of discussion of this issue, see the thread:

http://groups.google.com/group/sage-support/browse_thread/thread/5e650e35bd354c17/415d26fe242fd162

One thing to note about Mathematica: they used to have a RealOnly
package that would make (-1)^(1/3) be -1 (I believe this was in even the
last version of Mma). But even that package is now deprecated.

Jason

[William Stein]

When f is a Sage object and one does

sage: plot(f,...)

the plot function first tries f.plot(...).

I think I might not be opposed to f.plot(...) for symbolic objects being
a bit "smarter" in cases like x^(1/3) and drawing something that
seems like the best choice for "plotting symbolic expression f
on the interval (a,b)". For example, it could do something sensible
with powers like x^(1/3) that appear in the symbolic expression.

This is of course a VERY slippery slope -- the subtle bugs
and confusion that could arise from making the value of
(-3)^(1/3) different when plotting x^(1/3) from everything else
make me shudder.

-- William

[Carl Witty]

On Feb 19, 7:35 am, kcrisman <kcris...@gmail.com> wrote:
> If anyone is interested, here is a potential inconsistency:
>
> The first one works and yields something useful and manifestly real;
> the second one we already know yields Alex's slew of error messages:
>
> sage: show(plot(x^(1/2),-5,5))
> sage: show(plot(x^(1/3),-5,5))

I don't understand your question; when I try this example, I get
essentially identical error backtraces for both these examples.

Carl

[kcrisman]

> sage: show(plot(x^(1/2),-5,5))
> sage: show(plot(x^(1/3),-5,5))
>
> I don't understand your question; when I try this example, I get
> essentially identical error backtraces for both these examples.
>
> Carl

Huh! I would attach the graphic I get, but instead I'll post it to
the Files section. Just tried it again with 2.10.1, no problems with
the first, lots with the second.

-kcrisman

[William Stein]

On Feb 19, 2008 7:35 AM, kcrisman <kcri...@gmail.com> wrote:
> Sage will not catch on in the secondary and low-level postsecondary

> markets if this [ploting x^(1/3) from -5 to 5]

> doesn't work, so I hope that something like this might
> be possible.

Are you sure? You do realize that

Plot[x^(1/3), {x,-5,5}]

in Mathematica 6.0 also doesn't work (in the sense that it completely
ignores all the negative values of x).

You could argue probably successfully that Mathematica hasn't
caught on at all in "secondary and low-level postsecondary
markets"...

-- William

[Carl Witty]

OK, now I'm curious. I still get errors (with 2.10.1 and with
2.10.2.alpha0), as do a couple of people I asked to test this on IRC.

What OS/architecture are you using?

Do you get different results than this:
sage: import math
sage: math.sqrt(float(-1))
---------------------------------------------------------------------------
<type 'exceptions.ValueError'> Traceback (most recent call
last)

/home/cwitty/sage/<ipython console> in <module>()

<type 'exceptions.ValueError'>: math domain error

(The Python documentation says that the functions in the math module
are not portable in "exceptional" situations; perhaps we should try to
avoid these functions...)

Carl

[kcrisman]

Sorry I never attached the graphic - it's on my desktop, but the
Google Groups help leads me to believe someone has turned off file
posting on this group (on sage-support it seems to still be available,
though, given recent posts).

>
> OK, now I'm curious.  I still get errors (with 2.10.1 and with
> 2.10.2.alpha0), as do a couple of people I asked to test this on IRC.
>
> What OS/architecture are you using?
>

Mac OSX.4, PowerPC G4

> Do you get different results than this:
> sage: import math
> sage: math.sqrt(float(-1))
> ---------------------------------------------------------------------------
> <type 'exceptions.ValueError'>            Traceback (most recent call
> last)
>
> /home/cwitty/sage/<ipython console> in <module>()
>
> <type 'exceptions.ValueError'>: math domain error

Yes, I do!

sage: import math
sage: math.sqrt(float(-1))

nan

My computer's native Python installation also provides:

Python 2.3.5 (#1, Mar 20 2005, 20:38:20)
[GCC 3.3 20030304 (Apple Computer, Inc. build 1809)] on darwin
Type "help", "copyright", "credits" or "license" for more information.
>>> import math
>>> math.sqrt(float(-1))
nan

So presumably plot is somehow using this. Personally, I think this is
a feature, not a bug - if one specifies an implied real range for a
function, non-real outputs should just be ignored, which should avoid
the messing up of internal calculations that having a "real" method or
something would imply. For instance, on the same computer, the
following plots something nice and without complex branches or
something weird like that:

sage: show(plot(log(x),-5,5))

While looking into this I discovered:

sage: log(0)
---------------------------------------------------------------------------
<type 'exceptions.TypeError'>
[lots of stuff]
<type 'exceptions.TypeError'>: error evaluating "log(0)":
Error executing code in Maxima
CODE:
log(0);
Maxima ERROR:
log(0) has been generated.

Presumably log(0) should provide a more helpful (and much shorter)
error message. And what does this phrase, "log(0) has been
generated," mean?

- kcrisman

[Andy]

I jumped over from the sage-support group, where I had posted a
question about graphing x^(1/3). Alex Ghitza's response to me made
more sense when I read the background discussion here.

I am experiencing the same inconsistency revealed in the latter part
of this discussion. I am also running Sage locally on a MacOS X PPC
G4 machine, and it graphs plot(x^(1/2)) just as I expected. I was
stumped when I uploaded my worksheet to sagenb.org, and got the same
error message for the square root function, which I thought was
restricted in some way to cube roots. (I'm a novice.)

Ditto on the inconsistency of plotting log(x) and math.sqrt(float(-1))
-- which produces nan on my machine, and an error message on
sagenb.org. I thoroughly agree with kcrisman's reasoning that a
implied real range is a feature, not a bug, for plotting.

Andy