Download Discrete-Ordinate Method for Radiative Transfer: Theory and Numerical Implementation and more Study notes Environmental Science in PDF only on Docsity!
Lecture 20
Methods for solving the radiative transfer equation.
Part 3: Discrete-ordinate method.
Objectives:
- Discrete-ordinate method for the case of isotropic scattering.
2.Generalization of the discrete-ordinate method for inhomogeneous atmosphere.
- Numerical implementation of the discrete-ordinate method: DISORT
- Examples: Modeling with DISORT.
Required reading :
L02: 6.
Additional/Advanced reading :
Thomas G.E. and K. Stamnes, Radiative transfer in the atmosphere and ocean, 2000,
Chapter 8.1-8.
1. Discrete-ordinate method for the case of isotropic scattering. - A discrete-ordinate method has been developed by Chandrasekhar in about 1950
(see Chandrasekhar S., Radiative transfer, 1960, Dover Publications).
Recall the radiative transfer equation (Eq.17.7) for azimuthally independent diffuse
intensity:
( , )exp( / ) 4
0 0 0
0
1
1
0
−
I I P d F P d
dI
For isotropic scattering, the scattering phase function is 1. Hence we have
exp( / ) 4
0
0 0
1
1
0 τ μ π
ω τ μ μ
ω τ μ τ
τ μ
−
F
I I d d
dI [20.1]
Let’s apply the Gaussian quadratures to replace the integral in Eq.[20.1]
exp( / ) 4
0
0 0 0
=−
F
I a I d
dI n
j n
i j j
i i [20.2]
Inhomogeneous part
where i=-n, …,n (2n terms) and aj are the Gaussian weights (constants) and μ j are
quadrature angles (or points).
Eq.[20.2] is a system of 2n inhomogeneous differential equations:
Solution of Eq.[20.2] = general solution + particular solution
where the general solution is a solution of the homogeneous part of the Eq.[20.2]
Denoting Ii = Ii ( τ,μ i ) , the general solution of Eq.[20.2] can be found as
I (^) i = gi exp(− k τ ) [20.3]
Inserting Eq.[20.3] into Eq.[20.2], we obtain
=−
n
j n
g (^) i ik ajgj 2
( 1 )
ω 0 μ [20.4]
We can find gi in the form
g i = L /( 1 + μ ik )
where L is a constant to be determined. Substituting this expression for gi in Eq.[20.4],
we have
2 2 1
0
0 2 1 1
k
a
k
a
j
j
n
i j
j
n
j n^ μ
ω μ
ω
−
=− =
[20.5]
Eq.[20.5] gives 2n solutions for +/-kj (j=1,…,n).
Thus general solution is
exp( ) 1
j i j
j
j
i k k
L
I −
= ∑ [20.6]
where Lj are constants.
The particular solution can be found as
exp( / ) 4
0
0 0
i =^ hi −
F
I [20.7]
where hi are constants.
Inserting Eq.[20.7] into Eq.[20.2], we have
( 1 + / 0 )=^0 ∑ +
=−
n
j n
hi i ajhj
μ μ [20.8]
So we need to solve
( ) ( ) exp( / ) 4
0 0
0
1
1
0 0 ,
ϖ μ μ τ μ π
ω
ϖ μ μ τ μ μ
ω τ μ δ τ
τ μ μ
=
= −
P P F
I P P I d d
dI
o
m l
N
l m
m l
m l
m m l
m l
m l
N
l m
m
m
m
The general solution may be written
=−
n
j n
m j j
m j
m i j
m
I (τ,μ ) L φ (μ )exp( k τ)
m φ (^) j , , are coefficients to be determined. m k (^) j m Lj
The particular solution may be written
( τ ,μ)= (μ i )exp(−τ/ μ 0 )
m i
m Ip Z
( (^) i )
m Z μ is a function
( )
1
1 /
( ) ( ) ( ) 4
1 ( ) 0 0 0
0 0 0 0 0 i
m l
m l
N
l
m l i
m m m i m
m P
H H Z FP μ μ
ϖ ζ μ μ
μ μ ω μ π
− = −
The complete solution of the radiative transfer is
∑^ [20.14]
=−
n
j n
i
m m j j
m j
m i j
m
I (τ ,μ) L φ (μ )exp( k τ) Z (μ)exp( τ/ μ 0 )
i=-n,…,n
Let’s generalize the complete solution Eq.[20.14] of the radiative transfer for the
inhomogeneous atmosphere. The atmosphere can be divided into the N homogeneous
layers, each is characterized by a single scattering albedo, phase function, and optical
depth.
NOTE : If an atmospheric layer has gases, aerosols and/or clouds, one needs to calculate
the effective optical properties of this layer (see Eqs.[14.31]-[14.34]).
For l -th layer, we can write the solution using Eq.[20.14]. To simplify notations, let’s
consider the azimuthal independent case (i.e., m=0 ), so we have
∑^ [20.15] =−
n
j n
i
l l j j
l j
l i j
l I (τ ,μ ) L φ (μ )exp( k τ) Z (μ )exp( τ/ μ 0 )
Now, we need to match the boundary and continuity conditions between layers.
At the top of the atmosphere (TOA): no downward diffuse intensity
( 0 , ) 0 [20.16]
1 − = = i
l I μ
At the layer’s boundary : upward and downward intensities must be continuous
( , ) ( , ) [20.17]
1 l i
l l i
l I τ μ I τ μ
At the bottom of the atmosphere (assuming the Lamdertian surface):
( , ) [ (τ ) μ 0 0 exp(τ / μ 0 )] π
τ μ N N sur N i
l N
F F
r
I = + −
= ↓ [20.18]
Eqs.[20.16]-[20.18] provide necessary equations to find the unknown coefficients.
3. Numerical implementation of the discrete-ordinate method: DISORT
DISORT is a FORTRAN numerical code based on the discrete-ordinate method
developed by Stamnes, Wiscombe et al.
DISORT is openly available and has a good user-guide.
Some features:
DISORT applies to the inhomogeneous nonithothermal plane-parallel atmosphere.
A user may set-up any numbers of the plane-parallel layers.
Each layer must be characterized by the effective optical depth, single scattering
albedo and asymmetry parameter if the Henyey-Greenstein phase function is used.
- A user may use any phase function by providing the Legendre polynomial expansion
coefficients.
- A user selects a number of streams (keeping in mind that the computation time varies
as n 3 ).
A key problem is to obtain a solution for fluxes for strongly forward-peaked scattering.
DISORT allows predicting the intensity as a function of the direction and position at
any point in the atmosphere (i.e., not only at the boundaries of the layers).
