July 29th, 2013, 11:23 PM

How to Create a Multiple Point Vortice Vector Field
I want to create a vector field with multiple point vortices in a checkerboard arrangement with alternating rotational directions. I have written a code to create a single point vortex in the middle of my grid (see below) but would like to know how to go about creating this multiple point vector field.
Code:
from math import *
from pylab import *
# Set limits and number of points in grid
xmax = 10
xmin = xmax
NX = 30
ymax = 10
ymin = ymax
NY = 30
# Make grid and calculate vector components
x1 = linspace(xmin, xmax, NX)
y1 = linspace(ymin, ymax, NY)
X1, Y1 = meshgrid(x1, y1)
S21 = (X1)**2  (Y1)**2 # This is the radius squared
Bx1 = Y1/S21
By1 = +X1/S21
figure()
QP = quiver(X1,Y1,Bx1,By1)
quiverkey(QP, 0.85, 1.05, 1.0, '1 mT', labelpos='N')
# Set the left, right, bottom, top limits of axes
dx = (xmax  xmin)/(NX1)
dy = (ymax  ymin)/(NY1)
axis([xmindx, xmax+dx, ymindy, ymax+dy])
show()
Thanks in advance!
July 30th, 2013, 10:27 PM

Assuming superposition of linear systems, and assuming you just want to make a pretty picture, you could create the effect by shifting.
Code:
QP = quiver(X1,Y1,
Bx1(concatenate((Bx1[:,10:],Bx1[:,:10]),axis=1)+concatenate((Bx1[10:],Bx1[:10]))),
By1(concatenate((By1[:,10:],By1[:,:10]),axis=1)+concatenate((By1[10:],By1[:10]))))
For carefully positioned values, shift the axes, recompute the radius, and add all of them.
[code]
Code tags[/code] are essential for python code and Makefiles!
August 5th, 2013, 09:11 PM

Originally Posted by b49P23TIvg
Assuming superposition of linear systems, and assuming you just want to make a pretty picture, you could create the effect by shifting.
Code:
QP = quiver(X1,Y1,
Bx1(concatenate((Bx1[:,10:],Bx1[:,:10]),axis=1)+concatenate((Bx1[10:],Bx1[:10]))),
By1(concatenate((By1[:,10:],By1[:,:10]),axis=1)+concatenate((By1[10:],By1[:10]))))
For carefully positioned values, shift the axes, recompute the radius, and add all of them.
Thanks that works great! How can I keep adding vortices to continue the pattern? I've looked at the concatenate syntax and I don't really understand it.
August 13th, 2013, 09:48 PM

Are you still interested? I next solution will be a function accepting a list of centers for the vortices, and will work by recomputation not translation followed by summing. I can make the translation/sum work but I'd rather not. Even if it is cheaper.
[code]
Code tags[/code] are essential for python code and Makefiles!
August 15th, 2013, 07:13 AM

Originally Posted by b49P23TIvg
Are you still interested? I next solution will be a function accepting a list of centers for the vortices, and will work by recomputation not translation followed by summing. I can make the translation/sum work but I'd rather not. Even if it is cheaper.
Yeah mate i'm definitely still interested! as well as this, is there a simple way to output the magnitudes of all the points on the grid? say for the 10 by 10 grid, outputting the magnitude of each of the 100 points?
August 16th, 2013, 12:35 AM

Here we work in complex planeoften much easier for 2D problems (and I've now learned, was originally invented to handle planar geometry). The centers are given by a vector of centers x + y*1j. This plot nearly matches the previous response, except that shifting gives the effect of the solution in a repeated (coupled) grid of vortex centers. In this there are only the 3 centers.
Code:
from math import *
from pylab import *
TAU = 2*pi
interesting_but_wrong = False
# Set limits and number of points in grid
xmax = 10
xmin = xmax
NX = 30
ymax = 10
ymin = ymax
NY = 30
# Make grid and calculate vector components
x1 = linspace(xmin, xmax, NX)
y1 = linspace(ymin, ymax, NY)
X1, Y1 = meshgrid(x1, y1)
grid = X1 + 1j * Y1
j = 1j
t = (xmaxxmin)/6.0
T = t*j
centers = (7,0,7j)
U = j*zeros(grid.shape)
for (i, center) in enumerate(centers):
V = gridcenter
U += (V*(j*(12*(i&1))))/(abs(V)**2) # (j*(12*(i&1))) multiplies by +/ 1j matching the order of centers.
Bx1 = U.real
By1 = U.imag
figure()
QP = quiver(X1,Y1,Bx1,By1)
quiverkey(QP, 0.85, 1.05, 1.0, '1 mT', labelpos='N')
# Set the left, right, bottom, top limits of axes
dx = (xmax  xmin)/(NX1)
dy = (ymax  ymin)/(NY1)
axis([xmindx, xmax+dx, ymindy, ymax+dy])
print(abs(U)) # answering "is there a simple way to output the magnitudes of all the points on the grid? say for the 10 by 10 grid, outputting the magnitude of each of the 100 points?"
show()
Last edited by b49P23TIvg; August 16th, 2013 at 12:38 AM.
[code]
Code tags[/code] are essential for python code and Makefiles!