3
$\begingroup$

When I plot

r[x_, y_] := Sqrt[x^2 + y^2];
θ[x_, y_] := ArcTan[x, y];
a=3;

Q[k_, N_] := Module[{pts1}, pts1 = Table[{0, x}, {x, -2 N, 2 N, 0.5}];
StreamPlot[{-Cos[θ[x, y] k] a^(1/2)/r[x, y]^(
1/2) + (1 - Sin[θ[x, y] k]^2 a/r[x, y])^(1/2) + 
a^(1/2)/r[x, y]^(1/2) Cos[θ[x, y] k], 
a^(1/2)/r[x, y]^(1/2) Sin[θ[x, y] k]}, {x, -2 N, 
2 N}, {y, -2 N, 2 N}, PlotRange -> {Full, Full}, 
RegionFunction -> Function[{x, y}, x^2 + y^2 > a^2], 
StreamPoints -> {pts1, Automatic, Scaled[2]}, 
StreamStyle -> {Blue, Thick, "Line"}, PerformanceGoal -> "Quality", 
StreamScale -> Full]

Mathematica doesn't give any output.

Whereas, if I put the Abs around the square root bit, plotting something similar to what is shown here it works, even if I am not sure that those discontinuities in the gradients should appear (I mean the spikes near the x axes).

I report the code for the second case, the one that works, and the figure:

r[x_, y_] := Sqrt[x^2 + y^2];
θ[x_, y_] := ArcTan[x, y];
a = 3;
Q[k_, N_] := 
 Module[{pts1}, pts1 = Table[{0, x}, {x, -2 N, 2 N, 0.5}];
 StreamPlot[{-Cos[θ[x, y] k] a^(1/2)/r[x, y]^(1/2) + 
 Abs[(1 - Sin[θ[x, y] k]^2 a/r[x, y])]^(1/2) + 
 a^(1/2)/r[x, y]^(1/2) Cos[θ[x, y] k], 
 a^(1/2)/r[x, y]^(1/2) Sin[θ[x, y] k]}, {x, -2 N, 
 2 N}, {y, -2 N, 2 N}, PlotRange -> {Full, Full}, 
 RegionFunction -> Function[{x, y}, x^2 + y^2 > a^2], 
 StreamPoints -> {pts1, Automatic, Scaled[2]}, 
 StreamStyle -> {Blue, Thick, "Line"}, PerformanceGoal -> "Quality",
 StreamScale -> Full]]

 Q[1/2, 10]

enter image description here

What is wrong with that square root?

$\endgroup$
2
  • $\begingroup$ RegionPlot[(1 - Sin[\[Theta][x, y]/2]^2 a/r[x, y]) < 0, {x, -10, 10}, {y, -10, 10}] $\endgroup$ Feb 4, 2015 at 19:57
  • $\begingroup$ @belisarius That's not the problem, it's just that I am not expecting to see those points on the left of the empty circle. $\endgroup$ Feb 4, 2015 at 21:17

1 Answer 1

2
$\begingroup$

In order to examine the offending Sqrt, I added

EvaluationMonitor :> If[1 - Sin[θ[x, y] k]^2 a/r[x, y] < 0, 
  Print[x, "  ", y, "  ", r[x, y], "  ", θ[x, y], "  ", 1 - Sin[θ[x, y] k]^2 a/r[x, y]]]

to the first function Q in the Question. Interestingly, numerous instances of r[x, y] < a were printed. For instance

-2.27757  -1.41668  2.68223  -2.58514  -0.0341049

Of course, the Sqrt would produce complex numbers in those instances, and I imagine that StreamPlot could not cope with them when it began plotting. Putting Abs around the argument of the Sqrt, as in the second definition of Q, eliminates this issue, of course. To corroborate this assessment, I replaced

(1 - Sin[θ[x, y] k]^2 a/r[x, y])^(1/2)

by

(Max[1 - Sin[θ[x, y] k]^2 a/r[x, y], 0])^(1/2)

in the first definition of Q, and it then worked fine.

It is not obvious to me why StreamPlot computes and uses values outside the RegionFunction. Perhaps, it needs them in order to construct accurate stream lines near the edges of the RegionFunction.

By the way, "those discontinuities in the gradients" are caused by the discontinuity in Sin[θ[x, y] k] for k = 1/2 in

a^(1/2)/r[x, y]^(1/2) Sin[θ[x, y] k]

as can be seen from

Plot3D[a^(1/2)/r[x, y]^(1/2) Sin[θ[x, y] /2], {x, -20, 20}, {y, -20, 20}, 
 RegionFunction -> Function[{x, y, z}, x^2 + y^2 > a^2], PlotRange -> All]

vy in 3D

Update: More Strange Behavior

As an experiment, I changed the range of StreamPlot from {x, -2 N, 2 N}, {y, -2 N, 2 N} to {x, -2 N, 2 N}, {y, 0, 2 N}. The first definition of Q then produced a plot, although not what one might have expected.

restricted range

Two aspects of this plot are striking. Most streamlines end at x = 0, and streamlines also are plotted for -2 < y < 0. (The arrow heads appear in all cases that I have run. I use V 10.0.2.0.) Using EvaluationMonitor as above showed that StreamPlot continued to evaluate the functions outside RegionFunction, and additionally for y < 0.

On the other hand, the change of range in the second definition of Q produced the expected complete set of streamlines for y > 0 but also produced streamlines for -2 < y < 0, just as in the plot above.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.