BUG: Fix streamplot when velocity component is exactly zero. #1189

Merged
merged 4 commits into from Sep 2, 2012

Conversation

Projects
None yet
4 participants
Contributor

tonysyu commented Sep 1, 2012

When one component of the velocity field is zero, streamplot fails to draw streamlines. Here's a short example of the issue:

import matplotlib.pyplot as plt
import numpy as np

N = 10
Y, X = np.mgrid[:N, :N]
U = np.ones((N, N))
V = np.zeros((N, N))

plt.streamplot(X, Y, U, V)
plt.show()

Explanation of fix:
Negative and positive steps are treated differently by the code. Before, a zero gradient (cx, cy) was getting treated as a negative step. This resulted in a step size of -inf; treating zero as a positive step results in a step size of (positive) inf. This inf value gets filtered out since the integration routine chooses the smaller of the two (dsx, dsy) step sizes.

This pull request fails (merged 5639974 into 9448ff2).

Owner

efiring commented Sep 1, 2012

@tonysyu, I'm sorry, but I really don't like either the original form or the fix. One should not divide by zero. Relying on it to return inf is not good. At the very least, it raises a RuntimeWarning. Depending on the types of the arguments, it could fail completely. The zero-denominator case should be trapped and handled appropriatedly.

Contributor

tonysyu commented Sep 1, 2012

@efiring Yeah, it felt a little flakey to me too, but I was just trying to make a minimal change.

Do you think checking for zero and setting the step size to inf would be appropriate? This fix would still rely on min removing the inf later on; if both steps are inf then it would fail. This shouldn't happen, though, since _euler_step is only called when the local velocity causes a trajectory to go beyond the domain boundaries (which wouldn't happen if both velocity components are zero).

Owner

efiring commented Sep 2, 2012

@tonysyu, I think your suggestion above is OK. I don't see anything better, offhand.

Contributor

tonysyu commented Sep 2, 2012

I've added the checks to prevent divide-by-zero warnings. I could explicitly check if both values are zero, but like I said, this shouldn't happen.

I also made a couple of changes to suppress other warnings. I usually don't like changes that venture beyond the PR's topic, but it seemed minor enough. I can revert these, if necessary.

This pull request fails (merged 1d4682c into 9448ff2).

Member

dmcdougall commented Sep 2, 2012

@tonysyu Perhaps I'm missing something but I have a potentially silly question. If the velocity field is close to zero, passive tracers (basically streamlines) don't move in the fluid. Is there any reason you're setting the step-size to inf as opposed to just not doing an Euler step?

Contributor

tonysyu commented Sep 2, 2012

Not a silly question at all: the idea here is that only one of the two (x, y) values should be set to inf. E.g. in a horizontal flow, it'd take an infinite step to reach the upper or lower boundary, but the distance to the left or right boundary would be finite.

This routine shouldn't get called if both velocity components are zero, which is why I'm not inclined to add an explicit check for that case. As far as being "close to zero", I'm not sure; were you thinking of some sort of threshold?

Member

dmcdougall commented Sep 2, 2012

Well, I was more thinking that, in your example of a horizontal flow, you would simply ignore the vertical component.

Contributor

tonysyu commented Sep 2, 2012

That's essentially what happens. The inf stuff in this PR is used to calculate the step size along the stream line. Then, this step size gets scaled to x and y components (for horizontal flow, the y-scale is exactly zero). Is this equivalent to what you had in mind?

@efiring efiring added a commit that referenced this pull request Sep 2, 2012

@efiring efiring Merge pull request #1189 from tonysyu/streamplot-fix-zero-velocity
BUG: Fix streamplot when velocity component is exactly zero.
4c1d072

@efiring efiring merged commit 4c1d072 into matplotlib:master Sep 2, 2012

1 check failed

default The Travis build failed
Details
Member

dmcdougall commented Sep 3, 2012

By 'stepsize' I meant the Euler step in the Runge Kutta scheme for the odes:

\dot{x} = u(x)
\dot{y} = v(x)

I think we're talking about different things.

Contributor

tonysyu commented Sep 4, 2012

There's definitely a lot of room for misunderstanding since there are 4 step sizes: dx, dy, ds, and dt. By your last comment, I assume you mean dt. Basically, the scheme calculates the ds (or similarly dt, since this takes a single step based on the velocities at x, y) required to reach the closest horizontal and vertical boundaries and chooses the shorter of the two.

When you say "we're talking about different things", do you mean that your original suggestion no longer applies? Or am I missing something?

Member

dmcdougall commented Sep 4, 2012

My suggestion was to enforce that ds (or dy) is zero when the velocity at x (or y) is zero. Then you don't need to deal with inf.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment