There is a group of two-dimensional radiative transfer problems for which an exact solution can be found. These problems are connected with non-scattering media with the following types of boundary radiation: (1) cosine varying collimated radiation, (2) strip of collimated radiation, (3) cosine varying diffuse radiation, and (4) constant temperature strip. In these cases the two-dimensional problem can be reduced to one-dimensional integral equations by the method of separation of variables. These problems are considered in a series of papers by Breig and Crosbie (Breig, W.F. and Crosbie, A.L. Two-dimensional radiative equilibrium. J. Math. Anal. Appl., 1974, 46, 104-125.) where also a good review of literature on the subject is given. Their approach allowed to determine only the external radiation field.

Here we generalize the results of Breig and Crosbie by applying the method of approximating the kernel of the integral equation for the Sobolev resolvent function (which essentially is the regular part of the respective Green function) by a series of exponents. The resulting approximate equation has an exact solution which is also represented by a series of exponents. This allows us to define the auxiliary functions g and h through the resolvent function F and thus define all the relevant functions. 

The following is given in detail in: T. Viik, Temperature Distribution in a Semi-Infinite Atmosphere Subjected to Cosine Varying Collimated Radiation, Proc. Estonian Acad. Sci. Phys. Math., 49, 2000, pp.40-57.

2  SOLUTION OF THE EQUATION OF RADIATIVE TRANSFER

We are looking for the emissive power  in a homogeneous non-scattering plane-parallel two-dimensional gray atmosphere which is in local thermodynamic equilibrium. The radiative transfer in such an atmosphere is described by the following equation

cosq ¶I ¶tz +sinqsinf ¶I ¶ty +I= s p T4,

(1)

where I is the intensity; q, the polar angle measured from the inward normal to the atmosphere; f,azimuthal angle measured from the t_x-axis; s, the Stefan-Boltzmann constant; T, the temperature in the atmosphere; sT⁴,the emissive power. The optical depth t_z is measured downward from the boundary of the atmosphere and together with  t_x and t_y they form a right-hand rectangular co-ordinate system. We require that the energy is transferred only by radiation, i.e. there is no heat conduction or convection in the atmosphere.

Applying of integrating factor techniques to Eq. (1) we obtain the formal solution for the intensities of downward and upward moving radiation in the form

I+(ty,tz,m)=I0(ty+)exp(-tz/m)+ 1 p ó õ tz 0  sT4(ty¢,tz¢)exp(-(tz-tz¢)/m)dtz¢/m

(2)

and

I-(ty,tz,m)= 1 p ó õ ¥ tz  sT4(ty¢,tz¢)exp(-(tz¢-tz)/m)dtz¢/m,

(3)

where t_y⁺=t_y-t_ztanqsinf, t_y^¢=t_y+(t_z^¢-t_z)tanqsinf, m = cosq and I₀⁺ is the intensity incident on the boundary of the atmosphere [ ¹] .

As far as we require the atmosphere to be in radiative equilibrium we may write

4sT4(ty,tz)= ó õ 4p  Idw,

(4)

where w is the solid angle.

Substituting Eqs. (2) and (3) into Eq. (4) we obtain the equation for the emissive power

4sT4(ty,tz) = ó õ 2p  I0+(ty+)exp(-tz/m)dw+ 1 p ó õ 2p 0  ó õ 1 0  ó õ ¥ 0  sT4(ty¢,tz¢)exp(-| tz-tz¢| /m)dtz¢dm/mdf.       (5)

According to our assumption the incident intensity may be expressed as

I0+(ty+)=I0[ 1+eexp(ibty+)] d(m-m0)d(f),

(6)

where I₀ is a constant, (m₀=cosq₀,f) defines the direction of the incident collimated radiation, e is the amplitude of the cosine wave and d is the Dirac delta function. Boundary condition (6) means that the top of the atmosphere is illuminated strip-wise by a parallel beam at an angle q₀ while the strips are parallel to the x - axis and their widths are defined by the spatial frequency b as p/b in units of optical length t_y⁺. The illumination in the direction parallel to the y - axis varies according to the cosine law. Next we apply the concept of separation of variables to Eq. (5) by assuming that

sT4(ty,tz)= 1 4 I0[ Bb = 0(t,m0)+eBb(t,m0)exp(ibty)] ,

(7)

where B_b is the dimensionless emissive power and t = t_z. Using Eq. (7) in Eq. (5) gives us a simple integral equation for B_b in the form

Bb(t,m0)=exp(-t/m0)+ 1 2 ó õ ¥ 0  E1(t-t¢)Bb(t¢,m0)dt¢,

(8)

where the generalized exponential integral E₁ is defined as [ ²]

E1(t,b)= ó õ ¥ 1  exp(-| t| Ö   t2+b2)   dt Ö t2+b2 .

(9)

By substituting m₀ for (t²+b²)^-1/2 in Eq. (8) and multiplying both sides of it by dt/Ö{t²+b²}, and last, integrating from 1 to ¥, we arrive at the integral equation for the resolvent function F_b in the form

Fb(t)= 1 2 E1(t,b)+ 1 2 ó õ ¥ 0  E1(t-t¢)Fb(t¢)dt¢,

(10)

where

Fb(t)= 1 2 ó õ ¥ 1  Bb æ è t, Ö   t2+b2   ö ø dt Ö t2+b2 .

(11)

Next we introduce two functions h(t,m) and g(t,m) as follows [ ¹³]

hb(t,m)=1+ ó õ ¥ t  Fb(t)exp(-(t-t)/m)dt

(12)

and

gb(t,m)=exp(-t/m)+ ó õ t 0  Fb(t)exp(-(t-t)/m)dt.

(13)

In the following we need a system of equations which connect those two functions with each other

-m ¶hb(t,m) ¶t +hb(t,m) = mFb(t)+1, (14) m ¶gb(t,m) ¶t +gb(t,m) = mFb(t). (15)

Eqs. (14) and (15) can easily be found from Eqs. (12) and (13) by differentiating them with respect to t.

Sobolev  [ ¹⁸]  has shown that the solution of Eq. (8) may be written in the form

Bb(t,m0)=Bb(0,m0) é ë exp(-t/m0)+ ó õ t 0  Fb(t)exp(-(t-t)/m0)dt ù û ,

(16)

or, in our notation,

Bb(t,m0)=Bb(0,m0)gb( t,m0)

(17)

Formally this completes the solution of the problem to determine the temperature distribution in a semi-infinite atmosphere subjected to collimated cosine varying radiation.

Next we show how to find the emissive power at the boundary B_b(0,m₀) and the function g_b(t,m) at an arbitrary optical depth.

It is obvious that if b = 0 then Eq. (8) reduces to the equation describing radiation transfer in an one-dimensional medium which have been successfully solved by introducing the Sobolev's resolvent function [²] and then approximating it by a sum of exponents. Since Eq. (8) is linear and the kernel is a sum of exponents we may try to use the same technique.

First we change the variable u=(t²+b²)^-1/2 in Eq. (9) to reduce this formula to a more familiar form

E1(t,b)= ó õ p 0  exp(-t/u) Ö 1-b2u2 du u ,

(18)

where p=(1+b²)^-1/2.

To solve Eq. (10) we express the generalized exponent integral in Eq. (18) as a sum of exponents

E1(t,b)=2 N å k=1  wkmk-1Ykexp(-t/mk),

(19)

where the characteristic function is expressed as

Yk= 1 2 Ö 1-b2mk2 .

(20)

In Eq. (19) w_i and m_i are the weights and points of a Gauss quadrature rule in the interval ( 0,p) and N is the order of the quadrature [ ³] . The characteristic function Y (different from that which appears in the analysis by Breig and Crosbie [ ⁷] but nevertheless giving accurate results!) is not a polynomial but it has retained another important quality - it still is an even function of x.

If we have approximated the general exponential integral as a sum of exponents then Eq. (10) accepts an exact solution as a sum of exponents [ ¹³]

Fb(t)= N å i=1  aiexp(-sit).

(21)

In order to determine the coefficients a_i and s_i in Eq. (21) we use Eq. (21) in Eq. (10) and by equating the similar exponents we obtain the characteristic equation

1-2 N å i=1  wiY(mi,b) 1-mi2s2 =0,

(22)

and a linear algebraic system for coefficients a_i

N å k=1  ak 1-misk -mi-1=0,       i=1,...,N.

(23)

It is evident that Eq. (22) has exactly N pairs of non-zero solutions ±s_k if only b ¹ 0. If b = 0 then s₁=±0 is also a solution but as far as this takes us back to the thoroughly studied one-dimensional case, we shall not consider it here.

The roots of the characteristic equation satisfy the following inequalities

0 £ | s1| < mN-1 < | s2| < mN-1-1 < ... < | sN| < m1-1.

As far as the roots are bracketed we may use any of the well-recommended root-finding algorithm, e.g. Brent's method [ ¹⁵] .

In our approximation the functions h_b(t,m) and g_b(t,m)  (Eqs. (12) and (13)) may be written in the form

hb(t,m)=1+m N å i=1  aiexp(-sit) 1+sim

(24)

and

gb(t,m)=exp(-t/m)+m N å i=1  ai[exp(-sit)-exp(-t/m)] 1-sim .

(25)

It may happen that in certain cases we observe the apparent singularity at s_im = 1 but it can simply be removed by substituting the respective term in the sum for a_itexp(-t/m).

Since the formula for the emissive power at the boundary is given by Sobolev [ ¹⁸] in the form

Bb(0,m0)=1+ ó õ ¥ 0  Fb(t)exp(-t/m0)dt,

(26)

we have in our approximation

Bb(0,m0)=Hb(m0)=1+m N å i=1  ai 1+sim0 .

(27)

This concludes the solution of Eq. (8).

3  RADIATIVE FLUX

In this section we consider the formulation of the equations for the z-component of radiative flux in the atmosphere and respective calculations. According to Breig and Crosbie the z-component of radiative flux can be shown to satisfy the relationship

qz(ty,t)=I0Qb = 0( t,m0)+eI0Qb( t,m0) exp(ibty), (28)

where the dimensionless radiative flux is given by

Qb( t,m0) = m0exp( -t/m0) + 1 2 ó õ t 0  E2( t-t¢,b) Bb(t¢,m0)dt¢- 1 2 ó õ ¥ t  E2( t¢-t,b) Bb(t¢,m0)dt¢. (29)

In Eq. (29) the generalized second exponential integral can be written as

E2( t,b) = ó õ p 0  exp( -| t| /u) du ( 1-b2u2) 3/2 .

(30)

Substituting Eq. (30) into Eq. (29), changing the order of integration and taking into account Eqs. (2), (3), (14) and (15) we obtain

Qb( t,m0) = m0exp( -t/m0) + m0Hb( m0) ó õ p 0  uy1(m,b)du m0-u [ gb( t,m0) -gb( t,u) ] - m0Hb( m0) ó õ p 0  uy1(m,b)du m0+u [ gb( t,m0)+hb(t,u)-1] , (31)

where

y1(m,b)= 1 2( 1-b2u2) 3/2 .

(32)

The radiative flux at the boundary of an atmosphere is thus

Qb( 0,m0) = m0-m0Hb( m0) ó õ p 0  uy1(m,b)Hb(u)du m0+u .

(33)