The twostream equations with scattering

The two stream approximations to the full scattering equation are derived from Eq. 5.14 and Eq. 5.15 by constraining the angular distribution of the radiation in such a way as to allow all integrals appearing in these equations to be written in terms of either I+ + I_ or I+ — I_. In the resulting equations, flux in the upward stream is absorbed, or scattered into the downward stream, at a rate proportional to the upward stream intensity, and similarly for the downward stream. The two-stream approximations are an instance of what physicists euphemistically like to call "uncontrolled approximations," in that they are not actually exact in any useful limit but are nonetheless physically justifiable and perform reasonably well in comparison to more precise calculations. The two-stream approximiations have inevitable inaccuracies because it is not, in fact, possible to precisely determine the scattering or the absorption from knowledge of the upward and downward fluxes alone. The two-stream approximations can be thought of as the first term in a sequence of N-stream approximations which become exact as N gets large. Fortunately, N = 2 proves sufficiently accurate for most climate problems.

The general form of a two-stream approximation for diffuse radiation is where t* is the optical depth in the vertical direction (increasing upward) including scattering loss. The coefficients Yj depend on frequency, on the properties of the scatterers, and on the particular assumption about angular distribution of radiation that was made in order to derive an approximate two-stream form from the full angular-resolved equation. We recover the hemispherically isotropic Schwartzschild equations used in previous chapters by taking y2 = 7+ = 7- =0 and 71 = yb = 2. The terms proportional to yb represent the source due to thermal emission of radiation, while the terms proportional to 7+ and y_ represent the source of diffuse radiation caused by scattering of the direct beam. The direct beam is assumed to have flux Lq (generally the solar constant) in the direction of travel, and to travel at an angle Z relative to the vertical. There is no upward direct beam term because it is assumed that all direct beam flux scattered from the ground scatters into diffuse radiation.

The symmetry between the coefficients multiplying I+ and I_ is dictated by the requirement that the equations be invariant in form when one exchanges the upward and downward directions. 71 gives the rate at which flux is lost from the upward or downward radiation, while 72 gives the rate of conversion between upward and downward radiation by scattering. We can derive additional constraints on the 7j. Subtracting the two equations gives us the equation for net vertical flux dL(1+ - 1-) = -(7i - 72)(1+ + I_) + 2ybnB + +(7+ + y_)L0 exp(-(r^ - t*)/cos Z) (5.28)

First, we demand that in the absence of a direct-beam source, the fluxes reduce to black body radiation in the limit of an infinite isothermal medium. Since B is constant and L0 = 0 in this case, we may assume the derivative on the left hand side to vanish. Since I+ = 1— = nB for blackbody radiation, we find yb = 71 — 72. Comparing with Eq. 5.14 we also find that 7+ + 7— = vo. To further exploit Eq. 5.14 we must approximate f Idl as being proportional to I+ + I_. The constant of proportionality, which we shall call 27', depends on the angular distribution of radiation assumed. With this approximation it follows that 71 — 72 = 27'(1 — vo).

Exercise 5.5.1 As a check on the above reasoning, show that in the conservative scattering limit = 1 the sum of the diffuse vertical flux with the direct beam vertical flux is constant.

Next, we sum the equations for I+ and I_ to obtain

^(I+ + I_) = -(71 + 72)(I+ - I_) + (7+ - 7_)L0 exp(-(t^ - t*)/cos Z) (5.29)

This can be compared to the symmetric flux projection given by Eq. 5.15. Making an assumption about the angular distribution allows one to approximate JIH(cos 0) cos 0di as proportional to I+ +1_, and J IHdl and f IHdll each as being proportional to I+ -1_. In consequence, 71 + 72 = 27•(1-g^o), where 7 is related to the proportionality coefficient and g is a coefficient characterizing the asymmetry of the scattering. If H = cos 0, then g is in fact the asymmetry factor g defined by Eq. 5.17, but other forms of H yield somewhat different asymmetry factors, though these tend to be reasonably close to g. For example, with the form of H given by Eq. 5.20, we showed that g = 2g for phase functions that are truncated to their first three Fourier components. Finally, under circumstances when we can write G = gcos 0, it follows that 7+ - 7_ = -27vog. This relation holds exactly when H = cos 0, and imposing it for other forms of H introduces errors that are no worse than other errors that are inevitable in reducing the full scattering equation down to two streams.

The general form of the set of two-stream coefficients satisfying all the above constraints is then

The coefficients 7 and 7' are purely numerical factors that depend on the assumption about the angular distribution of radiation which is used to close the two-stream problem. All vertical dependence then comes in through vo, and possibly through g if the asymmetry properties of scattering

Hemi-isotropic

Eddington

2 V3

Table 5.3: Coefficients for various two-stream approximations particles vary with height. There are three common closures in use. The first is the hemispherically isotropic closure, which we used earlier in deriving the two-stream equations without scattering. In this closure, it is assumed that the flux is isotropic (i.e. I is constant) in each of the upward and downward hemispheres, but with a different value in each hemisphere. The hemi-isotropic closure is derived by using the weighting function H defined by Eq. 5.20 and making use of Eq. 5.21. Given the isotropy of the blackbody source term, it is generally believed that the hemi-isotropic approximation is most appropriate for thermal infrared problems, with or without scattering. Another widely used closure is the Eddington approximation. The Eddington closure is obtained by taking H = cos 0 and making use of Eq. 5.19. To complete the closure, /Icos2 is written in terms of I+ + I_ by assuming that the flux is truncated to the first two Fourier components, so I = a + b cos 0. This is probably the most widely used closure for dealing with solar radiation. It is generally believed that this closure is a good choice for dealing with both Rayleigh scattering and the highly forward-peaked scattering due to cloud particles, though the mathematical justification for this belief is not very firm. The quadrature approximation is similar, except that J" I cos2 is evaluated using a technique known as Gaussian quadrature, which yields a different proportionality constant from the Eddington closure. The defining coefficients for the three closures are given in Table 5.3.

When = 0 there is no scattering and so the upward and downward streams should become uncoupled. From Eq. 5.30 we see that this decoupling happens only if 7 = 7', a requirement that is satisfied for the hemi-isotropic and quadrature approximations but not for the Eddington approximation. It follows that the Eddington approximation can incur serious errors when the scattering is weak, though it can nonetheless outperform the other approximations when scattering is comparable to or dominant over absorption.

The two-stream equations form a coupled system of ordinary differential equations in two dependent variables. Therefore they require two boundary conditions. At the top of the atmosphere, there is generally no incoming diffuse radiation, so the boundary condition there is simply I_ =0. At the bottom boundary we require that the upward diffuse radiation be the sum of the upward emission from the ground with the reflected direct beam and downward diffuse radiation. In general, the direct beam reflection might in part yield a reflected direct beam (as in reflection from a mirror-like smooth surface), but in the following we'll assume that all reflection from the bottom boundary is diffuse. Thus, the boundary condition at t* =0 is

I+(0) = egnB(v,Ts) + agLe cosZexp(—t*/ cosZ) + agI_(0) (5.31)

where eg is the emissivity of the ground and ag is the albedo of the ground, both of which vary with v; Kirchoff's law implies that eg = (1 — ag) for any given frequency.