Opened 12 years ago

Last modified 10 days ago

#8450 needs_review defect

intermediate complex expression in real functions make many plot functions fail

Reported by: jason Owned by: was
Priority: minor Milestone: sage-9.5
Component: graphics Keywords:
Cc: mjo, egourgoulhon, gh-jungmath Merged in:
Authors: Michael Orlitzky Reviewers:
Report Upstream: N/A Work issues:
Branch: u/mkoeppe/ticket/8450 (Commits, GitHub, GitLab) Commit: d8e79086453876c2f1c43931590ffc94801fcaa1
Dependencies: #32234 Stopgaps:

Status badges

Description (last modified by vdelecroix)

All of the following plots fail

x, y = SR.var('x y')
contour_plot(abs(x+i*y), (x,-1,1), (y,-1,1))
density_plot(abs(x+i*y), (x,-1,1), (y,-1,1))
plot3d(abs(x+i*y), (x,-1,1),(y,-1,1))
streamline_plot(abs(x+i*y), (x,-1,1),(y,-1,1))

with

TypeError: unable to coerce to a real number

The culprit is the call to setup_for_eval_on_grid (from sage/plot/misc.py) that tries to compile the symbolic expression with fast_float. But since the expression involves an intermediate complex number the compilation fails. This can be tested with any of the two

fast_float(abs(x + i*y), x, y)
fast_callable(abs(x + i*y), vars=[x,y])

The function compilation succeeds if we ask for a complex function instead

fast_callable(abs(x + i*y), vars=[x,y], domain=complex)

See also this question on ask.sagemath.org.

Change History (29)

comment:1 Changed 11 years ago by jason

#5572 will help solve this.

comment:2 Changed 8 years ago by jdemeyer

  • Milestone changed from sage-5.11 to sage-5.12

comment:3 Changed 8 years ago by vbraun_spam

  • Milestone changed from sage-6.1 to sage-6.2

comment:4 Changed 7 years ago by vbraun_spam

  • Milestone changed from sage-6.2 to sage-6.3

comment:5 Changed 7 years ago by vbraun_spam

  • Milestone changed from sage-6.3 to sage-6.4

comment:6 Changed 2 years ago by vdelecroix

  • Description modified (diff)
  • Milestone changed from sage-6.4 to sage-8.8

comment:7 Changed 2 years ago by vdelecroix

  • Description modified (diff)
  • Summary changed from contour_plot chokes on function which involves imaginary numbers to intermediate complex expression in real functions make many plot functions fail

comment:8 Changed 2 years ago by vdelecroix

  • Description modified (diff)

comment:9 Changed 2 years ago by embray

  • Milestone sage-8.8 deleted

As the Sage-8.8 release milestone is pending, we should delete the sage-8.8 milestone for tickets that are not actively being worked on or that still require significant work to move forward. If you feel that this ticket should be included in the next Sage release at the soonest please set its milestone to the next release milestone (sage-8.9).

comment:10 Changed 2 months ago by mkoeppe

  • Cc mjo added
  • Dependencies set to #32234
  • Milestone set to sage-9.5
  • Status changed from new to needs_info

comment:11 follow-ups: Changed 2 months ago by mjo

  • Authors set to Michael Orlitzky
  • Branch set to u/mjo/ticket/8450
  • Commit set to a773fd3919b69e575d94cdf60135c419b5e08f24
  • Status changed from needs_info to needs_review

This is an "easy fix" but causes a cascade of other issues because removing domain=float turns a lot of things that used to be Nan into TypeError, ValueError, ZeroDivisionError, etc.

I'm still running a ptestlong on this, but everything under sage/plot now passes. Here are some thoughts:

  1. Do we really want to support passing (for example) the integer zero as a function to be plotted? We have doctests that check things like plot(ZZ(0), x, 0, 1). Supporting this requires special cases in the code, and IMO just perpetuates confusion about the difference between the integer 0 and the function const0.
  2. We now have "try to evaluate this as a real number, and return NaN (or skip it) if we can't" code in at least five places. Should this be made consistent, or (better yet) factored out? I've made it so that the doctests pass, but in some places we catch only e.g. TypeError, while in others we check for a nice long list.
  3. This needs a careful review, since it can change the appearance of plots. There was some other (now fixed) bug involved as you can see from the commit list, but for example, this fix accidentally broke list_plot.
  4. I don't really like using try/except in fast loops. Is there a better way to fix the doctest failures that the first commit caused and that the later ones fix?

comment:12 in reply to: ↑ 11 Changed 2 months ago by mkoeppe

Replying to mjo:

  1. Do we really want to support passing (for example) the integer zero as a function to be plotted? We have doctests that check things like plot(ZZ(0), x, 0, 1). Supporting this requires special cases in the code, and IMO just perpetuates confusion about the difference between the integer 0 and the function const0.

Allowing non-callables as input for plot is, I think, important to keep. Convenience for casual use is key. Given that plot takes variable names as part of its input, it's clear that there is some implicit construction of a symbolic or numerical function happening.

comment:13 in reply to: ↑ 11 ; follow-up: Changed 2 months ago by mkoeppe

Replying to mjo:

This is an "easy fix" but causes a cascade of other issues because removing domain=float turns a lot of things that used to be Nan into TypeError, ValueError, ZeroDivisionError, etc.

Sorry, where do you remove domain=float?

comment:14 in reply to: ↑ 13 Changed 2 months ago by mjo

Replying to mkoeppe:

Replying to mjo:

This is an "easy fix" but causes a cascade of other issues because removing domain=float turns a lot of things that used to be Nan into TypeError, ValueError, ZeroDivisionError, etc.

Sorry, where do you remove domain=float?

The first thing fast_float() tries is to use fast_callable() with domain=float:

try:
    return fast_callable(f, vars=vars, domain=float,
                         expect_one_var=expect_one_var)

We're now using fast_callable() directly, and omit the domain.

comment:15 follow-up: Changed 2 months ago by mkoeppe

Thanks.

I have not looked at the details, nor done any timing. But it seems to me that a better solution may be to try domain=float first, then domain=complex (or whatever is needed to make fast_callable use interp_cdf); and only in the end fall back to the general interpreter.

  1. We now have "try to evaluate this as a real number, and return NaN (or skip it) if we can't" code in at least five places. Should this be made consistent, or (better yet) factored out?

Yes, I guess there would be value in a version of the general interpreter that catches errors and replaces them by NaN returns

Last edited 2 months ago by mkoeppe (previous) (diff)

comment:16 Changed 2 months ago by mkoeppe

By the way, I just came across https://trac.sagemath.org/attachment/ticket/5572/improve_fast_callable.6.patch, which seems to have done a lot of what we're doing again now...

comment:17 in reply to: ↑ 15 Changed 2 months ago by mjo

Replying to mkoeppe:

Thanks.

I have not looked at the details, nor done any timing. But it seems to me that a better solution may be to try domain=float first, then domain=complex (or whatever is needed to make fast_callable use interp_cdf); and only in the end fall back to the general interpreter.

This is tough because there are so many plotting interfaces, and the fast-callables aren't created near the point of failure. In this case, setup_for_eval_on_grid() is preprocessing some of the plotting args and returning them to some other function that orchestrates the plotting. Then the failure occurs in another function that actually computes the plot points.

Some major refactoring would be needed to do this all more intelligently. (I think we would also need to generate some kind of custom exception to avoid wrapping ten layers of code in one big except TypeError block.)

comment:18 follow-up: Changed 2 months ago by mkoeppe

  • Priority changed from major to minor
  • Status changed from needs_review to needs_work

OK, thanks for the explanation. Then I would suggest that we declare the issue of this ticket as "minor" and set it aside for later consideration.

comment:19 Changed 2 months ago by git

  • Commit changed from a773fd3919b69e575d94cdf60135c419b5e08f24 to 3168bdc46b2b9cc16bb3bcc9333c95739e0e4791

Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:

006dd61Trac #8450: support "output" wrapper around fast_callable().
a746166Trac #8450: use fast_callable() in setup_for_eval_on_grid().
3168bdcTrac #8450: update the output from a few plotting doctests.

comment:20 in reply to: ↑ 18 Changed 2 months ago by mjo

  • Status changed from needs_work to needs_review

Replying to mkoeppe:

OK, thanks for the explanation. Then I would suggest that we declare the issue of this ticket as "minor" and set it aside for later consideration.

I actually hold a grudge against this bug from back when I was generating plots programmatically and the random appearance of 0.000000000001*I would kill the whole thing after an hour, so I haven't given up.

I tried a few more things, described in the first of the recent commit messages, and discovered problems with all of them. I did however eventually find a solution that is both unobtrusive and fast: wrap fast-callable evaluation in a class that can force the result to a specific output type. With domain=<something complex> and output=<something real>, we return NaN unless the imaginary part is trivial, in which case we return the real part. In all other situations we simply try to coerce the output to the given type.

This allows us to plot over CDF but convert the result to float (possibly NaN) only on the way out, without having to hack that conversion logic into a hundred plotting functions. It's probably as fast as any solution could be; complex (two-float) operations necessarily take longer than plain float computations, but the difference is small in my non-rigorous tests, less than 10%.

An output wrapper is not used unless it's needed, so this doesn't affect any other uses of fast_callable() in the tree.

comment:21 Changed 2 months ago by mkoeppe

  • Cc egourgoulhon gh-jungmath added

comment:22 follow-up: Changed 2 months ago by mkoeppe

This looks like a fine approach. I won't have time before mid August to look at it in detail. Just a quick comment:

+        # This tolerance was not chosen for any particular reason.
+        if abs(ipart) < 1e-8:
+            return self._output(rpart)
+        else:
+            return self._output("nan")

Probably best to make these tolerances user-configurable (as an optional init argument of OutputWrapper? and optional argument of fast_callable).

comment:23 in reply to: ↑ 22 Changed 2 months ago by mjo

  • Status changed from needs_review to needs_work

Replying to mkoeppe:

Probably best to make these tolerances user-configurable (as an optional init argument of OutputWrapper? and optional argument of fast_callable).

I'm glad you brought this up. I swept it under the rug when I had fifteen things on my internal problem stack, but it's a code smell: using a fixed tolerance isn't right, but adding the argument to fast_callable() and trying to explain how it will be used exposes too many implementation details.

I think I can simplify things even further.

comment:24 Changed 2 months ago by git

  • Commit changed from 3168bdc46b2b9cc16bb3bcc9333c95739e0e4791 to 2b88457367c28b5315e60d8223d0d79aa750f474

Branch pushed to git repo; I updated commit sha1. This was a forced push. Last 10 new commits:

58fdf1esrc/sage/plot/plot3d/parametric_surface.pyx: Remove handling of FastDoubleFunc
ec44c61src/sage/ext/fast_eval.pyx: Reduce/update documentation
4c93ae0Trac #32234: use fast_callable() for symbolic find_root().
56e9997Trac #32234: don't import fast_callable() from sage.ext.fast_eval.
18593ffTrac #32234: replace fast_float_arg() in plot3d.
5fc69caTrac #32234: fix symbolic find_local_minimum() docs.
cdef097Trac #32234: use domain=float for fast callables used by numpy.
95a887fTrac #8450: new fast-callable output wrapper for plotting.
f696662Trac #8450: use fast_callable() in setup_for_eval_on_grid().
2b88457Trac #8450: update the output from a few plotting doctests.

comment:25 Changed 2 months ago by mjo

  • Status changed from needs_work to needs_review

Much better now. There's no need to solve the output-conversion problem more generally; we need it only for plotting, and creating the wrapper under sage.plot.misc eliminates all of the problems with the previous iteration.

comment:26 Changed 2 months ago by git

  • Commit changed from 2b88457367c28b5315e60d8223d0d79aa750f474 to d8e79086453876c2f1c43931590ffc94801fcaa1

Branch pushed to git repo; I updated commit sha1. This was a forced push. Last 10 new commits:

e6f88d6src/sage/ext/fast_eval.pyx: Reduce/update documentation
1353bafTrac #32234: use fast_callable() for symbolic find_root().
c56fa05Trac #32234: don't import fast_callable() from sage.ext.fast_eval.
0fcec8eTrac #32234: replace fast_float_arg() in plot3d.
1479666Trac #32234: fix symbolic find_local_minimum() docs.
19aa975Trac #32234: use domain=float for fast callables used by numpy.
df062b3Trac #32234: replace two trivial fast-callables in plot3d().
411d34dTrac #8450: new fast-callable output wrapper for plotting.
6600158Trac #8450: use fast_callable() in setup_for_eval_on_grid().
d8e7908Trac #8450: update the output from a few plotting doctests.

comment:27 follow-up: Changed 6 weeks ago by mkoeppe

This is looking nice. Similar to what I said in comment:22, I think the imaginary tolerance should be customizable on the level of the plot functions -- in the same way that options such as detect_poles, adaptive_tolerance are offered.

comment:28 in reply to: ↑ 27 Changed 5 weeks ago by mjo

Replying to mkoeppe:

This is looking nice. Similar to what I said in comment:22, I think the imaginary tolerance should be customizable on the level of the plot functions -- in the same way that options such as detect_poles, adaptive_tolerance are offered.

Those only work for the plot() function, and not for anything else that calls setup_for_eval_on_grid(). But once #32234 gets reviewed, I guess I wouldn't mind adding imaginary_tolerance to plot.options. In that case, we should probably go back later and try to make the options to the various plotting functions consistent (in another ticket).

Last edited 5 weeks ago by slelievre (previous) (diff)

comment:29 Changed 10 days ago by mkoeppe

  • Branch changed from u/mjo/ticket/8450 to u/mkoeppe/ticket/8450
Note: See TracTickets for help on using tickets.