Gauss’ Law: Electrodynamically Correct?
G.R.Dixon, July 16, 2010
In this article the electric field flux of a relativistically oscillating point charge inside a surrounding spherical surface is investigated. The objective is to determine whether Gauss’ law is satisfied when the particle is at different epochs in its oscillation cycle.
Many elementary texts demonstrate that Gauss’ law is satisfied when the charge is at rest in various locations inside the spherical surface. This is shown to be a natural consequence of the charge’s 1/r2 electric field dependence. However, it is not obvious (at least to the author) that the same will be true when the charge oscillates with a maximum speed on the order of c. The electric field of such an oscillating charge deviates significantly from the field of a resting point charge.
Of course if a point charge moves in some prescribed way along the x axis, then its E field can always be computed at points on the spherical surface. Having done this, the flux can readily be derived. Fig. 1 depicts flux increments (d
FE) through area increments (rings) of a spherical surface of radius R=5 meters. The charge’s motion is
(1)
The time is t=0. theta is the azimuth angle. Note that q’s maximum speed is highly relativistic. (The data points were computed by the first Visual Basic routine in the Appendix to this article.)
Figure 1

d
FE , Oscillating Charge, t=0The net flux through the spherical surface is the sum of the values of d
FE in Fig. 1. Its value is computed to be 112994350.Now according to Gauss’ law,
o. (2)
In the present case q / eo = 112994350. Thus Gauss’ law is also satisfied by the oscillating point charge, despite the fact that its field is not electrostatic.
But what about other instants in time, when the charge is not at the center of the sphere? Figs. 2-4 repeat Fig. 1 for times t=T/8, t=T/4, and t=3T/8 (where T is the period of oscillation). In every case the computed flux of E differs from the "Coulomb flux" at most by 1 (i.e., less than 1.5E-6 percent) ... a difference possibly attributable to numerical error.
Figure 2

Oscillating Charge, t=T/8
Figure 3

Oscillating Charge, t=T/4
Figure 4

Oscillating Charge, t=3T/8
Since every periodic motion can be represented as a sum of SHMs, it might be wondered whether Gauss’ law is obeyed for a particle whose motion consists of such a sum. Fig. 5 depicts the motion of a particle with frequencies
w1 = .99c/4A, w2 = .99c/2A, and w3 = .99c/A. Fig. 6 depicts the incremental fluxes. (The curves were generated by the second routine in the Appendix.) Note that the sum of the incremental fluxes again nearly equals q / eo, indicating that here too Gauss’ law may be satisfied. (The discrepancy is only 1.2E-5 percent). Indeed, this result suggests that Gauss’ law should be satisfied for every periodic motion of q.Figure 5

Motion of Multiple Frequency Oscillating q
Figure 6

d
FE , Oscillating Charge, Multiple Frequencies, t=0It is perhaps advisable to conclude by briefly discussing what this study does not prove. In electrostatics, Gauss’ law holds for the flux through any surface and not only a spherical one. Thus the present study does not demonstrate the universal applicability of Gauss’ law in electrodynamic cases. Nevertheless it is instructive that the law does appear to apply for all motions within a spherical surface ... a result not always made clear in electromagnetic texts.
***Appendix***
Option Explicit
Private Sub cmdComputeEFlux_Click()
'*******************
'Compute the E Field Flux through a 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 r As Double = 5 'radius of surrounding surface
Const omega As Double = 0.99 * c / A 'max q speed is highly relativistic
Const steps As Long = 5000 'number of iterations
Const dtheta As Double = pi / steps 'aximuth angle increment
Const freq As Double = omega / (2 * pi)
Const tau As Double = 1 / freq
Const deltat As Double = tau / steps 'time between epochs
Dim j As Long 'loop counter
Dim t(steps) As Double 'time
Dim theta(steps) As Double 'aximuth angle
Dim dArea As Double 'Area increment of ring
Dim Px, Py As Double
Dim Ex, Ey As Double 'electric field components
Dim Enormal(steps) As Double
Dim nx, ny 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 vr As Double
Dim ar As Double
Dim Eflux As Double
Dim Flux(steps) As Double 'Flux increments
'Compute and plot x(t)
For j = 0 To steps - 1
t(j) = j * deltat
x(j) = A * Sin(omega * t(j))
Next j
Open "c:\\WINMCADC\Physics\GaussTest.PRN" For Output As #1
For j = 0 To steps - 1
Write #1, t(j), x(j)
Next j
Close
MsgBox ("Ready for plotting x(t)")
Eflux = 0
For j = 0 To steps - 1
theta(j) = (j + 0.5) * dtheta
nx = Cos(theta(j))
ny = Sin(theta(j))
Px = r * Cos(theta(j))
Py = r * Sin(theta(j))
dArea = 2 * pi * Py * r * dtheta
dtmin = 0
dtmax = 5 * (r + A)
Do
dt = (dtmin + dtmax) / 2
'Choose one of the following, depending when the flux is to be computed.
'tr = t(0) - dt
tr = t(steps / 8) - dt
'tr = t(steps / 4) - dt
'tr = t(3 * steps / 8) - 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 = q / (4 * pi * eps0) * dr / (drx * ux + dry * uy) ^ 3 * (ux * (c ^ 2 - vr ^ 2) + dry * (-uy * ar))
Ey = q / (4 * pi * eps0) * dr / (drx * ux + dry * uy) ^ 3 * (uy * (c ^ 2 - vr ^ 2) - drx * (-uy * ar))
Enormal(j) = Ex * nx + Ey * ny
Flux(j) = Enormal(j) * dArea
Eflux = Eflux + Flux(j)
Next j
Open "c:\\WINMCADC\Physics\GaussTest.PRN" For Output As #1
For j = 0 To steps - 1
Write #1, theta(j) * 360 / 2 / pi, Flux(j)
Next j
Close
MsgBox ("Ready for plotting")
MsgBox ("E Flux = " & Eflux & "; Gauss = " & q / eps0)
MsgBox ("Percent Difference = " & Abs(Eflux - q / eps0) / Eflux * 100)
Stop
End Sub
Private Sub cmdTheorem1_Click()
'*******************
'Show that Gauss law works for a point charge moving periodically, at multiple frequencies, when the flux is
'computed over a spherical surface.
'*******************
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
Const A As Double = 1 'amplitude of oscillation equals 1 meter
Const r As Double = 5 'radius of surrounding surface
Const omega3 As Double = 0.99 * c / A
Const omega1 As Double = omega3 / 4
Const omega2 As Double = omega3 / 2
Const steps As Long = 5000 'number of iterations
Const dtheta As Double = pi / steps
Const freq1 As Double = omega1 / (2 * pi)
Const tau1 As Double = 1 / freq1
Const deltat As Double = tau1 / steps
Dim j As Long
Dim x(steps) As Double
Dim t(steps) As Double
Dim theta(steps) As Double
Dim dArea As Double
Dim Px, Py As Double
Dim Ex(steps), Ey(steps) As Double 'electric field components
Dim ExTotal(steps), EyTotal(steps) As Double
Dim Enormal(steps) As Double
Dim nx, ny 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 xr1, xr2, xr3, xr As Double
Dim vr1, vr2, vr3, vr As Double
Dim ar1, ar2, ar3, ar As Double
Dim Eflux As Double
Dim Flux(steps) As Double 'Flux increments
Dim RestFlux As Double
'Compute and plot x(t)
For j = 0 To steps - 1
t(j) = j * deltat
x(j) = A * Sin(omega1 * t(j)) + A * Sin(omega2 * t(j)) + A * Sin(omega3 * t(j))
Next j
Open "c:\\WINMCADC\Physics\GaussTest.PRN" For Output As #1
For j = 0 To steps - 1
Write #1, t(j), x(j)
Next j
Close
MsgBox ("Ready for plotting x(t)")
Eflux = 0
For j = 0 To steps - 1
theta(j) = (j + 0.5) * dtheta
nx = Cos(theta(j))
ny = Sin(theta(j))
Px = r * Cos(theta(j))
Py = r * Sin(theta(j))
dArea = 2 * pi * Py * r * dtheta
dtmin = 0
dtmax = 5 * (r + A)
Do
dt = (dtmin + dtmax) / 2
tr = t(0) - dt
xr1 = A * Sin(omega1 * tr)
xr2 = A * Sin(omega2 * tr)
xr3 = A * Sin(omega3 * tr)
xr = xr1 + xr2 + xr3
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
vr1 = omega1 * A * Cos(omega1 * tr)
vr2 = omega2 * A * Cos(omega2 * tr)
vr3 = omega3 * A * Cos(omega3 * tr)
vr = vr1 + vr2 + vr3
ar1 = -(omega1 ^ 2) * A * Sin(omega1 * tr)
ar2 = -(omega2 ^ 2) * A * Sin(omega2 * tr)
ar3 = -(omega3 ^ 2) * A * Sin(omega3 * tr)
ar = ar1 + ar2 + ar3
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) * nx + Ey(j) * ny
Flux(j) = Enormal(j) * dArea
Eflux = Eflux + Flux(j)
Next j
Open "c:\\WINMCADC\Physics\GaussTest.PRN" For Output As #1
For j = 0 To steps - 1
Write #1, theta(j) * 360 / 2 / pi, Flux(j)
Next j
Close
MsgBox ("Ready for plotting")
MsgBox ("E Flux = " & Eflux & "; Gauss = " & q / eps0)
MsgBox ("Percent Difference = " & Abs(Eflux - q / eps0) / Eflux * 100)
Stop
End Sub