# paste this code at the end of VectorFieldPlot 2.4
# https://commons.wikimedia.org/wiki/User:Geek3/VectorFieldPlot
doc = FieldplotDocument('VFPt_flat_magnets_gap_potential+contour', commons=True,
width=800, height=800)
Bfield = Field({'coils':[ [0, 1, pi/2, 2, 0.5 ,1],
[0, -1, pi/2, 2, 0.5 ,1] ]})
Hfield = Field([ ['charged_disc', {'x0':-2, 'y0':-1.5, 'x1':2, 'y1':-1.5, 'Q':-1}],
['charged_disc', {'x0':-2, 'y0':-0.5, 'x1':2, 'y1':-0.5, 'Q':1}],
['charged_disc', {'x0':-2, 'y0':0.5, 'x1':2, 'y1':0.5, 'Q':-1}],
['charged_disc', {'x0':-2, 'y0':1.5, 'x1':2, 'y1':1.5, 'Q':1}] ])
doc.draw_magnets(Bfield)
U0 = Hfield.V([0., 1.5 + 0.02])
doc.draw_scalar_field(func=Hfield.V, cmap=doc.cmap_AqYlFs, vmin=-U0, vmax=U0)
U1 = Hfield.V([0., 1.5])
doc.draw_contours(func=Hfield.V, levels=sc.linspace(-U1, U1, 11)[1:-1])
nlines = 22
R0 = op.brentq(lambda x: Bfield.F([x, 0.])[1], 0, 3)
Sp = Startpath(Bfield, lambda t: sc.array([-R0 + 2. * R0 * t, 0.]))
xstart = [Sp.startpos((0.2+i) / (nlines-0.6))[0] for i in range(nlines)]
cond = lambda xy: fabs(xy[1]) < 1e-2 or fabs(xy[1]) > 1.4
for iline, x in enumerate(xstart):
line = FieldLine(Bfield, [x, 0.], directions='both', maxr=12)
doc.draw_line(line, linewidth=2.4, arrows_style={'potential':Hfield.V,
'at_potentials':[-0.3*U1, 0., 0.3*U1], 'condition_func':cond})
for x0, y0 in ((-1, -1), (-1, 1), (1, -1), (1, 1)):
line = FieldLine(Bfield, [2.3 * x0, 1. * y0], directions='both', maxr=5)
doc.draw_line(line, linewidth=2.4, arrows_style={'dist':2,
'offsets':{'start':1, 'end':0} })
doc.write()