Skip to content
This repository has been archived by the owner on Jan 30, 2023. It is now read-only.

Commit

Permalink
Trac #24536: factor out a FastCallableFloatWrapper class from plots.
Browse files Browse the repository at this point in the history
We have a few numerical methods that run into trouble in certain
situations:

  1. When the objective function is symbolic and contains complex
     sub-expressions, the real fast_callable() wrapper can't handle
     the (intermediate) complex results.

  2. When the function evaluations can technically be complex, but the
     imaginary parts are essentially numerical noise.

Our plotting infrastructure had the same problems, and we worked
around it by (a) using complex fast_callable() wrappers instead of
real ones, and (b) adding *another* wrapper around each fast-callable.
The so-called FastCallablePlotWrapper ignores imaginary noise, but
returns NaN if the result is truly complex -- the end result being
that complex results are "skipped" while plotting.

A similar solution should work for other numerical methods, but we
don't want to ignore complex results if we're trying to e.g. minimize
a function that is supposed to return real numbers. So instead of
using the existing FastCallablePlotWrapper, this commit factors out a
superclass called FastCallableFloatWrapper that works just like the
plot wrapper except that it raises an error when it sees a complex
result from the underlying fast_callable().
  • Loading branch information
orlitzky committed Mar 8, 2022
1 parent b21ca55 commit beca559
Show file tree
Hide file tree
Showing 2 changed files with 158 additions and 100 deletions.
145 changes: 145 additions & 0 deletions src/sage/ext/fast_callable.pyx
Expand Up @@ -2444,3 +2444,148 @@ cdef class Wrapper:
if isinstance(op, tuple) and op[0] == 'py_call':
py_calls.append(op[1])
return py_calls


class FastCallableFloatWrapper:
r"""
A class to alter the return types of the fast-callable functions.
When applying numerical routines (including plotting) to symbolic
expressions and functions, we generally first convert them to a
faster form with :func:`fast_callable`. That function takes a
``domain`` parameter that forces the end (and all intermediate)
results of evaluation to a specific type. Though usually always
want the end result to be of type ``float``, correctly choosing
the ``domain`` presents some problems:
* ``float`` is a bad choice because it's common for real
functions to have complex terms in them. Moreover precision
issues can produce terms like ``1.0 + 1e-12*I`` that are hard
to avoid if calling ``real()`` on everything is infeasible.
* ``complex`` has essentially the same problem as ``float``.
There are several symbolic functions like :func:`min_symbolic`,
:func:`max_symbolic`, and :func:`floor` that are unable to
operate on complex numbers.
* ``None`` leaves the types of the inputs/outputs alone, but due
to the lack of a specialized interpreter, slows down evaluation
by an unacceptable amount.
* ``CDF`` has none of the other issues, because ``CDF`` has its
own specialized interpreter, a lexicographic ordering (for
min/max), and supports :func:`floor`. However, most numerical
functions cannot handle complex numbers, so using ``CDF``
would require us to wrap every evaluation in a
``CDF``-to-``float`` conversion routine. That would slow
things down less than a domain of ``None`` would, but is
unattractive mainly because of how invasive it would be to
"fix" the output everywhere.
Creating a new fast-callable interpreter that has different input
and output types solves most of the problems with a ``CDF``
domain, but :func:`fast_callable` and the interpreter classes in
:mod:`sage.ext.interpreters` are not really written with that in
mind. The ``domain`` parameter to :func:`fast_callable`, for
example, is expecting a single Sage ring that corresponds to one
interpreter. You can make it accept, for example, a string like
"CDF-to-float", but the hacks required to make that work feel
wrong.
Thus we arrive at this solution: a class to wrap the result of
:func:`fast_callable`. Whenever we need to support intermediate
complex terms in a numerical routine, we can set ``domain=CDF``
while creating its fast-callable incarnation, and then wrap the
result in this class. The ``__call__`` method of this class then
ensures that the ``CDF`` output is converted to a ``float`` if
its imaginary part is within an acceptable tolerance.
EXAMPLES:
The ``float`` incarnation of "not a number" is returned instead
of an error being thrown if the answer is complex::
sage: from sage.ext.fast_callable import FastCallableFloatWrapper
sage: f = sqrt(x)
sage: ff = fast_callable(f, vars=[x], domain=CDF)
sage: fff = FastCallableFloatWrapper(ff, imag_tol=1e-8)
sage: fff(1)
1.0
sage: fff(-1)
Traceback (most recent call last):
...
ValueError: complex fast-callable function result
1.0*I for arguments (-1,)
"""
def __init__(self, ff, imag_tol):
r"""
Construct a ``FastCallableFloatWrapper``.
INPUT:
- ``ff`` -- a fast-callable wrapper over ``CDF``; an instance of
:class:`sage.ext.interpreters.Wrapper_cdf`, usually constructed
with :func:`fast_callable`.
- ``imag_tol`` -- float; how big of an imaginary part we're willing
to ignore before raising an error.
OUTPUT:
An instance of ``FastCallableFloatWrapper`` that can be
called just like ``ff``, but that always returns a ``float``
if no error is raised. A ``ValueError`` is raised if the
imaginary part of the result exceeds ``imag_tol``.
EXAMPLES:
The wrapper will ignore an imaginary part smaller in magnitude
than ``imag_tol``::
sage: from sage.ext.fast_callable import FastCallableFloatWrapper
sage: f = x
sage: ff = fast_callable(f, vars=[x], domain=CDF)
sage: fff = FastCallableFloatWrapper(ff, imag_tol=1e-8)
sage: fff(I*1e-9)
0.0
sage: fff = FastCallableFloatWrapper(ff, imag_tol=1e-12)
sage: fff(I*1e-9)
Traceback (most recent call last):
...
ValueError: complex fast-callable function result 1e-09*I for
arguments (1.00000000000000e-9*I,)
"""
self._ff = ff
self._imag_tol = imag_tol

def __call__(self, *args):
r"""
Evaluate the underlying fast-callable and convert the result to
``float``.
TESTS:
Evaluation either returns a ``float``, or raises a
``ValueError``::
sage: from sage.ext.fast_callable import FastCallableFloatWrapper
sage: f = x
sage: ff = fast_callable(f, vars=[x], domain=CDF)
sage: fff = FastCallableFloatWrapper(ff, imag_tol=0.1)
sage: try:
....: result = fff(CDF.random_element())
....: except ValueError:
....: result = float(0)
sage: type(result) is float
True
"""
z = self._ff(*args)

if abs(z.imag()) < self._imag_tol:
return float(z.real())
else:
raise ValueError(f"complex fast-callable function result {z} " +
f"for arguments {args}")
113 changes: 13 additions & 100 deletions src/sage/plot/misc.py
Expand Up @@ -13,6 +13,9 @@
# http://www.gnu.org/licenses/
#*****************************************************************************

from sage.ext.fast_callable import FastCallableFloatWrapper


def setup_for_eval_on_grid(funcs,
ranges,
plot_points=None,
Expand Down Expand Up @@ -444,62 +447,13 @@ def get_matplotlib_linestyle(linestyle, return_type):
(linestyle))


class FastCallablePlotWrapper:
class FastCallablePlotWrapper(FastCallableFloatWrapper):
r"""
A class to alter the return types of the fast-callable functions
used during plotting.
When plotting symbolic expressions and functions, we generally
first convert them to a faster form with :func:`fast_callable`.
That function takes a ``domain`` parameter that forces the end
(and all intermediate) results of evaluation to a specific type.
Though we always want the end result to be of type ``float``,
correctly choosing the ``domain`` presents some problems:
* ``float`` is a bad choice because it's common for real
functions to have complex terms in them. Moreover, when one is
generating plots programmatically, precision issues can
produce terms like ``1.0 + 1e-12*I`` that are hard to avoid if
calling ``real()`` on everything is infeasible.
* ``complex`` has essentially the same problem as ``float``.
There are several symbolic functions like :func:`min_symbolic`,
:func:`max_symbolic`, and :func:`floor` that are unable to
operate on complex numbers.
* ``None`` leaves the types of the inputs/outputs alone, but due
to the lack of a specialized interpreter, slows down plotting
by an unacceptable amount.
* ``CDF`` has none of the other issues, because ``CDF`` has its
own specialized interpreter, a lexicographic ordering (for
min/max), and supports :func:`floor`. However, none of the
plotting functions can handle complex numbers, so using
``CDF`` would require us to wrap every evaluation in a
``CDF``-to-``float`` conversion routine within the plotting
infrastructure. This slows things down less than a domain of
``None`` does, but is unattractive mainly because of how
invasive it would be to "fix" the output everywhere.
Creating a new fast-callable interpreter that has different input
and output types solves most of the problems with a ``CDF``
domain, but :func:`fast_callable` and the interpreter classes in
:mod:`sage.ext.interpreters` are not really written with that in
mind. The ``domain`` parameter to :func:`fast_callable`, for
example, is expecting a single Sage ring that corresponds to one
interpreter. You can make it accept, for example, a string like
"CDF-to-float", but the hacks required to make that work feel
wrong.
Thus we arrive at this solution: a class to wrap the result of
:func:`fast_callable`. Whenever we need to support intermediate
complex terms in a plot function, we can set ``domain=CDF`` while
creating its fast-callable incarnation, and then wrap the result in
this class. The ``__call__`` method of this class then ensures
that the ``CDF`` output is converted to a ``float``. Since
plotting tries to ignore unplottable points, this job is easier
than it would be in a more general context: we simply return
``nan`` whenever the result has a nontrivial imaginary part.
A fast-callable wrapper for plotting that returns ``nan`` instead
of raising an error whenever the imaginary tolerance is exceeded.
A detailed rationale for this can be found in the superclass
documentation.
EXAMPLES:
Expand All @@ -514,46 +468,7 @@ class FastCallablePlotWrapper:
1.0
sage: fff(-1)
nan
"""
def __init__(self, ff, imag_tol):
r"""
Construct a ``FastCallablePlotWrapper``.
INPUT:
- ``ff`` -- a fast-callable wrapper over ``CDF``; an instance of
:class:`sage.ext.interpreters.Wrapper_cdf`, usually constructed
with :func:`fast_callable`.
- ``imag_tol`` -- float; how big of an imaginary part we're willing
to ignore before returning ``nan``.
OUTPUT:
An instance of ``FastCallablePlotWrapper`` that can be called
just like ``ff``, but that always returns a ``float``, even if
it is ``nan``.
EXAMPLES:
The wrapper will ignore an imaginary part smaller in magnitude
than ``imag_tol``::
sage: from sage.plot.misc import FastCallablePlotWrapper
sage: f = x
sage: ff = fast_callable(f, vars=[x], domain=CDF)
sage: fff = FastCallablePlotWrapper(ff, imag_tol=1e-8)
sage: fff(I*1e-9)
0.0
sage: fff = FastCallablePlotWrapper(ff, imag_tol=1e-12)
sage: fff(I*1e-9)
nan
"""
self._ff = ff
self._imag_tol = imag_tol

def __call__(self, *args):
r"""
Evaluate the underlying fast-callable and convert the result to
Expand All @@ -566,14 +481,12 @@ def __call__(self, *args):
sage: from sage.plot.misc import FastCallablePlotWrapper
sage: f = x
sage: ff = fast_callable(f, vars=[x], domain=CDF)
sage: fff = FastCallablePlotWrapper(ff, imag_tol=1e-8)
sage: fff = FastCallablePlotWrapper(ff, imag_tol=0.1)
sage: type(fff(CDF.random_element())) is float
True
"""
z = self._ff(*args)

if abs(z.imag()) < self._imag_tol:
return float(z.real())
else:
try:
return super().__call__(*args)
except ValueError:
return float("nan")

0 comments on commit beca559

Please sign in to comment.