# relativisticDynamicsVersion3.py

from pylab import *

from time import time

start = time()

c = 2.998e8 # Speed of light in m/s

m = 0.938e9 # Mass in eV/c^2

Efield = 5e6 # Electric field in Volts per meter

x = 0 # Position in meters

v = 0 # Velocity in meters/second

t = 0 # Time in seconds

dt = 1e-8 # Time STEP in seconds

lorentz = 1 / sqrt(1 - v**2 / c**2) ## ADDED ##

E = lorentz * m ## MODIFIED ##

# Create arrays using initial values

tArray = array(t)

xArray = array(x)

vArray = array(v)

EArray = array(E)

lorentzArray = array(lorentz) ## ADDED ##

while t < 1e-5:

# The dynamics:

lorentz = 1 / sqrt(1 - v**2 / c**2) ## ADDED ##

a = Efield*c**2/(lorentz**3*m) ## MODIFIED ##

t = t + dt

x = x + v*dt

v = v + a*dt

E = lorentz * m ## MODIFIED ##

# Append the new values onto arrays

tArray = append(tArray, t)

xArray = append(xArray, x)

vArray = append(vArray, v)

EArray = append(EArray, E)

lorentzArray = append(lorentzArray, lorentz) ## ADDED ##

end = time()

print('Elapsed time:', end-start)

# The following lines are added for Exercise 7

KE_numerical = lorentz * m - m

KE_analytical = Efield*x

print('Numerical KE: ', KE_numerical)

print('Analytical KE:', KE_analytical)

print('Percent Difference:', 100*abs(KE_numerical-KE_analytical)/KE_analytical)

# Create plots

figure(1)

plot(tArray, vArray, linewidth=2)

xlabel('Time (sec)')

ylabel('Speed (m/s)')

grid(True)

ylim([0, 3.1e8])

savefig('ultrarelativistic_v_vs_t.png')

show()

figure(2)

plot(tArray, xArray, linewidth=2)

xlabel('Time (sec)')

ylabel('Position (m)')

grid(True)

savefig('ultrarelativistic_x_vs_t.png')

show()

figure(3)

plot(tArray, EArray, linewidth=2)

xlabel('Time (sec)')

ylabel('Energy (eV)')

grid(True)

savefig('ultrarelativistic_E_vs_t.png')

show()