#To run this program install Python and Vpython
#see: http://vpython.org/
#for Windows download and install Python and Vpython from:
# http://vpython.org/contents/download_windows.html
#Copy this program and paste it in IDLE (for example), save it as
# Synchrotron.py and run it.
#
'''
Simple Synchrotron model
The arrow represent the accelerated Electric field
Curved and straight parts drawn with coils so that the electron can be seen
In curved parts the particle trajectory is curved to magnetic field.
'''
#Author: Manuel Jose Paez
#Instituto de Fisica
#Universidad de Antioquia
#Medellin, Colombia
#mpaez@fisica.udea.edu.co
from visual.graph import * #graphics and math classes
scene=display(x=0,y=0,width=600,height=600,range=200,
title='Synchrotron. B: cyan arrows, E: electric field white arrow')
f=frame()
electron=sphere(pos=(0,150,0),radius=5,color=color.green) #the electron
xf=zeros(5,Float)
yf=zeros(5,Float)
xc=zeros(5,Float)
yc=zeros(5,Float)
xvert=zeros((13),Float) #vertices x-component
yvert=zeros((13),Float) #vertices y-component
r=150.0
Efield=arrow()
B=zeros((9,2),Float)
def magnet(angle):
mag1=curve(frame=f,range=200,color=color.yellow,radius=1)
mag2=curve(frame=f,range=200,color=color.yellow,radius=1)
mag3=curve(frame=f,range=200,color=color.yellow,radius=1)
mag4=curve(frame=f,range=200,color=color.yellow,radius=1)
mag5=curve(frame=f,range=200,color=color.yellow,radius=1)
mag6=curve(frame=f,range=200,color=color.yellow,radius=1)
for i in range(0,5):
th=math.pi/24.0
xf[i]=170.0*math.sin(angle+i*th)
yf[i]=170.0*math.cos(angle+i*th)
xc[i]=130.0*math.sin(angle+i*th)
yc[i]=130.0*math.cos(angle+i*th)
mag1.pos=([(xf[0],yf[0],30),(xf[1],yf[1],30),(xf[2],yf[2],30),(xf[3],yf[3],30),(xf[4],yf[4],30),
(xc[4],yc[4],30),(xc[3],yc[3],30),(xc[2],yc[2],30),(xc[1],yc[1],30),(xc[0],yc[0],30),(xf[0],yf[0],30)])
mag2.pos=([(xf[0],yf[0],-30),(xf[1],yf[1],-30),(xf[2],yf[2],-30),(xf[3],yf[3],-30),(xf[4],yf[4],-30),
(xc[4],yc[4],-30),(xc[3],yc[3],-30),(xc[2],yc[2],-30),(xc[1],yc[1],-30),(xc[0],yc[0],-30),(xf[0],yf[0],-30)])
mag3.pos=([(xf[0],yf[0],-30),(xf[0],yf[0],30)])
mag4.pos=([(xf[4],yf[4],-30),(xf[4],yf[4],30)])
mag5.pos=([(xc[0],yc[0],-30),(xc[0],yc[0],30)])
mag6.pos=([(xc[4],yc[4],-30),(xc[4],yc[4],30)])
def pipestraight(angle): #straight parts of the accelerator
pio6=math.pi/6.
xpos=r*math.sin(angle)
ypos=r*math.cos(angle)
if angle<11.*math.pi/6.0:
pipe=helix(pos=(xpos,ypos,0),radius=15,axis=(r*(math.sin(pio6+angle)-math.sin(angle)),
+r*(math.cos(pio6+angle)-math.cos(angle)),0),coils=78*0.25, thickness=2,color=color.orange)
else:
pipe=helix(pos=(xpos,ypos,0),radius=15,axis=(r*(math.sin(pio6+angle)-math.sin(angle)),
+r*(math.cos(pio6+angle)-math.cos(angle)),0),coils=78*0.2, thickness=2,color=color.orange)
def pipecurved(angle): #Where the magnets are located
th=math.pi/36.0
for i in range (0,6):
xpos=150.*math.sin(angle+i*th)
ypos=150.*math.cos(angle+i*th)
sect=helix(pos=(xpos,ypos,0), radius=15,thickness=2,color=color.orange,coils=3,
axis=(150*(math.sin(angle+(i+1)*th)-math.sin(angle+i*th)),
150*(math.cos(angle+(i+1)*th)-math.cos(angle+i*th)),0))
def structure(): #ensemble all parts to form sinchrotron
for i in range (0,6):
magnet(i*math.pi/3.0)
pipecurved(i*math.pi/3.)
for i in range (1,7):
pipestraight((2*i-1)*math.pi/6.)
def setB(): #nine arrows to represent B field in each curved path
B[0,0]=150.0 #radius
B[0,1]=15.0 #angle
B[1,0]=150.0
B[1,1]=7.5
B[2,0]=150.0
B[2,1]=22.5
B[3,0]=140.0
B[3,1]=15.0
B[4,0]=140.0
B[4,1]=7.5
B[5,0]=140.0
B[5,1]=22.5
B[6,0]=160.0
B[6,1]=15.0
B[7,0]=160.0
B[7,1]=7.5
B[8,0]=160.0
B[8,1]=22.5
def Bfield(radi,angl):
radians=angl*math.pi/180.0
for i in range (0,7):
x=radi*math.sin(radians+i*math.pi/3.0)
y=radi*math.cos(radians+i*math.pi/3.0)
#d=frame() #to make cilinder with cone one unit
cylinder(pos=(x,y,30),radius=2,axis=(0,0,-50),color=color.cyan)
#cone(pos=(x,y,0),radius=4,axis=(0,0,-10),color=color.cyan)
cone(frame=f,pos=(x,y,-20),radius=4,axis=(0,0,-15),color=color.cyan)
#d.pos=vector(x,40,z+15)
def vertices():
for i in range (0,13):
xvert[i]= r*math.sin(i*math.pi/6.0)
yvert[i]=r*math.cos(i*math.pi/6.0)
def orbit(n):
angle=0.0
Eangle=0.0
th=math.pi/(6.0*n)
ar=0.7*cos(Eangle)
fac=1
for j in range(1,7):
for i in range (1,n+1): #curved parts
#print j, angle,Eangle*180/math.pi
rate(25+5*j) #to show it faster
electron.pos=vector(r*math.sin(angle+th*i),r*math.cos(angle+th*i))
Eangle +=30.0*math.pi/180./n
ar=0.7*cos(Eangle)
if ar<0:# position of the arrow depends if E positive or negative
fac=2.0
else:
fac=1.0
Efield.pos=(r*math.sin(2*math.pi-fac*math.pi/18.),
r*math.cos(2*math.pi-fac*math.pi/18.),0)
Efield.axis=(-ar*(xvert[12]-xvert[11]),-ar*(yvert[12]-yvert[11]),0)
xincr=(xvert[2*j]-xvert[2*j-1])/n #number of steps in straigh part
yincr=(yvert[2*j]-yvert[2*j-1])/n #in x and y, changes with j
m=n #number of steps for electron
if j==6:
m=n-2 #particle is accelerated in this straight section 6
for i in range(1,m+1): #straight parts
rate(25+5*j) #faster rate each whole turn
electron.pos=vector(xvert[2*j-1]+i*xincr,yvert[2*j-1]+i*yincr,0)
Eangle +=30.0*math.pi/180.0/m #simulates increase in E frequency
ar=0.7*cos(Eangle)
if ar<0:
fac=2.0
else:
fac=1.0
Efield.pos=(r*math.sin(2*math.pi-fac*math.pi/18.),
r*math.cos(2*math.pi-fac*math.pi/18.),0)
Efield.axis=(-ar*(xvert[12]-xvert[11]),-ar*(yvert[12]-yvert[11]),0)
angle +=math.pi/3.0
structure()
vertices()
setB()
kB=0
for i in range (20,3, -2):
Br=B[kB,0] #radius to place B
Ban=B[kB,1] #angle to put it
Bfield(Br,Ban)
orbit(i)
kB +=1