Charge Seesaw

Consider a moment at time $t$ where the rod makes an angle of $\theta$ with the horizontal. At this configuration, the coordinates of each of the masses are:

$\vec{r}_1 = \left(\begin{matrix} \cos{\theta} \\ \sin{\theta}\end{matrix}\right)$ $\vec{r}_2 = \left(\begin{matrix} -\cos{\theta} \\ -\sin{\theta}\end{matrix}\right)$

Velocities are:

$\vec{v}_1 = \dot{\theta}\left(\begin{matrix} -\sin{\theta} \\ \cos{\theta}\end{matrix}\right)$ $\vec{v}_2 = \dot{\theta}\left(\begin{matrix} \sin{\theta} \\ -\cos{\theta}\end{matrix}\right)$

The total electrostatic potential energy of this system of three charges in this configuration is:

$\mathcal{V} = \frac{1}{\sqrt{(1-\cos{\theta})^2 + (1+\sin{\theta})^2}} + \frac{1}{\sqrt{(1+\cos{\theta})^2 + (1-\sin{\theta})^2}} + \frac{1}{2}$

The total kinetic energy of this system is:

$\mathcal{T} = \frac{m}{2} \lvert \vec{v}_1 \rvert^2 + \frac{m}{2} \lvert \vec{v}_2 \rvert^2$ $\implies \mathcal{T} = \dot{\theta}^2$

It can be seen that: $\theta(0) = \dot{\theta}(0) = 0$ . Applying the conservation of energy equation leads to:

$\mathcal{T} + \mathcal{V} = \mathcal{T}_{initial} + \mathcal{V}_{initial}$

This gives:

$\dot{\theta}^2 = \frac{1}{2} + 1 + \frac{1}{\sqrt{5}} - \mathcal{V}$ $\dot{\theta}^2 = \frac{3}{2} + \frac{1}{\sqrt{5}} - \mathcal{V}$

As time progresses, $\theta$ will increase since the charge on the right is initially closer to the fixed charge and will experience a greater upwards repulsive force that sets the system into motion.

$\implies \dot{\theta} = \left(\frac{3}{2} + \frac{1}{\sqrt{5}} - \mathcal{V}\right)^{1/2}$

Now, with a little bit of visualization, one can imagine that as the rod rotates by 90 degrees, it's configuration will be similar to it's initial configuration in the sense that the potential energy of the system at that instant will be similar to the initial potential energy. If the rod were to rotate beyond 90 degrees, the potential energy must increase as the movable charge comes closer to the fixed charge. An increase of PE would violate the energy conservation equation, which is a physical impossibility. Thus, the rod rotates from zero to 90 degrees and back to complete one full oscillation. Moving from zero to 90 degrees will be done in half a time period. Therefore, this half time period can be computed by separating and integrating the above equation as follows:

$\int_{0}^{\pi/2} \frac{d\theta}{\left(\frac{3}{2} + \frac{1}{\sqrt{5}} - \mathcal{V}\right)^{1/2}} = \int_{0}^{T/2} dt$

The evaluation beyond this point is done numerically.

$\boxed{T \approx 8.689466511664962}$

@Karan Chatrath has shown how to use energy methods to solve this problem. This method relies on recognizing that the trajectory is determined by the location of the "equipotential" points. Another way is to model the torques and run a time-domain simulation. I have included the code for reference. The code looks for the time separation between two local maxima of the rod angle. The time-domain solution shows the rod angle oscillating between $0$ and $90$ degrees, just as the equipotential conditions suggest.

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105

import math

dt = 10.0**(-5.0)

m = 1.0
q = 1.0
R = 1.0
k = 1.0

x3 = 1.0
y3 = -1.0

I = 2.0*m*(R**2.0)

####################################

t = 0.0
count = 0

theta = 0.0
thetad = 0.0
thetadd = 0.0

max1 = 0.0
max2 = 0.0

while t <= 50.0:

    theta = theta + thetad*dt
    thetad = thetad + thetadd*dt

    x1 = R*math.cos(theta)
    y1 = R*math.sin(theta)

    x2 = R*math.cos(theta + math.pi)
    y2 = R*math.sin(theta + math.pi)

    D31x = x1 - x3
    D31y = y1 - y3

    D32x = x2 - x3
    D32y = y2 - y3

    D31 = math.sqrt(D31x**2.0 + D31y**2.0)
    D32 = math.sqrt(D32x**2.0 + D32y**2.0)

    u31x = D31x/D31
    u31y = D31y/D31

    u32x = D32x/D32
    u32y = D32y/D32

    F31 = k*q*q/(D31**2.0)
    F32 = k*q*q/(D32**2.0)

    F31x = F31*u31x
    F31y = F31*u31y

    F32x = F32*u32x
    F32y = F32*u32y

    T31 = x1*F31y - y1*F31x
    T32 = x2*F32y - y2*F32x
    T = T31 + T32

    thetadd = T/I

    if (t > 0.0) and (t < 10.0) and (theta > max1):
        max1 = theta
        t1 = t

    if (t > 10.0) and (t < 20.0) and (theta > max2):
        max2 = theta
        t2 = t

    t = t + dt
    count = count + 1

    theta_deg = 180.0*theta/math.pi

    #if count % 1000 == 0:
        #print t,theta_deg

####################################

print ""
print ""

T = t2 - t1

print dt
print T

#>>> 
#0.001
#8.683
#>>> ================================ RESTART ================================
#>>> 
#0.0001
#8.68889999998
#>>> ================================ RESTART ================================
#>>> 
#1e-05
#8.68940999967
#>>>