The diffusion, reacceleration, convection and loss terms in eq. (1) can all be finite-differenced for each dimension (R,z,p) in the form

∂ψi	=	ψit+Δt − ψit	=	α1ψi-1t+Δt  −  α2ψit+Δt  +  α3ψi+1t+Δt	+	qi
∂t		Δt		Δt

where all terms are functions of (R,z,p). In the Crank-Nicholson implicit method (Press et al., 1992) the updating scheme is

ψit+Δt =  ψit +  α1ψi-1t+Δt  −  α2ψit+Δt  +  α3ψi+1t+Δt  + qiΔt

Thus the tridiagonal system of equations,

−  α1ψi-1t+Δt  +  (1 + α2)ψit+Δt -  α3ψi+1t+Δt  = ψit + qiΔt

Process   Coordinate  α1 /Δt α2 /Δt α3 /Δt  Diffusion  R  Dxx 2Ri – ΔR 2Ri (ΔR)2   Dxx 2Ri 2Ri (ΔR)2   Dxx  2Ri + ΔR  2Ri (ΔR)2    z  Dxx /(Δz)2   2Dxx /(Δz)2   Dxx /(Δz)2   Convectiona  z>0(V>0)  V(zi-1 )/Δz   V(zi )/Δz  0 z<0(V<0) 0  − V(zi )/Δz  − V(zi+1 )/Δz   p (  dV  >0 )  dz  0  1  dV  pi  3  dz  Pii+1  1  dV  pi+1  3  dz  Pii+1  Diffusive  reaccelerationa  p   −  Dpp,i − Dpp,i-1  + Pii−12 +   2  (   Dpp,i  +   Dpp,i−1 ) Pii−1  Pi−1i+1 pi–1   −  Dpp,i − Dpp,i−1  + Pii−12 +  2Dpp,i (    1  +  1 ) +  Pi−1i+1  Pii+1 Pii−1  +  2Dpp,i Pii−1 pi 2Dpp,i+1 Pi−1i+1 Pii+1    Energy lossa p  0 −    1  dpi Pii+1 dt −    1  dpi+1 Pii+1 dt  Fragmentation  R,z,p  0  1/3 τf  0  Radioactive  decay R,z,p  0  1/3 τr  0

a Pji ≡  pi – pj

As an example, the coefficients for the radial term are derived here.

1

∂

(

RDxx

∂ψ

)

=

2

Dxx

{

Ri+1

ψi+1– ψi

–  Ri-1

ψi– ψi-1

}

R

∂R

∂R

Ri

Ri+1– Ri-1

Ri+1– Ri

Ri– Ri-1

Setting  R_i+1 − R_i  = R_i − R_i-1  = ΔR, one can obtain the following expressions in terms of our standard form (eq. [1])

α1/Δt =  Dxx   2Ri – ΔR  2Ri (ΔR)2

α2/Δt = Dxx  2Ri  2Ri (ΔR)2

α3/Δt =  Dxx   2Ri + ΔR  2Ri (ΔR)2

In terms of 3-D momentum phase-space density ƒ(p) the diffusive reacceleration equation is

∂ƒ(p)

= divp •

[

Dppgradp ƒ(p)

]

=

1

∂

[

p2Dpp

∂ƒ(p)

]

∂t

p2

∂p

∂p

The distribution is asumed isotropic so ƒ(p) = ƒ(p) where p = |p|. First we rewrite the equation in terms of ψ(p) = 4πp²ƒ(p) instead of ƒ(p) and expand the inner differential:

∂ψ

=

∂

[

p2Dpp

∂

ψ

]

=

∂

Dpp

[

∂ψ

−

2ψ

]

∂t

∂p

∂p

p2

∂p

∂p

p

The differencing scheme is then

2

[

Dpp,i+1

(

ψi+1– ψi

–

2ψi+1

)

–

Dpp,i-1

(

ψi– ψi-1

–

2 ψi-1

)]

pi+1– pi-1

pi+1– pi

pi+1

pi– pi-1

pi-1

In terms of our standard form (eq.[1]) the coefficients for reacceleration are

α1/Δt =	2 Dpp,i-1		(	1	+	2	)
	pi+1– pi-1			pi– pi-1		pi-1

α2/Δt =	2		(	Dpp,i+1	+	Dpp,i-1	)
	pi+1– pi-1			pi+1– pi		pi– pi-1

α3/Δt =	2 Dpp,i+1		(	1	+	2	)
	pi+1– pi-1			pi+1– pi		pi+1

One more scheme (#2) comes from further detalization

	dψ	=	∂ Dpp	∂ψ	+	Dpp	∂2ψ	– 2	∂	Dppψ
	dt		∂p	∂p			∂p2		∂p	p

Here it is

α1/Δt =  –	Dpp,i–  Dpp,i-1		+	2 Dpp,i	+	2 Dpp,i-1
	( pi– pi-1 ) 2			( pi+1– pi-1)  ( pi– pi-1)		( pi– pi-1) pi-1

α2/Δt =  –	Dpp,i–  Dpp,i-1		+	2 Dpp,i	(	1	+	1	)	+	2 Dpp,i
	( pi– pi-1 ) 2			pi+1– pi-1		pi+1– pi		pi– pi-1			( pi– pi-1) pi

α3/Δt =	2 Dpp,i+1
	( pi+1– pi-1)  ( pi+1– pi)

|  Home  |  Manual  |  Models  |  Theory  |  Forms  |  Results  |  News  |  References  |  Contact us  |