,

,

W1, the work per cycle expended to drive q1, is positive. But if q1 is decreased, say to q1 = .0001 q2, then the work per cycle to drive q1 is negative. There must be at least one value of q1 in the range .0001 q2 < q1 < q2 such that W1 = 0. At this value of q1 no net work per cycle is expended to make q1 move with the specified motion. Indeed if one constructs a surface around the vibrating q1 (but excluding q2), then the net flux of field energy per cycle time through this surface is zero. Fig. 1 plots W1 vs. q1 over the range .075 q2 < q1 < .077 q2.

Private Sub cmdRangeOfq1_Click()

'****************

'Compute W1, the work per cycle expended to drive the charge q1 sinusoidally, for

'a range of values of q1. Output W1 vs. q1 for plotting purposes.

'****************

'Physical and mathematical constants

Const c As Double = 299792000# 'Speed of light

Const epsilon0 As Double = 0.00000000000885 'Permittivity constant

Const pi As Double = 3.14159265358979

Const Steps As Long = 100

Const q1max As Double = 0.077

Const q1min As Double = 0.075

Const deltaq As Double = (q1max - q1min) / Steps

Const q2 As Double = 1

Const Amp As Double = 1 'Amplitude of each oscillation

Const omega As Double = 0.01 * c / Amp 'Angular frequency of each oscillation

Const freq As Double = omega / (2 * pi) 'frequency of each oscillation

Const tau As Double = 1 / freq 'period of each oscillation

Const lambda As Double = c * tau 'wavelength

Const deltat As Double = tau / Steps 'time interval between epochs

Const L As Double = 1 * lambda 'separation

'Variables

Dim W1(Steps), W2 As Double

Dim x1(Steps), x2(Steps) As Double 'current coordinates of q1 and q2

Dim v(Steps), gamma(Steps), a(Steps), dadt(Steps) As Double

Dim BigIndex, Index As Long 'Loop index

Dim t(Steps) As Double 'Current time

Dim tr As Double 'Retarded time

Dim dt As Double 't - tr

Dim dtmin As Double 'Minimum value for dt

Dim dtmax As Double 'Maximum value for dt

Dim xr2 As Double 'q2 Retarded positions (on x-axis)

Dim vrx2 As Double 'Retarded velocities

Dim arx2 As Double 'Retarded accelerations

Dim Drx2 As Double 'component of vector Dr

Dim Dr2 As Double 'magnitude of vector Dr

Dim ux2 As Double 'component of vector u

Dim Ex2(Steps) As Double 'x-component of electric field vector

Dim Fx1(Steps) As Double 'Interactive force on q1

Dim FRadReact1(Steps) As Double 'Radiation Reaction force on q1

Dim FAgent1(Steps) As Double 'Agent counteraction

Dim P1(Steps) As Double 'Agent expended power

Dim q1(Steps) As Double 'Values of q1

'Compute values for W1.

For BigIndex = 0 To Steps - 1

W1(BigIndex) = 0

Next BigIndex

For BigIndex = 0 To Steps - 1

'Set value of q1.

q1(BigIndex) = q1min + BigIndex * deltaq

Debug.Print q1(BigIndex)

For Index = 0 To Steps - 1

'Compute the current positions and derivatives for q1.

t(Index) = Index * deltat

x1(Index) = -L / 2 + Amp * Sin(omega * t(Index))

v(Index) = omega * Amp * Cos(omega * t(Index))

gamma(Index) = 1 / Sqr(1 - v(Index) ^ 2 / c ^ 2)

a(Index) = -omega ^ 2 * Amp * Sin(omega * t(Index))

dadt(Index) = -omega ^ 3 * Amp * Cos(omega * t(Index))

'Compute the radiation reaction force on q1.

FRadReact1(Index) = q1(BigIndex) ^ 2 / (6 * pi * epsilon0 * c ^ 3) * gamma(Index) ^ 4 * (dadt(Index) + 3 * gamma(Index) ^ 2 * v(Index) * a(Index) ^ 2 / c ^ 2)

'Compute the needed retarded quantities for q2.

dtmin = 0

dtmax = 5 * L / c

Do

dt = (dtmax + dtmin) / 2

tr = t(Index) - dt

xr2 = L / 2 + Amp * Sin(omega * tr)

Drx2 = x1(Index) - xr2

Dr2 = Abs(Drx2)

If Abs(c * dt - Dr2) < 2 ^ (-30) Then Exit Do

If c * dt - Dr2 > 0 Then

dtmax = dt

Else

dtmin = dt

End If

Loop

vrx2 = omega * Amp * Cos(omega * tr)

arx2 = -(omega ^ 2) * Amp * Sin(omega * tr)

'Compute the components of vector u.

ux2 = c * Drx2 / Dr2 - vrx2

'Compute the electric field components of q2 at q1

Ex2(Index) = q2 / (4 * pi * epsilon0) * Dr2 / (Drx2 * ux2) ^ 3 * (ux2 * (c ^ 2 - vrx2 ^ 2))

Next Index

'Now compute the electric (interactive) force on q1.

For Index = 0 To Steps - 1

Fx1(Index) = q1(BigIndex) * Ex2(Index)

Next Index

'Then compute W1.

For Index = 0 To Steps - 1

FAgent1(Index) = -FRadReact1(Index) - Fx1(Index)

P1(Index) = FAgent1(Index) * v(Index)

Next Index

For Index = 0 To Steps - 1

W1(BigIndex) = W1(BigIndex) + P1(Index) * deltat

Next Index

Next BigIndex

'Output values of q1 and W1 for plotting purposes.

Open "c:\WINMCAD\ZeroRadiator.prn" For Output As #1

For BigIndex = 0 To Steps - 1

Write #1, q1(BigIndex), W1(BigIndex)

Next BigIndex

MsgBox ("Ready for plotting")

Close

Stop

End Sub