I am trying to replicate some standard structures in the study of ellipse fields in Mathematica, and I'm having much more trouble than is really reasonable, so I'd like to put this on the floor here in the hope that someone here can point out some simpler method that I'm currently not seeing.
More specifically, I would like to have simple and reliable ways to produce streamline plots like the following,
from the paper
Geometry of phase and polarization singularities, illustrated by edge diffraction and the tides. M.V. Berry. M.S. Soskin, M.V. Vasnetsov (Eds.), Second International Conference on Singular Optics, Proc. SPIE, 4403, p. 1 (2001). Author eprint.
Basically, the deal is the following: you have some complex-valued vector field, say, $$ \mathbf F(x,y) = (1+iby, i(1-x)) $$ (where $b$ is a parameter, the main interesting points of which are $b\in\{-1,1,3\}$), and you see it as encoding an ellipse at each point much in the same way as for an electromagnetic wave (i.e. the ellipse is traced out by $\mathrm{Re}(\mathbf F(x,y)e^{-it})$ for $t\in[0,2\pi]$); you then look for the streamlines of the major and minor axis of these ellipses.
The major and minor axis of each ellipse can be obtained as in this question as the real and imaginary parts of $$ \mathbf A+i\mathbf B = \frac{\sqrt{\mathbf E^*\cdot \mathbf E^*}}{\left|\sqrt{\mathbf E^*\cdot \mathbf E^*}\right|}\mathbf E, $$ and as a naive start, one can do fairly well by putting these in explicitly and asking for a StreamPlot:
Block[{F, A, B, R = 0.6},
F[x_, y_] = {1 + I b y, I (1 - x)};
A[x_, y_] =
1/Abs[Sqrt[F[x, y]\[Conjugate].F[x, y]\[Conjugate]]] Re[
Sqrt[F[x, y]\[Conjugate].F[x, y]\[Conjugate]] F[x, y]];
B[x_, y_] =
1/Abs[Sqrt[F[x, y]\[Conjugate].F[x, y]\[Conjugate]]] Im[
Sqrt[F[x, y]\[Conjugate].F[x, y]\[Conjugate]] F[x, y]];
Table[
StreamPlot[
A[x, y]
, {x, -R, R}, {y, -R, R}
, ImageSize -> 400
]
, {b, {-1, 1, 5}}]
]
However, I have a couple of problems here:
I cannot control where the streamlines start and stop, and the defaults do some pretty ugly spacing and they start and eliminate streamlines in mid-flight, which I would like to avoid. For some reason, if I try to give it an explicit
StreamPointsspecification, Mathematica just hangs forever.The stream plots have a discontinuity on the negative $x$ axis, where the major axis $\mathbf A(x,y)$ discontinuously changes sign. This discontinuity is both:
- spurious, because ellipses are symmetric and it doesn't matter which way you choose the axis, but also
- unavoidable, because of the way these ellipse points have been set up: if you track the ellipse around the origin, it will come back to itself while acquiring a 180° rotation; this won't affect the ellipse but it means that $\mathbf A(x,y)$ must change direction at some point.
I would like to have a solid construction that allows me to get past both of these limitations, but after several attacks I can't really get there. I'm OK with NDSolveing my way through this, but it keeps choking on the discontinuity and with the accumulation point at the origin.
(Also, I am wary of solutions that place explicit stock in the symmetry of the problem, by e.g. splitting up the plane into different regions. It's not a deal killer, but I would like to be able to handle cases where the branch cut falls along some non-symmetrically-placed line.)
Does anyone see how this can be done?
