Skip to content

Instantly share code, notes, and snippets.

@orzelc
Created January 18, 2014 18:51
Show Gist options
  • Save orzelc/8494567 to your computer and use it in GitHub Desktop.
Save orzelc/8494567 to your computer and use it in GitHub Desktop.
VPython mock-up of charged tape problem as a mass-on-a-spring, because the equations are easier to deal with.
from visual import *
# Constants and other definitions
e=1.6e-19
k=9e9
g=9.80
length=0.21
height=0.3
sep=0.1524
mass=0.0003
charge=2.0e-8
scale=1
spring=100
drag=0.001
# Pivot points for the pendulums
pivot1=sphere(pos=vector(sep/2,height,0), radius=0.02, color=color.white)
pivot2=sphere(pos=vector(-sep/2,height,0), radius=0.02, color=color.white)
# The suspended point charges that stand in for charged tapes
mass1=sphere(pos=vector(sep/2, height-length,0), radius=0.02, color=color.green)
mass1.m=mass
mass1.q=charge
mass1.p=vector(0,0,0)
mass1.pos.y=mass1.pos.y-mass1.m*g/spring #stretch to account for the mass
mass2=sphere(pos=vector(-sep/2,height-length,0), radius=0.02, color=color.yellow)
mass2.m=mass
mass2.q=charge
mass2.p=vector(0,0,0)
mass2.pos.y=mass2.pos.y-mass2.m*g/spring #stretch to account for the mass
# Strings for the pendulum, to make it look pretty
string1=cylinder(pos=pivot1.pos, axis=mass1.pos-pivot1.pos, radius=0.005, color=color.white)
string2=cylinder(pos=pivot2.pos, axis=mass2.pos-pivot2.pos, radius=0.005, color=color.white)
# Gravitational forces
fgrav1=mass1.m*g*vector(0,-1,0)
fgrav2=mass2.m*g*vector(0,-1,0)
# Treating the strings as springs to get the right tension force
stretch1=string1.axis.mag-length
stretch2=string2.axis.mag-length
ft1=-(spring*stretch1)*(string1.axis/string1.axis.mag)
ft2=-(spring*stretch2)*(string2.axis/string2.axis.mag)
# Electrostatic force between the charged balls
r=mass1.pos-mass2.pos
rhat=r/r.mag
felec1=-(k*mass1.q*mass2.q/r.mag2)*rhat
felec2=-felec1
# Arrows showing the force
fgarr1=arrow(pos=mass1.pos, axis=scale*fgrav1, color=color.blue)
fgarr2=arrow(pos=mass2.pos, axis=scale*fgrav2, color=color.blue)
ftarr1=arrow(pos=mass1.pos, axis=scale*ft1, color=color.magenta)
ftarr2=arrow(pos=mass2.pos, axis=scale*ft2, color=color.magenta)
felarr1=arrow(pos=mass1.pos, axis=scale*felec1, color=color.cyan)
felarr2=arrow(pos=mass2.pos, axis=scale*felec2, color=color.cyan)
# Loop parameters
tmax=5
t=0
dt=0.001
while t<tmax:
rate(1000)
#print(t)
#Recalculate forces
stretch1=string1.axis.mag-length
stretch2=string2.axis.mag-length
ft1=-(spring*stretch1)*(string1.axis/string1.axis.mag)
ft2=-(spring*stretch2)*(string2.axis/string2.axis.mag)
r=mass1.pos-mass2.pos
rhat=r/r.mag
felec1=-(k*mass1.q*mass2.q/r.mag2)*rhat
felec2=-felec1
#Drag force to make the system settle near equlibirum point
fdrag1=-drag*mass1.p/mass1.m
fdrag2=-drag*mass2.p/mass2.m
fnet1=fgrav1+ft1+felec1+fdrag1
fnet2=fgrav2+ft2+felec2+fdrag2
# Update momentum and position of the charged balls
mass1.p=mass1.p+fnet1*dt
mass1.pos=mass1.pos+mass1.p*dt/mass1.m
mass2.p=mass2.p+fnet2*dt
mass2.pos=mass2.pos+mass2.p*dt/mass2.m
# Update all the miscellaneous graphics
string1.axis=mass1.pos-pivot1.pos
string2.axis=mass2.pos-pivot2.pos
fgarr1.pos=mass1.pos
fgarr2.pos=mass2.pos
ftarr1.pos=mass1.pos
ftarr1.axis=scale*ft1
ftarr2.pos=mass2.pos
ftarr2.axis=scale*ft2
felarr1.pos=mass1.pos
felarr1.axis=scale*felec1
felarr2.pos=mass2.pos
felarr2.axis=scale*felec2
t=t+dt
# If the balls are in contact, stop the loop.
r=mass1.pos-mass2.pos
if r.mag<(mass1.radius+mass2.radius):
print(t)
t=tmax+dt
# Print out some parameters for data recording
angle=abs(asin((mass2.pos.x-pivot2.pos.x)/string2.axis.mag))
print(t/tmax,r.mag/(mass1.radius+mass1.radius),(r.x)/(sep),angle)
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment