Opened 8 years ago

Closed 2 years ago

#10775 closed enhancement (fixed)

streamline plots

Reported by: jason Owned by: jason, was
Priority: major Milestone: sage-7.4
Component: graphics Keywords:
Cc: kcrisman, alauve, jakobkroeker, paulmasson, jhonrubia6 Merged in:
Authors: Paul Masson Reviewers: Travis Scrimshaw
Report Upstream: Completely fixed; Fix reported upstream Work issues:
Branch: 640cb6a (Commits) Commit: 640cb6aad1a6b6380adabcb2c4d24d10eca3e718
Dependencies: Stopgaps:

Description

Recently there was a thread on the matplotlib list where someone posted matplotlib code for a streamline plot. It looks very nice. If the code doesn't get included in matplotlib, we could still ask about including it in Sage to do streamline plots.

Mailing list message: http://permalink.gmane.org/gmane.comp.python.matplotlib.general/26362

Mailing list thread: http://comments.gmane.org/gmane.comp.python.matplotlib.general/26354

Example code: http://web.mit.edu/speth/Public/streamlines.py

Example plot: http://web.mit.edu/speth/Public/streamlines.png

Attachments (3)

streamlines.py (6.6 KB) - added by jason 8 years ago.
streamplot.py (9.5 KB) - added by jason 8 years ago.
DE Solutions (1).sws (2.7 KB) - added by jason 6 years ago.

Download all attachments as: .zip

Change History (43)

Changed 8 years ago by jason

Changed 8 years ago by jason

comment:2 Changed 8 years ago by jason

I've attached both code files so that we have them even if the above links go down for whatever reason.

comment:4 Changed 7 years ago by kcrisman

  • Report Upstream changed from N/A to Fixed upstream, but not in a stable release.

Yo, this is now in mpl!

I'm not sure how to figure out which release it's in, if it's in a release yet. It looks like it will probably be in 1.1.1, as 1.1.0 has been out for longer than this was in.

comment:5 Changed 7 years ago by jason

See also http://sage.cs.drake.edu/home/pub/60/ for another example.

Last edited 6 years ago by jason (previous) (diff)

Changed 6 years ago by jason

comment:6 Changed 6 years ago by jason

I attached the example I pointed out in the sage.cs.drake.edu server.

comment:7 Changed 6 years ago by kcrisman

#13693 is for the appropriate mpl upgrade for this.

comment:8 Changed 6 years ago by jdemeyer

  • Milestone changed from sage-5.11 to sage-5.12

comment:9 Changed 5 years ago by vbraun_spam

  • Milestone changed from sage-6.1 to sage-6.2

comment:10 Changed 5 years ago by vbraun_spam

  • Milestone changed from sage-6.2 to sage-6.3

comment:11 Changed 5 years ago by vbraun_spam

  • Milestone changed from sage-6.3 to sage-6.4

comment:12 Changed 4 years ago by kcrisman

Here are some upstream examples.


Here is another example from somewhere in this discussion, allowing phase plots. Not directly related but looks nice with such things.

# based on code from oddrobot, http://sagenb.org/home/pub/1532/
from sage.calculus.desolvers import desolve_system_rk4
x,y,t=var('x y t')
class DESolution:
    def __init__(self,system,time_range,initial,stepsize=0.05):
        
        self.tvar=time_range[0]
        self._times=srange(time_range[1],time_range[2],stepsize)
        self.vars=[v for v,_ in initial]
        self.dim=len(self.vars)
        self._soln=desolve_odeint(system, ics=[v for _,v in initial], times=self._times, dvars=self.vars, ivar=self.tvar)
        
    def phase_plot(self,vars=None,color='blue',**kwargs):
        # find which indices the specified variables are
        if vars is not None:
            vars_index=[self.vars.index(v) for v in vars]
        elif self.dim<=3:
            vars_index=range(self.dim)
        else:
            vars_index=range(2)
            
        p=line(self._soln[:,vars_index],color=color,**kwargs)
        # add an arrow head showing which way we are going around the phase line
        half=int(self._soln.shape[0]/2)
        p+=arrow(self._soln[half,vars_index], self._soln[half+1,vars_index],color=color)
        if len(vars_index)==2:
            p.axes_labels([str(self.vars[v]) for v in vars_index])
        return p
        
    def coordinates(self,colors=None,**kwargs):
        if colors is None:
            colors=rainbow(len(self.vars))
        p=Graphics()
        p+=sum(line(zip(self._times,self._soln[:,i]), color=colors[i], legend_label=str(self.vars[i]),**kwargs) for i in range(self.dim))
        return p

comment:13 Changed 2 years ago by tscrim

  • Cc alauve jakobkroeker paulmasson jhonrubia6 added
  • Milestone changed from sage-6.4 to sage-7.4
  • Report Upstream changed from Fixed upstream, but not in a stable release. to Completely fixed; Fix reported upstream

I was hoping we could resurrect this ticket as I am teaching a class on ODE's this semester and would like to plot streamlines on top of the slope fields. Does someone have some time to work on this? Alas, I don't know how much time I might have to do so.

So streamlines are included in matplotlib: http://matplotlib.org/examples/images_contours_and_fields/streamplot_demo_features.html

There is also another plot package I came across that is open-source python that does streamlines: https://plot.ly/python/streamline-plots/

comment:14 Changed 2 years ago by paulmasson

Since I spent too much time today reading graphics code, I can take a crack at this. It looks like a new class is needed, modeled on PlotField but calling the pyplot method streamplot instead of quiver. It shouldn't be too difficult but won't happen overnight.

Travis, if you see a way to avoid a new class, then let me know before I dive in tomorrow.

comment:15 Changed 2 years ago by paulmasson

To clarify: by new class I mean an individual file like scatter_plot.py, which calls the pyplot method scatter.

Last edited 2 years ago by paulmasson (previous) (diff)

comment:16 Changed 2 years ago by kcrisman

Yes, this is exactly how you should do it. You could even keep it in the same module as the vector field plots, you can see that e.g. contour and density plots are in the same file.

comment:17 Changed 2 years ago by tscrim

Thank you for taking point on this Paul!

I agree with Karl-Dieter, you should definitely create a new class. However, contour and density plots are in separate files. Although, I have no preference whether or not you put in a separate file; it fits very well with plot_field.py.

comment:18 Changed 2 years ago by paulmasson

  • Branch set to u/paulmasson/10775

comment:19 Changed 2 years ago by paulmasson

  • Authors set to Paul Masson
  • Commit set to 9576ba7abe47606fdf7667ea88e4bf8cf6dbefad
  • Status changed from new to needs_review

This was surprisingly easy to get working! I copied over plot_field.py and left variable names the same for consistency. Then it was a matter of rearranging how the data arrays were assembled and changing the pyplot function name.

Doctests and useful examples have been updated. I left the _allowed_options untouched because they might be relevant.


New commits:

852c719Initial commit
e6094a3Initial working commit
188b3c7Update doctests
4c2ff90Add streamline_plot
9576ba7Update documentation

comment:20 Changed 2 years ago by tscrim

Sorry to ask for a bit more. Could we add a way to input an initial point, allow input as a single function to treat it as a slope field, and an option to display the vector/slope field? If this is too much trouble, I am okay with positive reviewing this ticket as-is.

comment:21 follow-up: Changed 2 years ago by paulmasson

Travis, if by "initial point" you mean initial conditions, aren't those just part of the function arguments given to the plot? What more do you intend with this feature?

The idea of a convenience function that switches between streamline and vector/slope plots is interesting, but shouldn't be part of streamline_plot itself. That one should stand alone, partly because the streamplot function it uses behaves differently from the quiver function used by the others. Some of the examples from plot_field.py had to be left out of this new file because streamplot returned errors for them.

Would this convenience function be used regularly by others, or just you for one class? What's the rule of thumb for adding those sorts of features to Sage?

If you still want an additional convenience function added to this file, give me a name for it and its options.

comment:22 in reply to: ↑ 21 Changed 2 years ago by tscrim

Replying to paulmasson:

Travis, if by "initial point" you mean initial conditions, aren't those just part of the function arguments given to the plot?

Yes, but how I am reading the code, it looks like it does streamlines starting at each integer point.

What more do you intend with this feature?

Plot a solution for a given IVP.

The idea of a convenience function that switches between streamline and vector/slope plots is interesting, but shouldn't be part of streamline_plot itself. That one should stand alone, partly because the streamplot function it uses behaves differently from the quiver function used by the others. Some of the examples from plot_field.py had to be left out of this new file because streamplot returned errors for them.

If you look at the slope plot, it redirects to the vector plot with special options. We should be able to automate this by looking at the input, a tuple for a vector plot and otherwise consider it a slope plot.

Would this convenience function be used regularly by others, or just you for one class? What's the rule of thumb for adding those sorts of features to Sage?

I would be surprised if this did not get significant uptake in a variety of other classes. One of my first forays into Python was plotting a single streamline for the Lorenz attractor. However, there are a number of things in Sage that are likely only used by a few people, so the rule is it gets added if someone is willing to positively review it.

If you still want an additional convenience function added to this file, give me a name for it and its options.

Granted, I would eventually like a 3d version, but that can most certainly be a followup ticket.

I should have time next week to actually directly work on this (i.e. write code).

comment:23 Changed 2 years ago by git

  • Commit changed from 9576ba7abe47606fdf7667ea88e4bf8cf6dbefad to 9201066f0ce6d2aab310a875b00a2eeeb4ec75e3

Branch pushed to git repo; I updated commit sha1. New commits:

77651dfRemove unused keyword arguments
9201066Add density option

comment:24 follow-up: Changed 2 years ago by paulmasson

Travis, I went through the existing options, removed those that are not accepted by streamplot and added one that is, density. This value must be a single float or a pair of floats in a list, although I haven't added any validation for input values. It will at least allow you to see more streamlines and verify that they do not always start at integers.

There is another streamplot option I though might be interesting called start_points, but all it does is restrict streamlines to only those passing through an input list of points. You can test it by adding an option start_points=[[0,.1],[0,-.1]] or the like. Doesn't look too interesting to me.

comment:25 in reply to: ↑ 24 Changed 2 years ago by tscrim

  • Branch changed from u/paulmasson/10775 to u/tscrim/streamline_plot-10775
  • Commit changed from 9201066f0ce6d2aab310a875b00a2eeeb4ec75e3 to 584486bf675c0c4b813ca2a75f745369f55f65a6
  • Reviewers set to Travis Scrimshaw

Replying to paulmasson:

Travis, I went through the existing options, removed those that are not accepted by streamplot and added one that is, density. This value must be a single float or a pair of floats in a list, although I haven't added any validation for input values. It will at least allow you to see more streamlines and verify that they do not always start at integers.

I added some input checking, but density will still need to be documented.

There is another streamplot option I though might be interesting called start_points, but all it does is restrict streamlines to only those passing through an input list of points. You can test it by adding an option start_points=[[0,.1],[0,-.1]] or the like. Doesn't look too interesting to me.

It is useful for showing numerical solutions to initial value problems (IVPs) of (nonlinear) ODEs. Also useful for seeing how stable solutions are. Although it is very annoying that it can't plot start points that are too close together:

sage: x,y = var('x,y')
sage: f = (x + y) / sqrt(x^2 + y^2)
sage: pts = [[-1., -2.+val/8] for val in range(10)]
sage: streamline_plot((cos(x), sin(x+y)), (x,-3,3), (y,-3,3), start_points=pts)
/home/travis/sage-build/local/lib/python2.7/site-packages/sage/repl/rich_output/display_manager.py:590: RichReprWarning: Exception in _rich_repr_ while displaying object: 
  RichReprWarning,
Graphics object consisting of 1 graphics primitive

Although this is not something that would (nor should) block this getting a positive review.


New commits:

c443cc4Merge branch 'u/paulmasson/10775' of git://trac.sagemath.org/sage into u/paulmasson/10775
584486bAdding start points option and slope field input.

comment:26 Changed 2 years ago by paulmasson

  • Branch changed from u/tscrim/streamline_plot-10775 to u/paulmasson/streamline_plot-10775

comment:27 Changed 2 years ago by paulmasson

  • Commit changed from 584486bf675c0c4b813ca2a75f745369f55f65a6 to 9f3a5c25cd924b53bb8c4cb58fe6bd992440e40a

Ah, so that's what you had in mind for the additions.

The documentation would not build without some fixes. Sphinx plots do not allow ^ for exponentiation, needs to be **, and variable declarations do not carry over between them. Documentation now builds.


New commits:

9f3a5c2Fix Sphinx plots

comment:28 Changed 2 years ago by paulmasson

Don't know if you noticed, Travis, but in the example with start_points the streamlines come close to the start points but do not pass through them. I thought it might be because you didn't process the start points as in this example

http://matplotlib.org/examples/images_contours_and_fields/streamplot_demo_start_points.html

but even with setting up a numpy array as indicated and transposing it, the streamlines still do not pass through the start points. In fact the resulting graph is exactly the same as not bothering to process the points, so apparently that isn't necessary.

Problem with streamplot? Thought you should know before any students ask what's going on.

comment:29 Changed 2 years ago by tscrim

Thanks for fixing that Paul. I guess we should add a note and example about the behavior of start_points. Also, do you want me to document the density option?

Actually, I think I should not be normalizing the input when we want to consider it as a slope field, which makes sense when considering plotting the tangent lines. I believe it just adds a unnecessary computation that cancels when considering splitting dy / dx to (dx, dy).

comment:30 Changed 2 years ago by git

  • Commit changed from 9f3a5c25cd924b53bb8c4cb58fe6bd992440e40a to 478b1e19adc03befdc57d983aea8a09fc6e3e4dd

Branch pushed to git repo; I updated commit sha1. New commits:

78c53a0Improve documentation and example
478b1e1Code block for start_points

comment:31 Changed 2 years ago by paulmasson

Documentation for implemented plot options now included. I also updated your language for start_points using that from the matplotlib website and added the note.

I'm including a commented code block for the official way to handle start_points in the event it is useful later.

comment:32 follow-up: Changed 2 years ago by tscrim

I disagree with adding start_points to the @options, where if that is changed, it is perpetuated in all future plots, which is not what we want. I will revert this back when I make my next round of changes.

Does the "official way" for start_points make any difference?

comment:33 in reply to: ↑ 32 Changed 2 years ago by paulmasson

Replying to tscrim:

Does the "official way" for start_points make any difference?

Not that I can tell. xpos_array and ypos_array need to be explicit numpy arrays to avoid errors. The documentation says that start_points is also supposed to be a numpy array, so I added this block of commented code in case that is ever enforced, but you can remove it if you feel that is more appropriate.

comment:34 Changed 2 years ago by tscrim

Let's use the official way. Also (not that it makes any real difference), you don't need to pop off the option:

    if 'start_points' in options:
        xstart_array, ystart_array = [], []
        for point in options['start_points']:
            xstart_array.append(point[0])
            ystart_array.append(point[1])
        options['start_points'] = numpy.array([xstart_array, ystart_array]).T

comment:35 Changed 2 years ago by git

  • Commit changed from 478b1e19adc03befdc57d983aea8a09fc6e3e4dd to 79e879535dc2096ed3e15a77f4ff276253e72e44

Branch pushed to git repo; I updated commit sha1. New commits:

79e8795Update start_points handling

comment:36 Changed 2 years ago by tscrim

  • Status changed from needs_review to positive_review

Thank you for all your work on this.

comment:37 Changed 2 years ago by fbissey

Not sure what Volker will find in his current run but I am failing to build the documentation in sage-on-gentoo

[plotting ] /scratch2/portage/sci-mathematics/sage-9999/work/sage-9999/src-python2_7/build/lib/sage/plot/streamline_plot.py:docstring of sage.plot.streamline_plot.streamline_plot:35: WARNING: Exception occurred in plotting streamline_plot-1
[plotting ] from /scratch2/portage/sci-mathematics/sage-9999/work/sage-9999/src-python2_7/doc/en/reference/plotting/sage/plot/streamline_plot.rst:
[plotting ] Traceback (most recent call last):
[plotting ] File "/usr/lib64/python2.7/site-packages/matplotlib/sphinxext/plot_directive.py", line 517, in run_code
[plotting ] six.exec_(code, ns)
[plotting ] File "/usr/lib64/python2.7/site-packages/matplotlib/externals/six.py", line 672, in exec_
[plotting ] exec("""exec _code_ in _globs_, _locs_""")
[plotting ] File "<string>", line 1, in <module>
[plotting ] File "<string>", line 2, in <module>
[plotting ] File "/scratch2/portage/sci-mathematics/sage-9999/work/sage-9999/src-python2_7/build/lib/sage/misc/decorators.py", line 554, in wrapper
[plotting ] return func(*args, **options)
[plotting ] File "/scratch2/portage/sci-mathematics/sage-9999/work/sage-9999/src-python2_7/build/lib/sage/plot/streamline_plot.py", line 316, in streamline_plot
[plotting ] for point in options['start_points']:
[plotting ] TypeError: 'NoneType' object is not iterable

I have this several times. I note that the bots have a similar problem in their logs.

comment:38 Changed 2 years ago by tscrim

  • Branch changed from u/paulmasson/streamline_plot-10775 to u/tscrim/streamline_plot-10775
  • Commit changed from 79e879535dc2096ed3e15a77f4ff276253e72e44 to 640cb6aad1a6b6380adabcb2c4d24d10eca3e718
  • Status changed from positive_review to needs_review

I had thought the first part of comment:32 was done by Paul. I removed start_points from the @options and did some little doc tweaking. Paul, can you make a quick check of my changes?


New commits:

1e2b2d0Merge branch 'u/paulmasson/10775' of git://trac.sagemath.org/sage into u/tscrim/streamline_plot-10775
fa8cda9Merge branch 'u/paulmasson/streamline_plot-10775' of git://trac.sagemath.org/sage into u/tscrim/streamline_plot-10775
640cb6aSome last little touchups.

comment:39 Changed 2 years ago by paulmasson

  • Status changed from needs_review to positive_review

Documentation now builds. Setting back to positive review.

comment:40 Changed 2 years ago by vbraun

  • Branch changed from u/tscrim/streamline_plot-10775 to 640cb6aad1a6b6380adabcb2c4d24d10eca3e718
  • Resolution set to fixed
  • Status changed from positive_review to closed
Note: See TracTickets for help on using tickets.