1/* -*-ePiX-*- */
2#include "epix.h"
3using namespace ePiX;
4
5// number of field lines per charge
6const int N(13);
7
8// charge magnitudes and locations
9const int charge1(1);
10const int charge2(1);
11
12const P Q1( 1,0);
13const P Q2(-1,0);
14
15// inverse-square electric field from a point charge at the origin
16P unit_charge(const P& arg)
17{
18  return recip(arg|arg)*arg;
19}
20
21// electric field simulated by superimposing individual fields, then
22// re-scaling so that field -> 0 near the charge location (improves
23// plot quality:)
24P E(double x, double y)
25{
26  P temp(x,y);
27  // superposition of charge fields
28  P E_temp(charge1*unit_charge(temp-Q1) + charge2*unit_charge(temp-Q2));
29
30  // re-scale
31  return (1.0/(E_temp|E_temp))*E_temp;
32}
33
34P potential(double x, double y)
35{
36  // J rotates a vector by 1/4 turn; parallel to equipotentials
37  return J(E(x,y));
38}
39
40const double MAX(3);
41
42int main()
43{
44  picture(P(-MAX,-MAX), P(MAX,MAX), "4x4in");
45
46  begin();
47  set_crop();
48  degrees(); // change angle mode
49
50  // plot field lines
51  blue();
52  for (int i=0; i < N; ++i)
53    {
54      // initial points trace a small circle about Q1 or Q2,
55      ode_plot(E, Q1+polar(0.05, i*360.0/N), 10, 120);
56      ode_plot(E, Q2-polar(0.05, i*360.0/N), 10, 120);
57
58      // location of arrowhead
59      P pt(flow(E, Q2-polar(0.05, i*360.0/N), 3, 12));
60    }
61
62  green();
63  for (int i=-10; i < 10; ++i)
64    {
65      ode_plot(potential, Q1+polar(0.25*pow(0.8, i), 0), 2*M_PI, 120);
66      ode_plot(potential, Q2-polar(0.25*pow(0.8, i), 0), 2*M_PI, 120);
67    }
68
69  dot_size(6);
70  circ(Q1);
71  circ(Q2);
72
73  magenta();
74  label(Q1, "$+$");
75  label(Q2, "$+$");
76
77  blue();
78  for (int i=0; i < N; ++i)
79    {
80      P pt(flow(E, Q2-polar(0.05, i*360.0/N), 3.5, 12));
81      arrow(pt, pt+0.01*E(pt.x1(), pt.x2()));
82      arrow(-pt, -pt+0.01*E(-pt.x1(), -pt.x2()));
83    }
84
85  end();
86}
87