Make plot API to work directly with axis

[sylee957]

I think that a lot of plotting API can be redone by making functions that take matplotlib Axes and operates in-place to Axes to append the data points and options to work with existing matplotlib instsance, rather than building the plot from scratch.

And we may need stuff more corresponding to the matplotlib object-oriented API than pyplot API, https://matplotlib.org/api/index.html#the-object-oriented-api.

However, this can still be tricky to supersede all the options that matplotlib have. For example, plotting API can often have parametric color specifications. And if that is available in some plotting API, it should have match with sympy functions or geometric primitives.

from sympy import *
import matplotlib.pyplot as plt
import numpy as np

x = Symbol('x')

def mpl_axes_plot(axes, expr, xrange, **kwargs):
    x, start, stop = xrange
    func = lambdify(x, expr)
    x = np.linspace(start, stop, 1000)
    y = func(x)
    return axes.plot(x, y, **kwargs)


fig, ax = plt.subplots()
mpl_axes_plot(ax, sin(x), (x, -5, 5))
mpl_axes_plot(ax, cos(x), (x, -5, 5))

plt.show()

from sympy import *
import matplotlib.pyplot as plt
import numpy as np

x = Symbol('x')

def mpl_axes_step(axes, expr, xrange, **kwargs):
    x, start, stop, step = xrange
    func = lambdify(x, expr)
    x = np.arange(start, stop, step)
    y = func(x)
    return axes.step(x, y, **kwargs)

fig, ax = plt.subplots()
mpl_axes_step(ax, sin(x), (x, -5, 5, 0.5))
mpl_axes_step(ax, cos(x), (x, -5, 5, 0.5))

plt.show()

from sympy import *
import matplotlib.pyplot as plt
import numpy as np

x = Symbol('x')

def mpl_axes_fill_between(axes, expr1, expr2, xrange, **kwargs):
    x, start, stop = xrange
    func1 = lambdify(x, expr1)
    func2 = lambdify(x, expr2)
    x = np.linspace(start, stop, 1000)
    y1 = func1(x)
    y2 = func2(x)
    return axes.fill_between(x, y1, y2, **kwargs)

fig, ax = plt.subplots()
mpl_axes_fill_between(ax, sin(x), cos(x), (x, -5, 5))

plt.show()

from sympy import *
import matplotlib.pyplot as plt
import numpy as np

x = Symbol('x')

def mpl_axes_stem(axes, expr, xrange, **kwargs):
    x, start, stop, step = xrange
    func = lambdify(x, expr)
    x = np.arange(start, stop, step)
    y = func(x)
    return axes.stem(x, y, use_line_collection=True, **kwargs)

fig, ax = plt.subplots()
mpl_axes_stem(ax, sinc(x), (x, -10, 10, 0.5))

plt.show()

from sympy import *
import matplotlib.pyplot as plt
import numpy as np

x = Symbol('x')
y = Symbol('y')

def mpl_axes_quiver(axes, vect, xrange, yrange, **kwargs):
    expr1, expr2 = vect
    x, xstart, xstop = xrange
    y, ystart, ystop = yrange
    func1 = lambdify([x, y], expr1)
    func2 = lambdify([x, y], expr2)
    
    x = np.linspace(xstart, xstop, 10)
    y = np.linspace(ystart, ystop, 10)
    xx, yy = np.meshgrid(x, y)
    u = func1(xx, yy)
    v = func2(xx, yy)
    axes.quiver(xx, yy, u, v, **kwargs)

fig, ax = plt.subplots()
mpl_axes_quiver(ax, [-x / sqrt(x**2 + y**2 + 10), y / sqrt(x**2 + y**2 + 10)], (x, -5, 5), (y, -5, 5))

plt.show()

from sympy import *
import matplotlib.pyplot as plt
import numpy as np

x = Symbol('x')
y = Symbol('y')

def mpl_axes_streamplot(axes, vect, xrange, yrange, **kwargs):
    expr1, expr2 = vect
    x, xstart, xstop = xrange
    y, ystart, ystop = yrange
    func1 = lambdify([x, y], expr1)
    func2 = lambdify([x, y], expr2)
    
    x = np.linspace(xstart, xstop, 20)
    y = np.linspace(ystart, ystop, 20)
    xx, yy = np.meshgrid(x, y)
    u = func1(xx, yy)
    v = func2(xx, yy)
    axes.streamplot(xx, yy, u, v, **kwargs)

fig, ax = plt.subplots()
mpl_axes_streamplot(ax, [-x / sqrt(x**2 + y**2 + 10), y / sqrt(x**2 + y**2 + 10)], (x, -5, 5), (y, -5, 5))

plt.show()

from sympy import *
import matplotlib.pyplot as plt
import numpy as np

x = Symbol('x')
y = Symbol('y')

def mpl_axes_contour(axes, expr, xrange, yrange, **kwargs):
    x, xstart, xstop = xrange
    y, ystart, ystop = yrange
    func = lambdify([x, y], expr)
    
    x = np.linspace(xstart, xstop, 100)
    y = np.linspace(ystart, ystop, 100)
    xx, yy = np.meshgrid(x, y)
    z = func(xx, yy)
    axes.contour(xx, yy, z, **kwargs)

fig, ax = plt.subplots()
mpl_axes_contour(ax, sin(x**2 + y**2) / (x**2 + y**2), (x, -5, 5), (y, -5, 5))

plt.show()

from sympy import *
import matplotlib.pyplot as plt
import numpy as np

x = Symbol('x')
y = Symbol('y')

def mpl_axes_matshow(axes, expr, xrange, yrange, **kwargs):
    x, xstart, xstop = xrange
    y, ystart, ystop = yrange
    func = lambdify([x, y], expr)
    
    x = np.linspace(xstart, xstop, 100)
    y = np.linspace(ystart, ystop, 100)
    xx, yy = np.meshgrid(x, y)
    z = func(xx, yy)
    axes.matshow(z, **kwargs)

fig, ax = plt.subplots()
mpl_axes_matshow(ax, sin(x**2 + y**2) / (x**2 + y**2), (x, -5, 5), (y, -5, 5))

plt.show()

from mpl_toolkits.mplot3d import axes3d
import matplotlib.pyplot as plt
import numpy as np
from sympy import *

fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')

x = Symbol('x')
y = Symbol('y')

def mpl_axes3d_plot_wireframe(axes, expr, xrange, yrange, **kwargs):
    x, xstart, xstop = xrange
    y, ystart, ystop = yrange
    func = lambdify([x, y], expr)
    
    x = np.linspace(xstart, xstop, 20)
    y = np.linspace(ystart, ystop, 20)
    xx, yy = np.meshgrid(x, y)
    z = func(xx, yy)
    axes.plot_wireframe(xx, yy, z, **kwargs)

mpl_axes3d_plot_wireframe(ax, sin(sqrt(x**2 + y**2)), (x, -5, 5), (y, -5, 5))

plt.show()

from mpl_toolkits.mplot3d import axes3d
import matplotlib.pyplot as plt
import numpy as np
from sympy import *

fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')

x = Symbol('x')
y = Symbol('y')

def mpl_axes3d_plot_surface(axes, expr, xrange, yrange, **kwargs):
    x, xstart, xstop = xrange
    y, ystart, ystop = yrange
    func = lambdify([x, y], expr)
    
    x = np.linspace(xstart, xstop, 100)
    y = np.linspace(ystart, ystop, 100)
    xx, yy = np.meshgrid(x, y)
    z = func(xx, yy)
    axes.plot_surface(xx, yy, z, **kwargs)

mpl_axes3d_plot_surface(ax, sin(sqrt(x**2 + y**2)), (x, -5, 5), (y, -5, 5))

plt.show()

from mpl_toolkits.mplot3d import axes3d
import matplotlib.pyplot as plt
import numpy as np
from sympy import *
​
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
​
x = Symbol('x')
y = Symbol('y')
​
def mpl_axes3d_scatter(axes, expr, xrange, yrange, **kwargs):
    x, xstart, xstop = xrange
    y, ystart, ystop = yrange
    func = lambdify([x, y], expr)
    
    x = np.linspace(xstart, xstop, 20)
    y = np.linspace(ystart, ystop, 20)
    xx, yy = np.meshgrid(x, y)
    z = func(xx, yy)
    axes.scatter(xx, yy, z, **kwargs)
​
mpl_axes3d_scatter(ax, sin(sqrt(x**2 + y**2)), (x, -5, 5), (y, -5, 5))
​
plt.show()

from mpl_toolkits.mplot3d import axes3d
import matplotlib.pyplot as plt
import numpy as np
from sympy import *

fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')

x = Symbol('x')
y = Symbol('y')

def mpl_axes3d_plot(axes, vect, xrange, **kwargs):
    x, xstart, xstop = xrange
    expr1, expr2 = vect
    func1 = lambdify(x, expr1)
    func2 = lambdify(x, expr2)
    
    x = np.linspace(xstart, xstop, 100)
    y = func1(x)
    z = func2(x)
    axes.plot(x, y, z, **kwargs)

mpl_axes3d_plot(ax, (cos(x), sin(x)), (x, -5, 5))

plt.show()

from mpl_toolkits.mplot3d import axes3d
import matplotlib.pyplot as plt
import numpy as np
from sympy import *

fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')

x = Symbol('x')
y = Symbol('y')
z = Symbol('z')

def mpl_axes3d_quiver(axes, vect, xrange, yrange, zrange, **kwargs):
    x, xstart, xstop = xrange
    y, ystart, ystop = yrange
    z, zstart, zstop = zrange
    
    expr1, expr2, expr3 = vect
    func1 = lambdify([x, y, z], expr1)
    func2 = lambdify([x, y, z], expr2)
    func3 = lambdify([x, y, z], expr3)
    
    x = np.linspace(xstart, xstop, 5)
    y = np.linspace(ystart, ystop, 5)
    z = np.linspace(zstart, zstop, 5)
    xx, yy, zz = np.meshgrid(x, y, z)
    u = func1(xx, yy, zz)
    v = func2(xx, yy, zz)
    w = func3(xx, yy, zz)
    axes.quiver(xx, yy, zz, u, v, w, normalize=True, **kwargs)

mpl_axes3d_quiver(ax, [-y, x, 0], (x, -5, 5), (y, -5, 5), (z, -5, 5))

plt.show()

[oscarbenjamin]

I think the best thing to do would be to make functions that do the boilerplate for lambdify etc in the examples above and produce arguments that are in the suitable form for passing to matplotlib. Then we can just document to users how to pass the results on to matplotlib (or other libraries) and we don't need any specific code for each type of plot except where we have a convenience function like sympy.plot.

[moorepants]

These are nice @sylee957. As a start, all of these could be added to the documentation as is.

I agree that operating on mpl's Axes objects a good approach (and is standard in most of libraries).

As for the future of sympy.plot(), someone needs to lay out a plan to morph it into something better. I like the idea of supporting different backends and think that should remain, but maybe in a different form. It'd be nice to be able to do sympy.plot(expr, backend='plotly'), sympy.plot(expr, backend='bokeh'), sympy.plot(expr, backend='altair'), sympy.plot(expr, backend='matplotlib'), etc.

[moorepants]

But, maybe, as @oscarbenjamin mentions, SymPy should only output numerical data for given expressions that is in a form for a particular plotting library to consume. Then we'd not have any plotting code present, only data generation code.

[Smit-create]

I like the idea of supporting different backends and think that should remain, but maybe in a different form

I think this could be a good idea to have it. We can make it simpler by just passing the expr to sympy.plot, then using backend option in lambdify to convert it to the required form. This will also lead to a clearer and simple code base of plotting.

[sylee957]

I also agree that having the old plot module can be good because it is more easier to use with some restriction of advanced features. But the features can operate on top of the low-level axex manipulation APIs

But the problem of the backend approach is that it is difficult to share options between plots. For example, even if we only have TextBackend and matplotlib, but they share almost no options between, and currently even the data arrays are not shared efficiently.

Of course, if the original author of the plotting module had intended to develop more advanced features like drag-able plot or animation, there can be a good point of keeping sympy expression in a plot object.

[oscarbenjamin]

The backend approach is a terrible idea as it means that sympy's plotting module has to grow all the API for all of the things that the backends can do. Users and contributors end up wasting their time trying to do things that are already implemented in matplotlib.

@sylee957 I think your suggestion here is good and should be used in the current plotting API. Also I think that the documentation and tutorial for plotting should make it very clear that sympy's plotting is just using matplotlib and should show users how to combine sympy with the matplotlib API before mentioning the convenience routines like sympy.plot. I think users would also benefit from examples showing how to use other libraries e.g. mayavi etc.

What I don't think we should do is extend the sympy plotting API apart from providing convenient utilities that users can use to prepare the data for input to other plotting libraries.

[moorepants]

The backend approach is a terrible idea

I feel like this kind of statement is certainly what a developer of SymPy might think, not so sure about users of SymPy though.

[moorepants]

But the problem of the backend approach is that it is difficult to share options between plots.

I don't think we'd need to give access to many options. The shared backend approach could work in a very constrained way, but provide the user with an underlying object or data structure that can be taken into a specific plotting library when they need something outside of our constraints.

For example, if the matplotlib axes object(s) were accessible from our Plot() object in a very simple way, it would open up a lot. You can could create the base MPL plot with a simple command and then operate on the axes with MPL. Seems like the best of both worlds.

[oscarbenjamin]

For example, if the matplotlib axes object(s) were accessible from our Plot() object in a very simple way, it would open up a lot. You can could create the base MPL plot with a simple command and then operate on the axes with MPL. Seems like the best of both worlds.

I think this is a good idea. It does call into question the purpose of the Plot instance in the first place though.

The backend approach is a terrible idea

I feel like this kind of statement is certainly what a developer of SymPy might think, not so sure about users of SymPy though.

I think that we can make things easier and better for users by making it clear how to use sympy with matplotlib. The OP idea from @sylee957 as well as yours above about providing access to the axes object would both help with that. The most useful thing for users is just being clear that sympy's plotting functionality uses matplotlib which already has many features and plenty of documentation and both of your ideas help to do that.

Building a duplicate API on the other hand does not benefit users since the API will always be incomplete and inconsistent with the API that it wraps and its relationship to matplotlib will remain unclear to users.

I also think that mission creep is an important problem that sympy needs to worry about. The sympy codebase is already enormous and something like plotting can suck a lot of contributor and user time reimplementing features that already exist in related libraries when for most users all that's really needed is a few utility functions and better documentation.

[moorepants]

I think this is a good idea. It does call into question the purpose of the Plot instance in the first place though.

I think it is there to take away the boiler plate lambdify and plot setup for typical plots. The other use it currently seems to have is handling some of the numerical issues associated with plotting complex functions (singularities for example). Some functions require non-trival amounts of NumPy and MPL code to get correct. Overall I think there is a place for a symbolic based plot command that "just works". Other symbolic libraries have plot commands. For example, it's nice in Sage that plotting is simple. One use case is using SymPy to teach math. It is a burden to have to explain numpy and matplotlib when you just want a student to see a plot of their math expression. I use SymPy's plot function for this often.

I think that we can make things easier and better for users by making it clear how to use sympy with matplotlib.

Agreed.

I also think that mission creep is an important problem that sympy needs to worry about.

That calls for a separate thread 😄

[oscarbenjamin]

I think this is a good idea. It does call into question the purpose of the Plot instance in the first place though.

I think it is there to take away the boiler plate lambdify and plot setup for typical plots.

That's what the plot function does. What I mean is what is the purpose of the Plot instance itself? I think that the MPL axes instance is more useful.

The other use it currently seems to have is handling some of the numerical issues associated with plotting complex functions (singularities for example). Some functions require non-trival amounts of NumPy and MPL code to get correct.

Where possible that should be improved in other parts of sympy to make it easier rather than being part of the plotting module.

Overall I think there is a place for a symbolic based plot command that "just works". Other symbolic libraries have plot commands.

I agree that it is good to have a simple plot command. I just want it to remain simple and the documentation should make it clear that it is deliberately simple and is just a wrapper for another library so that it is clear to everyone how to access the more advanced features that are already available from the other library. There is definite benefit for sympy to implement domain-specific plotting features that are not implemented somewhere else but we need to be careful not to end up in a situation of reimplementing things that already exist elsewhere.

An example of something that is nontrivial and makes sense for sympy but is not implemented directly in matplotlib is domain colouring:
https://en.wikipedia.org/wiki/Domain_coloring

[sylee957]

Although the frontend-backend approach used in Plot hadn't proved much usefulness, I don't think that there is much evidence that the design was fundamentally flawed.
The issues like #19051 are unrelated to this, and neither your comments about the documentation issue is significantly relevant for this because I think that it can be resolved by improving the docs for what we currently have.

And I'm afraid that changing the return types of plot is not trivial and can cause a lot of deprecation issues, and it can confuse users in other way if we cannot provide much reason about how such features were just experimental and it had to be deprecated after some trials.

I think that frontend-backend approach is valid in the sense that there are some distinguishable 'major' features from minor options for plots.
For example, everyone would design plots with major options like title, labels, ticks, colors, etc.
So having Plot.title, Plot.xlabel as a unified control for such options regardless of the implementation details of matplotlib, plotly, ... can still be a good idea.

The concern only happens when too much trivial or implementation-specific options are added for plots, and if we are unable to control how such keywords mangles each other, but I think that the contemporary CAS softwares can be reviewed to investigate whether addition of some plotting options are really necessary or not.