Testing Gauss’ Law for a Relativistically Oscillating
Point Charge in an Egg
G.R.Dixon, July 22, 2010
In a previous article it was demonstrated that Gauss’ law is obeyed for a relativistically oscillating point charge within an ellipsoid. In this article it is demonstrated that the same may be true when the charge oscillates within an egg-shaped surface. The flux of E through this surface is computed at time t=0, which is when the particle passes through the origin. Once again the particle’s motion is
.
The egg-shaped cross section is formed by joining a circle (left "half") to an ellipse (right "half"). Fig. 1 plots the top half of the egg cross section.
Figure 1

Egg Cross Section
Fig. 2 plots the normal component of E at points on the egg.
Figure 2

Enormal vs x
The flux of Enormal through the ellipsoid surface is computed to be 1.130155E8. On the other hand, q/eo1.129943E8. The percent difference (1.9E-2%) might be ascribable to numerical error, or it might denote a small but significant failure of Gauss’ law for the egg-shaped surface. The program used to compute the data in this article is included in the Appendix.
***Appendix***
Private Sub cmdEgg_Click()
'*******************
'Compute the E Field Flux through an egg-shaped surface containing
'an oscillating charge.
'*******************
Const c As Double = 300000000# 'speed of light
Const eps0 As Double = 0.00000000000885 'permittivity constant
Const pi As Double = 3.141592654
Const q As Double = 0.001 'oscillating charge equals .001 coulomb
Const A As Double = 1 'amplitude of oscillation equals 1 meter
Const omega As Double = 0.99 * c / A 'max q speed is highly relativistic
Const steps As Long = 5000 'number of iterations
Dim deltax As Double
Dim i, j, jend As Long 'loop counter
Dim dArea(steps) As Double 'Area increment of ring
Dim Px, Py As Double
Dim Ex(steps), Ey(steps) As Double 'electric field components
Dim Enormal(steps) As Double
Dim dtmin, dtmax, dt As Double
Dim ux, uy As Double
Dim drx, dry, dr As Double
Dim tr As Double 'retarded time
Dim x(steps), xr As Double
Dim y(steps) As Double
Dim dydx(steps) As Double
'Dim r(steps) As Double
Dim ds(steps) As Double
Dim vr As Double
Dim ar As Double
Dim Eflux As Double
Dim Flux(steps) As Double 'Flux increments
Dim normx(steps) As Double
Dim normy(steps) As Double
Dim deltay(steps) As Double
Dim slopenormal(steps) As Double
Dim phi(steps) As Double
Eflux = 0
deltax = (1 + Sqr(2)) / steps
'MsgBox ("deltax = " & deltax)
j = 0
x(j) = -1
Do While x(j) < 0
j = j + 1
x(j) = -1 + j * deltax
y(j) = Sqr(1 - x(j) ^ 2)
Loop
jend = j
'MsgBox ("jend = " & jend)
For j = jend To steps - 1
x(j) = 0 + (j - jend) * deltax
y(j) = Sqr(1 - x(j) ^ 2 / 2)
Next j
Open "c:\\WINMCADC\Physics\GaussTest.PRN" For Output As #1
For j = 0 To steps - 1
Write #1, x(j), y(j)
Next j
Close
MsgBox ("Ready to plot x and y")
For j = 0 To jend - 1
If x(j) <> -1 Then
deltay(j) = -x(j) / Sqr(1 - x(j) ^ 2) * deltax
Else
deltay(j) = deltax
End If
ds(j) = Sqr(deltax ^ 2 + deltay(j) ^ 2)
dArea(j) = 2 * pi * y(j) * ds(j)
Next j
For j = jend To steps - 1
If x(j) <> Sqr(2) Then
deltay(j) = 1 / 2 / Sqr(1 - x(j) ^ 2 / 2) * x(j) * deltax
Else
deltay(j) = deltax
End If
ds(j) = Sqr(deltax ^ 2 + deltay(j) ^ 2)
dArea(j) = 2 * pi * y(j) * ds(j)
Next j
For j = 0 To jend - 1
dydx(j) = deltay(j) / deltax
If dydx(j) <> 0 Then
slopenormal(j) = Abs(-1 / dydx(j))
phi(j) = Atn(slopenormal(j))
normx(j) = -Cos(phi(j))
normy(j) = Sin(phi(j))
Else
normx(j) = 0
normy(j) = 1
End If
Px = x(j)
Py = y(j)
dtmin = 0
dtmax = 10 * A
Do
dt = (dtmin + dtmax) / 2
tr = -dt
xr = A * Sin(omega * tr)
drx = Px - xr
dry = Py
dr = Sqr(drx ^ 2 + dry ^ 2)
If Abs(c * dt - dr) < 2 ^ (-30) Then Exit Do
If c * dt - dr > 0 Then
dtmax = dt
Else
dtmin = dt
End If
Loop
vr = omega * A * Cos(omega * tr)
ar = -(omega ^ 2) * A * Sin(omega * tr)
ux = c * drx / dr - vr
uy = c * dry / dr
Ex(j) = q / (4 * pi * eps0) * dr / (drx * ux + dry * uy) ^ 3 * (ux * (c ^ 2 - vr ^ 2) + dry * (-uy * ar))
Ey(j) = q / (4 * pi * eps0) * dr / (drx * ux + dry * uy) ^ 3 * (uy * (c ^ 2 - vr ^ 2) - drx * (-uy * ar))
Enormal(j) = Ex(j) * normx(j) + Ey(j) * normy(j)
Flux(j) = Enormal(j) * dArea(j)
Eflux = Eflux + Flux(j)
Next j
For j = jend To steps - 1
If j = jend Then
dydx(j) = 0
Else
dydx(j) = deltay(j) / deltax
End If
If dydx(j) <> 0 Then
slopenormal(j) = Abs(-1 / dydx(j))
phi(j) = Atn(slopenormal(j))
normx(j) = Cos(phi(j))
normy(j) = Sin(phi(j))
Else
normx(j) = 0
normy(j) = 1
End If
Px = x(j)
Py = y(j)
dtmin = 0
dtmax = 10 * A
Do
dt = (dtmin + dtmax) / 2
tr = -dt
xr = A * Sin(omega * tr)
drx = Px - xr
dry = Py
dr = Sqr(drx ^ 2 + dry ^ 2)
If Abs(c * dt - dr) < 2 ^ (-30) Then Exit Do
If c * dt - dr > 0 Then
dtmax = dt
Else
dtmin = dt
End If
Loop
vr = omega * A * Cos(omega * tr)
ar = -(omega ^ 2) * A * Sin(omega * tr)
ux = c * drx / dr - vr
uy = c * dry / dr
Ex(j) = q / (4 * pi * eps0) * dr / (drx * ux + dry * uy) ^ 3 * (ux * (c ^ 2 - vr ^ 2) + dry * (-uy * ar))
Ey(j) = q / (4 * pi * eps0) * dr / (drx * ux + dry * uy) ^ 3 * (uy * (c ^ 2 - vr ^ 2) - drx * (-uy * ar))
Enormal(j) = Ex(j) * normx(j) + Ey(j) * normy(j)
Flux(j) = Enormal(j) * dArea(j)
Eflux = Eflux + Flux(j)
Next j
Open "c:\\WINMCADC\Physics\GaussTest.PRN" For Output As #1
For j = 0 To steps - 1
Write #1, x(j), Enormal(j)
Next j
Close
MsgBox ("Ready to plot Enormal")
MsgBox ("E Flux = " & Eflux & "; Gauss = " & q / eps0)
MsgBox ("Percent Difference = " & Abs(Eflux - q / eps0) / Eflux * 100)
Stop
End Sub