<< Click to Display Table of Contents >> second_order_time |
{ SECOND_ORDER_TIME.PDE
This example shows the integration of Bessel's Equation as a test of the
time integration capabilities of FlexPDE.
Bessel's Equation for order zero can be written as
t^2*dtt(w) + t*dt(w) + t^2*w = 0
Dividing by t^2 and avoiding the pole at t=0, we can write
dtt(w) + dt(w)/t + w = 0
FlexPDE cannot directly integrate second order time equations, so we define an
auxiliary variable v=dt(w) and write a coupled pair of equations
dt(v) + v/t + w = 0
dt(w) = v
We use a dummy spatial grid of two cells and solve the equation at each node.
You can try varying the value given for ERRLIM to see how it behaves.
}
title "Integration of Bessel's Equation"
select
ngrid=1
errlim=1e-4 { increase accuracy to prevent accumulation of errors }
Variables
v (threshold=0.1)
w (threshold=0.1)
definitions
L = sqrt(2)
t0 = 0.001 { Start integration at t=0.001 }
Initial values { Initialize to known values at t=t0 }
w = 1-2.25*(t0/3)^2
v = -0.5*t0 + 0.5625*t0*(t0/3)^2
equations
v: dt(v) +v/t + w = 0
w: dt(w) = v
boundaries
region 1
start(-L,-L) line to (L,-L) to (L,L) to (-L,L) to close
time 0.001 to 4*pi { Exclude t=0 }
plots
for t=0.01 by 0.01 to 0.1 by 0.1 to 1 by 1 to endtime
history(w,bessj(0,t)) at (0,0) as "W(t) and BESSJ0(t)"
history(w-bessj(0,t)) at (0,0) as "Absolute Error"
history(v,-bessj(1,t)) at (0,0) as "V(t) and dt(BESSJ0(t))"
history(v+bessj(1,t)) at (0,0) as "Slope Error"
history(deltat)
end