# Contour Plot not working with pseudo-piecewise function

Hi, I'm fairly new to Sage and don't know if this is a bug or just my own limited understanding. I want to plot a function 'f(g1, g2)' but only in some relevant area of the plot. All other points should display/evaluate to zero. This should work fine in Octave or numpy but Sage for some reason gives me an error-message which I fail to understand.

Function assignment like so:

def f(g1, g2):
if g1*g2 > 0.0:
if abs(g2) < 1/abs(g1):
return sqrt(0.5*((1+g1*g2)*abs(g1+g2))/(sqrt((1.0-g1*g2)^3)*g1*g2))
else:
return 0


Then the Plot:

var('g1 g2')
contour_plot(f(g1,g2), (g1,-4,4), (g2,-4,4),plot_points=500, cmap='afmhot_r', colorbar=True, contours=srange(1,6,0.5))


and at last the error-message:

ValueError: zero-size array to ufunc.reduce without identity


Strangely, if I assign a single number to 'contours' (or leave it out to get the default value) the error-message does not show but instead I get an empty/white plot.

If I do not use a conditional/if-statement in my function-assignment the contour-plot shows just fine (except that it didn't do what I wanted). Does anyone know how to fix this or why it doesn't work or at least why I'm being silly for even trying it that way?

edit retag close merge delete

Sort by » oldest newest most voted I can't reproduce the error (4.7.1), but I can see at least one problem.

The way the scoping works, f(g1, g2) is evaluated before the plot is made, so your plot will be zero (or maybe None) everywhere. Stick a "print g1,g2" inside the function to convince yourself this is true. You probably want "f" instead of "f(g1,g2)".

The following seems to work for me. Note that I've also changed it so that it returns 0 instead of None sometimes like your original code did-- if that's what you want, switch back.

var("g1 g2")
def f(g1, g2):
if g1*g2 > 0.0:
if abs(g2) < 1/abs(g1):
return sqrt(0.5*((1+g1*g2)*abs(g1+g2))/(sqrt((1.0-g1*g2)^3)*g1*g2))
return 0

z=contour_plot(f, (g1,-4,4), (g2,-4,4),plot_points=500, cmap='afmhot_r', colorbar=True, contours=srange(1,6,0.5))
show(z) more

Awesome, thanks for the quick reply. I had never considered the scope being an issue here.