Reciprocity relations in radiative transfer by spherical clouds

Journal of Quantitative Spectroscopy & Radiative Transfer 64 (2000) 151}172

Reciprocity relations in radiative transfer by spherical clouds H.C. van de Hulst Leiden Observatory, Postbus 9513, 2300 RA Leiden, Netherlands Received 1 August 1998

Abstract Time reversal i.e., the replacement of sources by detectors and conversely, leads in any con"guration of light scattering particles to mutually reciprocal problems, which have essentially the same solution. This principle is applied to a set of problems with spherically symmetric radiation "elds in a homogeneous spherical cloud. Internal or external sources and detectors are permitted. The individual scatterers may have arbitrary albedo a and an arbitrary phase function. The optical depth R along the radius of the cloud is also arbitrary. Each function belonging to a particular source}detector combination is de"ned as a quantity measurable in a "ctitious experiment. The principle of detailed balancing then is used to establish their reciprocity relations. Additional inter-relations are based on the conservation of energy and on volume integration. The relations permit surprising checks on analytical, numerical, or graphical results cited from papers treating very diverse problems. A set of examples for isotropic scattering is examined in detail and shows excellent agreement throughout. Among the possible extensions is a theorem on a homogeneous cloud of arbitrary shape with arbitrary phase function; this theorem relates the radiation "eld arising from uniform external illumination to the radiation "eld arising from uniformly embedded sources.  1999 Elsevier Science Ltd. All rights reserved.

1. Reciprocity as a tool Reciprocity in radiative transfer is a well-known concept (e.g. [1, Section 2.7; 2, Chapter 3]) but an underestimated tool. This paper addresses a variety of radiative transfer problems in a homogeneous spherical cloud. It uses the reciprocity principle to "nd rigorous cross-checks between the solutions of those problems. The idea is to maintain the con"guration of scattering and absorbing particles but to replace all (isotropic) sources by (isotropic) detectors and conversely. This e!ects a time reversal of all relevant light paths and de"nes two mutually reciprocal problems, which should have essentially the same answer. This &principle' is forged into a precise tool by the method of detailed balancing. The intriguing aspect of exploring these relations is that the quantities to be equated have often been computed for use in quite di!erent "elds of science and published in di!erent journals. Typical 0022-4073/99/$ - see front matter  1999 Elsevier Science Ltd. All rights reserved. PII: S 0 0 2 2 - 4 0 7 3 ( 9 9 ) 0 0 1 1 0 - 7

152

H.C. van de Hulst / Journal of Quantitative Spectroscopy & Radiative Transfer 64 (2000) 151}172

questions arising in problems with outside sources are: How much radiation penetrates to the center of the cloud? How much of it is sent back into space? The inspiration for such problems often comes from astronomy: a planet illuminated by the sun, an interstellar dust cloud exposed to the radiation from one or many stars, etc. Another class of problems assumes radiation sources inside the cloud. Frequently made assumptions are that there is just one source placed at the center of the cloud, or that isotropic sources are uniformly distributed through the cloud. Again, we seek the fate of the radiation from these sources. In some papers this is narrowed down to the question how much of it escapes. Such problems have served as practising problems for heat isolation and reactor physics but have also numerous applications in astronomy. In fact, the escape of radiation from a star, the grandmother of all radiative transfer problems, falls in this class. The tools used in this paper are simple but unconventional. The &classical' start with the equation of radiative transfer is not made. Instead, each function representing a physical quantity is de"ned, together with its normalization, as the outcome of a speci"c (idealized) physical experiment. We are convinced, that the properties of such a function do not change if we choose the lengthier but more familiar de"nition as the solution of a mathematical equation with given boundary conditions. The examples shown in Section 4 support this conviction. The paper is organized as follows. In Section 2 we select the set of problems we shall work with. The various source}detector con"gurations are arranged in a 5;5 checkerboard diagram (Fig. 1) and for each con"guration a function is de"ned expressing the net response of detectors to sources. Section 3 lists the interrelations between these functions. These are of four kinds. The reciprocity relations, the central topic of this paper, are proven in Section 3.1. Sections 3.2 and 3.3 treat the fairly trivial relations based on averaging and on integrals of the (external) spherical re#ection function. The not so obvious complementary relations, based on the conservation of energy are derived in Section 3.4. The correctness of these results is further supported by Section 4, which presents a number of analytical and numerical checks, based on examples cited from a variety of published papers. A brief listing of possible extensions (Section 5) and some "nal remarks (Section 6) end the paper.

2. A set of spherical problems 2.1. Terminology Since the terminology may vary with author or "eld of application, we brie#y list below the terms as we shall use them, mostly following Chandrasekhar [3]: E E E E E

Cloud: The con"guration of scattering particles. Position: Any point in space, usually inside the cloud. Volume element: A small volume around this position. Direction: Any direction in the full solid angle 4p. Intensity: Energy streaming at any position in any direction per unit solid angle per unit area in a wavelength interval chosen narrow enough to ensure that the &basic quantities' (see below) are constant in this interval. For simplicity, we assume that the intensity, incident and scattered, is unpolarized.

H.C. van de Hulst / Journal of Quantitative Spectroscopy & Radiative Transfer 64 (2000) 151}172

153

E Average intensity: A function of position, obtained by averaging the intensity over the full solid angle 4p. The average intensity is proportional to radiation density. If, in addition, a volume average over all positions in the cloud is taken, this will be explicitly stated. E Extinction coe$cient: i, the weakening of a unidirectional ray per unit length. The extinction coe$cient is the sum of the scattering coe$cient, ai, and the absorption coe$cient, (1!a)i, where a is the albedo (see below). E Albedo: Ratio scatt. coe!./ext. coe!. (04a41). E Scattering angle: a, the angle between an incident ray and any chosen direction of scattering. E Phase function: '(a), the function describing the intensity scattered in any direction and normalized to average 1 over the full solid angle 4p. E Basic quantities characterizing the scattering medium: We use this term for the three quantities i, a and '(a). Note 1. In a cloud with oriented particles, i and a depend on direction as well as position and the de"nition of '(a) as function of one angle no longer su$ces. We assume that such orientation e!ects are not present. Note 2. In the de"nition of intensity, the units have intentionally been left unde"ned. The reader may choose these at will. See introduction to Section 2.3. 2.2. Assumptions The set of conceivable transfer problems arising from all kinds of assumptions about the distributions of scatterers, sources, and detectors, is unlimited, and so is, in principle, the set of relevant reciprocity relations. After scanning the literature for papers in which a problem had been solved with su$cient precision to be eligible for a reciprocity check, we chose a limited set for a full exploration. This set of transfer problems is de"ned by the following assumptions: The scattering cloud: We assume the cloud to be spherical and homogeneous, with the &basic quantities' (see Section 2.1) constant throughtout the cloud. The radius is R, making the optical length along the radius R"iR.

(1)

No restriction is placed on the value of R. In the following, we shall avoid frequent multiplications by i by choosing i\ as the unit of length, so that optical length"geometrical length. Sources: All sources are isotropic. We consider a variety of arrangements of the sources; see Fig. 1. The basic choice is between line 1, where outside sources form a uniformly radiating wall surrounding the cloud at large distance, or line 2, where inside sources homogeneously cover a thin spherical shell with (geometrical and optical) radius r , with 04r 4R. The choices on the   remaining lines of Fig. 1 follow from line 2. Line 3 results from the speci"cation r "R, and  corresponds to a shell of sources covering the surface of the cloud. Line 5, by the speci"cation r "0, represents an isotropic source at the cloud center. Line 4 corresponds to sources homogene ously distributed through the volume of the cloud and has been placed as an intermediate case between the extremes of lines 3 and 5. It involves volume averaging; if f (r) is a function de"ned

154

H.C. van de Hulst / Journal of Quantitative Spectroscopy & Radiative Transfer 64 (2000) 151}172

Fig. 1. Diagram showing the selected set of problems and the function relevant to each particular combination of sources and detectors.

along the radius, then its volume average is



1 f 2"3R\

0 f (r)r dr. 

(2)

Detectors: All detectors are isotropic. The various con"gurations mimic those of the sources, but may be made fully independently, so that each square of the checker-board arrangement of Fig. 1 represents a separate problem. The basic choices are in columns a and b, and the derived choices in columns c, d, and e. This completes the set of assumptions but a few comments may be made. Only the four squares a1, b1, a2, and b2, in Fig. 1 are basic. The remaining 21 squares are derived. It is advisable to keep this diagram at hand in reading this paper just like it is easier to replay a chess game with a diagram at hand. Quantities found in squares which have mirror positions with respect to the main diagonal are related by the reciprocity principle (Section 3.1). In addition, we shall show (Section 3.4) that conservations of energy leads to &complementarity relations' between squares in the same column on lines 1 and 4 and between squares on the same line in columns a and d. The net e!ect is that up to four squares in Fig. 1 refer to directly related physical problems. The assumption of spherical symmetry may often be relaxed without a!ecting the quantities de"ned or the validity of their interrelations. For instance, it su$ces to assume that either the sources or the detectors have a spherically symmetric distribution. Also, the &surrounding wall' of line 1 or column a may have any form and be at any distance.

H.C. van de Hulst / Journal of Quantitative Spectroscopy & Radiative Transfer 64 (2000) 151}172

155

2.3. Dexnitions of the physical quantities One physical quantity in each square of Fig. 1 characterizes the (average) response of the detectors to the sources in that particular con"guration. We shall now proceed to de"ne these quantities, but a few general remarks in preparation to those de"nitions will be helpful. (1) Each quantity is de"ned as the dimensionless ratio of a numerator and a denominator of the same physical dimension. This renders the choice of units irrelevant, provided the same choice is employed in the numerator and in the denominator. As a consequence, the wording may be varied in wide limits. For instance, the word &radiation' in these de"nitions may also be replaced by &energy' or &number of photons'. Likewise, &average intensity' may be replaced by &radiation density' or by &number density of photons'. (2) In the interest of a simple notation, we anticipate the validity of reciprocity by employing the same notation for the quantities de"ned in squares in a reciprocal position, i.e., mirrored with respect to the main diagonal of Fig. 1. The need for separate, independent de"nitions for these two squares is not a!ected by this decision. We formulate these independent de"nitions below and postpone to Section 3.1 the proofs that both de"nitions lead indeed to identical mathematical functions. (3) Having given the de"nitions for the four basic squares of Fig. 1, we have to "nd the corresponding de"nitions for the derived squares. The substitutions r "R, or r"R, present no  problem. But the substitutions r "0, or r"0, lead to degeneracies, which we have to remove by  introducing separate de"nitions. 2.3.1. Dexnitions for the basic squares Square a1: Sources: a uniform, distant sky. radiation re#ected back to the sky . A"&Bond albedo'" radiation incident on the cloud

(3)

Square a2: Sources: a shell of isotropic sources at optical distance r from the center of the cloud.  radiation escaping to the sky . (4) P(r )"&escape probability'"  radiation emitted by the sources Square b1: Sources: (as in a1) a uniform, distant sky. average intensity at point r in cloud . P(r)"&gain'" average intensity anywhere in absence of cloud Square b2: Sources: (as in a2) shell of sources at r . We de"ne  aver. int. in cloud at distance r from center r . >(r, r )" )  r aver. int. at center of same shell source if cloud is removed 

(5)

(6)

2.3.2. Dexnitions for the derived squares We shall ignore squares c3 and e5, where >(r, r ) is in"nite. The remaining 19 squares are treated  line by line.

156

H.C. van de Hulst / Journal of Quantitative Spectroscopy & Radiative Transfer 64 (2000) 151}172

Line 1: Sources as in a1. No peculiarities. 1P2 is obtained by averaging according to Eq. (2). Line 2: A shell of sources, as in a2 and b2. The substitution in square c2 o!ers no problem. In square d2 we de"ne the new quantity R aver. int. with the shell source, averaged over cloud volume ) . (7) S(r )"  aver. int. at center of same shell without cloud 3r  The numerator of Eq. (7) is obtained by taking at each position the average intensity (see Section 2.1) and then averaging over the cloud volume by the recipe of Eq. (2). The denominator is the same as in Eq. (6). The front factor, which is dimensionless, since we work with optical lengths only, has been chosen to give the relation to P(r), Eq. (22) below, a simple form. On the "nal square of this line, e2, we need again a new quantity, for >(r, r ) becomes 0. We de"ne  aver. intensity at center with shell source at r  . (8) Z(r )"  aver. intensity at center of same shell without cloud Here the denominator is again the same as in the preceding de"nitions, and the numerator is an obvious choice. Line 3: No comments needed. Line 4: The sources are uniformly distributed through the cloud. Square a4 o!ers no problem, see Eq. (13) below. Square b4 requires a new de"nition. Call the emission (arbitrary units) by the embedded sources in a volume element d< inside the cloud per unit solid angle j d<. By de"nition,  j does not vary with direction or position. Each volume element also produces scattered radiation,  which in the same units may be written as j d<, where j is a function of direction and position.   We average j over all directions in a point at distance r from the center, and thus obtain the  function j (r). Now we are ready to de"ne  1 j (r) (9) S(r)" )  , a j  where the front factor ensures that squares b4 and d2 contain the same function (proof in Section 3.1). The functions in the remaining squares of line 4 follow from Eq. (9) by the substitutions r"R or r"0 (sq. c4 and e4), or by volume averaging according to Eq. (2) (sq. d4). Line 5: An isotropic source at the center of the cloud. The quantity in square a5 is obtained by putting r "0 in Eq. (4). In square b5 we enter:  aver. intensity at r in cloud with central source . (10) Z(r)" aver. intensity at same distance without cloud This same equation serves for square c5 by the substitution r"R. Square d5 requires a limiting procedure, for we wish to make sure that the S(0) entered here is the limit for r "0 of S(r ) de"ned   by Eq. (7). We start by stating that one isotropic source gives at distance l in vacuum the average intensity Al\, where A is a constant which we need not specify since it will cancel out. The correctness of this equation may be seen, e.g., by imagining the source to have a non-zero diameter. The result then is an intensity, which does not vary with l inside a solid angle which goes as l\, from which the equation for the average intensity follows at once. A number of C identical sources

H.C. van de Hulst / Journal of Quantitative Spectroscopy & Radiative Transfer 64 (2000) 151}172

157

give at distance l the average intensity CAl\. Applying this to Eq. (7) with l"r we see that the  r in the denominator of the front factor cancels the r\ in the remaining denominator. Hence, if we   let a shell containing C sources gradually shrink, the full denominator of (7) reaches at r "0 the  limit 3CA. The numerator of (7) applied to the shrinking source leads directly to the numerator of Eq. (11), below. In the denominator we have to choose a new reference value. The obvious choice, shown in (11), assumes the value CAR\. The front factor of (11) now automatically follows and the full de"nition becomes R aver. int. with central source, averaged over cloud volume . S(0)" ) aver. int. at distance R from same source in vacuum 3

(11)

3. Interrelations Any equation in this section is rigorous, satis"ed by the functions entered in Fig. 1, and valid for arbitrary cloud radius R, arbitrary but constant albedo a, and arbitrary but constant phase function. These relations fall, by physical content and by mathematical form, into distinct classes as follows. 1. 2. 3. 4.

Reciprocity relations connecting mirror squares across the main diagonal of Fig. 1. Integral relations derived from averaging over the cloud volume. Expressions in terms of bimoments of the spherical re#ection function. Complementarity relations between column a and column d, or lines 1 and 4 of Fig. 1.

3.1. Reciprocity relations This subsection forms the core of the present paper. It contains mostly proofs and only one new equation, Eq. (12). We have taken a mortgage on the reciprocity relations by already entering in mirror squares of Fig. 1 functions denoted by the same symbol. This mortgage will now be paid o! by proving that this was permitted. We use the time-honored method of detailed balancing: we equip the cloud with two sets of test spheres, which are su$ciently small not to bother the scattering process. First, set I is used as sources, and set II as detectors. In the reciprocal experiment, this choice is reversed. If we now assume that all emitting and absorbing test spheres are at the same temperature, detailed balancing should exist between the energy transferred in the "rst experiment from set I to set II, and in the reciprocal experiment from set II to set I. As the case may be, a set may be replaced by a single test sphere, or by a black wall surrounding the cloud. This is a heavy hammer to hit this nail but it works well. I feel excused from justifying the method in more detail, for it has been textbook material in astrophysics for at least 70 years, with the even older Kirchho! relation between emission and absorption as the most famous example. More theoretical background may be found in [4] and more about its application to multiple light scattering in Chapter 3 of Ref. [2]. We may prescribe any amount of detail in these "ctitious experiments. For instance, sources and detectors sensitive to radiation of a particular polarization, in a particular direction, or in a particular wavelength band, are not excluded in principle.

158

H.C. van de Hulst / Journal of Quantitative Spectroscopy & Radiative Transfer 64 (2000) 151}172

Such re"nements are not considered in this paper. This permits us to symbolize the sources and detectors by isotropic emitters and absorbers. We imagine each test sphere to be a small black sphere of radius c and temperature ¹, corresponding to Planck intensity B. Both c and B have to be adjusted to the chosen length unit i\, but this o!ers no problem, since they cancel out in each proof. As an emitter, one test sphere emits the energy pcB per unit time per unit solid angle and 4ppcB in all directions together. As an absorber, the same test sphere absorbs per unit time the energy 4ppcJ, where J is the local average intensity. We are now ready for the detailed proofs, in each case based on the requirement that the energy transferred from set I to set II in the "rst experiment equals the energy transferred from set II to set I in the reciprocal experiment. Squares a2 and b1, function P(R): We balance the energy exchange between a test sphere at distance r from the cloud center and a black outside wall. Both sides of the equation become 4pBpcP(r), where at one side P(r) is as de"ned in a2 and on the other side as de"ned in b1. This proves the identity. This is, in fact, a special case of the more general &gain' concept discussed in Ref. [2, p. 26]. It is also found implicitly in earlier texts (e.g. [5, p. 168]). Square b2, function >(r, r ): The statement to be proved is that this function is symmetric:  >(r , r)">(r, r ).  

(12)

In this proof and the ones which follow, we shall consistently write J for the average intensity, which depends on position only. Evidently, J has a di!erent form for each di!erent "ctitious experiment. Imagine at r a shell covered with N test spheres per unit area (set I) and at r a shell  covered with M test spheres per unit area (set II). Set I would produce in vacuum at the center J"NpcB. Hence, by De"nition (6), in the cloud at r it creates J"NpcB(r /r)>(r, r ). The   spheres of set II absorb together the energy 4prM4pJpc. So the stream from set I to set II is proportional to rr >(r, r ). This must equal the corresponding stream from set II to set I, which   proves Eq. (12). Squares b4 and d2, function S(r): Assume: set I is a shell with N test spheres per unit area at r. Set II is a homogeneous distribution of test spheres through the cloud with M spheres per unit volume. Set I, acting as sources, would give at the center in vacuum J"NpcB. Hence, by de"nition (7) it produces, averaged through the cloud J"NpcB(3r/R)S(r). This is picked up by each sphere of set II, making the full energy transfer I P II equal to 4pJpc(4pR/3)M. If set II is the emitter, the embedded sources emit per unit volume per unit solid angle j "MpcB. Observing that in  de"nition (9) we may write j "aJ, we "nd the J at position r to be J"MpcBS(r). This must be  multiplied by 4ppc4prN in order to "nd the net stream II P I. Equating the energy streams in the two experiments, we see that S(r) at both sides is the same function. Squares b5 and e2, function Z(r): Take for set I the shell at r from the preceding example and for set II one test sphere at the center. Set I as emitter produces at the center, by de"nition (8), J"NpcBZ(r). The energy transferred to set II is 4ppc times this value. The central sphere, set II, acting as emitter, is seen in free space from distance r under the solid angle pcr\. Since in free space the intensity is 0 in the remaining solid angle, J(r) in free space follows at once and by Eq. (10) we obtain J(r)"(pcB/4pr)Z(r) in the cloud. Multiplication with 4ppc4prN then gives the energy #owing from set II to set I. Equating the streams in both experiments, we "nd that the function Z(r) at both sides must be the same.

H.C. van de Hulst / Journal of Quantitative Spectroscopy & Radiative Transfer 64 (2000) 151}172

159

Squares d5 and e4, function S(0): A separate proof for this case is not really necessary, for Eq. (11) has been chosen to guarantee that S(0) is the limit of S(r ) for r "0. Hence, the balance proven for   S(r) is also valid in the limit r"0. However, we may also leave this consideration aside and perform the check as in the other cases. The result (which we do not write out) is again a full balance. 3.2. Integral relations Each of the functions appearing in column d is, in principle, obtained by volume averaging from the corresponding function in column b. This means an integration over r. The complete set of these integral relations is

  

0 3 P(r)r dr, 1P2" ) R 

(13)

0 3 1S2" ) S(r)r dr, R 

(14)

0 1 >(r, r )r dr, S(r )" )   r  

(15)



S(0)"

0

Z(r) dr.

(16)



Here (13) and (14) are straight volume averages of the form (2), owing to the fact that both P(r) and S(r) are de"ned as the average intensity divided by a quantity not dependent on r. Eq. (15) is found by "rst expressing the average intensity at r from (6) in terms of >(r, r ) and then applying the  averaging procedure of Eq. (2). In the resulting integral we recognize Eq. (7). Quite similarly, we may express by (10) the average intensity at r produced by a central source, and then average by Eq. (2). The factor R/r arising from the di!erent denominators in (10) and (11) cancels the r in the integrand, leading to (16). We further note that de"nitions (6) and (8) show that r >(r, r )  . Z(r )"lim   r P

(17)

Many further checks may be thought of. For instance, relation (17) with r and r interchanged, may  be obtained from (6) and(10). And applying the limiting process (17) to Eq. (15), we recover Eq. (16). 3.3. Relations involving the spherical reyection function It was shown earlier [6] that a spherical re#ection function R(k, k ) may be de"ned for  a homogeneous spherical cloud with an arbitrary optical diameter and an arbitrary phase function and albedo. Like all other quantities discussed in Ref. [6], the so-called re#ection includes the zero-order light exiting from the sphere after no scattering at all. Two bimoments of this function

160

H.C. van de Hulst / Journal of Quantitative Spectroscopy & Radiative Transfer 64 (2000) 151}172

appear in the context of the present paper. They are de"ned by

 

 

R(k, k )2k2k dk dk ,        R(k, k )2k dk dk . AH"      The A is identical to the Bond albedo de"ned in Eq. (3). The AH appears in the relation A"

(18) (19)

(20) P(R)"#AH.   The obvious meaning of Eq. (20) is that a shell of isotropic sources placed at the surface emits half of its photons into space, while the other half is re#ected by the cloud according to Eq. (19). No further relations of this kind exist, for all other squares of Fig. 1 involve either sources or detectors inside the cloud. We further recall from Ref. [6] that in the case of isotropic scattering all bimoments, including A and AH, are expressible in products of the moments of the X and > functions of a plane slab with optical thickness 2R. A referee to that paper asked for a clari"cation of the precise relation between two such results appearing in the literature. This question triggered the more systematic review of reciprocity relations reported in the present paper. 3.4. Complementarity relations No matter how the sources are arranged, their energy must end up escaping to the outside world, or being absorbed inside the cloud. This gives on each line of Fig. 1 a relation between the quantities found in column a and column d on that line. Since these relations are based on an inescapable either}or situation we shall call them complementarity relations. The energy absorbed per unit volume per unit time is 4p(1!a)J. Hence, these complementarity relations all involve the factor 1!a. The full set of these relations is 1!A"(4R/3)(1!a)1P2,

(21)

1!P(r)"(1!a)S(r),

(22)

1!1P2"(1!a)1S2.

(23)

Two more lines of Fig. 1 are covered by taking r"0 or R in Eq. (22). Eq. (21) is most easily derived by placing the cloud inside a thermal enclosure with uniform intensity B. The energy originating from the wall and absorbed in the entire cloud per unit time then is (1!A)4ppRB. But we may also calculate the average energy absorbed in each volume element and then integrate, which yields (1!a)4p(4pR/3)B1P2. Equating these expressions we "nd Eq. (21). The factor 4R/3 is seen to arise as the ratio of the volume to the projected area, which may be geometrically interpreted as the average chord length through the cloud. We brie#y digress for an important generalization. Eq. (21) remains valid for a convex cloud of arbitrary shape if the average chord length is de"ned as 4;(volume)/(total surface area). This simple relation has quite a history in the early theories on neutron transport. Case and Zweifel [1, pp. 51}55] refer to it as &Theorem 3'.

H.C. van de Hulst / Journal of Quantitative Spectroscopy & Radiative Transfer 64 (2000) 151}172

161

The other complementarity relations (22) and (23) can be derived in a similar manner by starting from a shell of sources at r or from a homogeneous distribution of primary sources in the cloud. They may have similar roots in history. Conservative scattering (a"1) yields by (21)}(23) the trivial result A"P(r)"1P2"1. However, curious relations are found if we divide by the factor 1!a before taking the limit a"1. We thus obtain






3 dA 1" 4R da S(r)"

,

(24)

?

dP(r) da

,

(25)

?

d1P2 1S2" da




3 dA " 8R da

. (26) ? ? These relations are believed to be new and should, like all others in this section, be valid for arbitrary phase functions. 4. Checks and examples The fun of this paper was in the checking. Some analytical checks are reported in Sections 4.1 and 4.2. In search of opportunities for numerical checks, we found a handful of papers in the past 30 years in which solutions of problems of radiative transfer in a homogeneous sphere had been presented in tables or graphs, regrettably all based on isotropic scattering. More such papers may exist but the wide range of professional journals which may contain such papers makes the search di$cult. By good fortune the set of parameters (R, a) chosen in such examples sometimes coincides. This enables us to make checks as accurate as the tables or graphs permit. Up to four squares of Fig. 1 may thus be connected by simple multipliers in virtue of reciprocity and complementarity. The more dissimilar the problems posed and the methods used in the cited papers, the more gratifying it was to "nd an accurate match. Additional checks are based on numerical integration [Eqs. (13)}(16)] or di!erentiation [Eqs. (24)}(26)]. 4.1. Limiting forms for a"0 All de"nitions of Section 2.2, include the zero-order scattering. This arises from the straight propagation of the light from sources to detectors, with or without a certain extinction. In the limit a"0 only this component remains. We label these limits with the su$x 0 and list them below for all quantities de"ned earlier. The physics is simple, but the spherical shape leads for some functions to expressions involving exponential integrals up to order 3. We write b"2R in Eqs. (30), (33), and (34). > (r, r )"(1/2)[E ("r !r")!E (r #r)],       1!S (r)"P (r)"(1/2r)[E (R!r)!E (R#r)#RE (R!r)!RE (R#r)],      

(27) (28)

162

H.C. van de Hulst / Journal of Quantitative Spectroscopy & Radiative Transfer 64 (2000) 151}172

Fig. 2. Finding the zero-order radiation "eld from a thin shell of sources: (a) at a point inside shell, (b) at a point outside shell.

leading in the two limits to: 1!S (0)"P (0)"exp(!R),   1!S (R)"P (R)"(1/2b)[1#b!exp(!b)],   and upon volume averaging to

(29)

1!1S 2"1P 2"(3/4R)(1!A ).    Furthermore,

(31)

Z (r)"exp(!r),  2 2 A " [1!exp(!b)]! exp(!b),  b b

(30)

(32) (33)

AH"(1/b)[1!exp(!b)]. (34)  Eqs. (27)}(34) may be derived in many ways. All of them follow, in principle, from the mapping theorem in Ref. [6]. Eqs. (28), (33), and (34) are cited directly from that paper. The simpler ones, Eqs. (29) and (32) can be read easily from the de"ning equations in Section 2.2. We sketch in addition a derivation of (27) from scratch. In Fig. 2A we seek the intensity in a point interior to a shell source (r(r ). The zero-order intensity at this point is  I(r, a)"(p/"cos b")exp(!s), where p is a constant and a, b, and s are de"ned in Fig. 2A. If, however, a point outside the shell is chosen with r'r (Fig. 2B), the line of sight from the measuring point may miss the shell, be  tangent to it, or cut it in two places, which contribute to the intensity in that direction with di!erent extinction factors and with di!erent b, but with the same "cos b". The integration of these intensities over all directions, necessary to "nd the numerator of Eq. (6), leads to complicated integrals, which remarkably reduce in both cases to Eq. (27). We are now ready for the crucial question. Do the forms (27)}(34) satisfy (as they should) all interrelations (12)}(17) and(20)}(23), derived for arbitrary a in Sections 3.1}3.4? Performing these

H.C. van de Hulst / Journal of Quantitative Spectroscopy & Radiative Transfer 64 (2000) 151}172

163

substitutions leads in some cases to rather complicated integrals, which however, may be reduced and yield in each case a full identity. This completes the check. 4.2. Expressions in terms of functions dexned for a slab A very extensive literature (see, e.g., Ref. [2]), with both analytical and numerical results, exists on radiative transfer in a plane-parallel slab. Exclusively for isotropic scattering, we may use a classical mapping theorem (Ref. [6, Section 3]) to make these results useful in solving spherical problems. All transfer problems with spherical symmetry in a homogeneous spherical cloud of radius R, diameter 2R"b, with isotropic scattering can be rigorously mapped on equivalent problems for a slab with optical thickness b. Since a detailed discussion falls outside the scope of this paper, we state without proof one gratifying but hardly surprising result. This is the link between the function >(r, r ) for isotropic scattering to Sobolev's [7, p. 99] &resolvent function'  !(q,q) for a slab with thickness 2R: 1 >(r, r )" [!(R!r, R!r )!!(R#r, R!r )].    a

(35)

Since no published tables of !(q, q) are available, Eq. (35) is not directly useful for numerical work. Indirectly it is, for Sobolev also shows that !(q, q) can be expressed by an integral containing only the function of one variable '(q)"!(q, 0)"(a/2)g , (36) \ which is a function elaborately studied by Sobolev [5,7] and tabulated by Nagirner [8] and without the factor a/2 in Ref. [2, Table 17]. At r "R, Eq. (35) becomes  1 (37) >(r, R)" ['(R!r)!'(R#r)]. a The numerical checks, which now follow, all refer to isotropic scattering. We have divided the subject matter into Section 4.3, dealing with functions of r but limited to R"1, and Section 4.4, dealing with functions not dependent on r, with a variety of values of R. 4.3. The internal radiation xeld for R"1 and isotropic scattering Each of the four functions of r only listed in column b of Fig. 1 represents the internal radiation density along a radius of the cloud, from r"0 (center) to r"R (surface). Their accurate de"nitions are in Section 2.3. We have collected numerical data on these functions for the most often chosen practicing example, R"1, i.e., for a cloud with optical path 1 along the radius and 2 along the diameter. This collection is shown in Fig. 3. Most of these data have been copied from the literature, applying appropriate multiplying factors as derived above. Only few have been newly computed or interpolated. This "gure has been presented earlier [9] without a detailed explanation. We drew data from the following papers. P(r) at a"0.9 was taken as o/2 from Siewert and Thomas [10, Table 1]. At a"0.3 we had accurate values only at the end points (Section 4.4) and

164

H.C. van de Hulst / Journal of Quantitative Spectroscopy & Radiative Transfer 64 (2000) 151}172

Fig. 3. Average intensity (in panel 3 multiplied by r) in a spherical cloud with four di!erent assumptions about the primary sources of radiation. Optical depth along the radius is 1 in all panels.

interpolated at intermediate r in the style of Fig. 4, below. At a"0 it was computed from (28). The function 1#S(r) at a"1 is found in the upper curve of Heaslet and Warming [11, Fig. 3], and also in [12, Fig. 1]. The value of S(r) at all other a was derived from P(r) by (22). Z(r) at a"0 was computed from (32). At a"0.3 and 0.9, it is identical to the functions tabulated by Siewert and Grandjean [13, Table 3], also cited in [14, Tables 1 and 2]. No values of Z(r) for a"1 were found in the literature, except at the two ends (Section 4.4). This curve was newly computed by the mapping theorem. Finally, >(r, R) was computed from (35) for a'0. Values of '(q) from Ref. [8] and of g (q) from Table 17 of Ref. [2], were used in this computation. An interpolation was made \ at a"0.3. The value at a"0 followed from (27). All four functions shown in Fig. 3 are interwoven by their values at the end points r"0 and R, and by certain integrals over r. These interrelations, some of which cannot easily be guessed, are documented in Section 4.4 to which we refer for full detail. Let us now look at the physical content of the curves of Fig. 3. Each function has a dual meaning. For clarity we discuss only the meaning where r is the position of the detector (corresponding to column b of Fig. 1). Three of the four graphs show direct plots of the internal average intensity: P(r) as created by external uniform radiation, S(r) as created by uniform internal sources, and >(r, R)/r by sources uniformly covering the surface. Each of these is #at at the center, as expected. It is obvious that S(r) must decrease towards the surface, where extra losses occur. The two other functions decrease towards the center, except for P(r) at a"1, because the radiation has to penetrate from outside. The function Z(r) is the correction factor to be applied to the radiation density of a bare central source, which by itself goes as r\. Two e!ects compete in Z(r). Exponential attenuation, seen in its pure form at a"0, decreases the radiation density toward the cloud's surface; multiple scattering increases it. The net e!ect is that the curves for larger a go through a maximum. This function has been placed in the third panel of Fig. 3 in order to illustrate clearly the equality of the values at r"R in the third panel to those at r"0 in the fourth panel by Eq. (17).

H.C. van de Hulst / Journal of Quantitative Spectroscopy & Radiative Transfer 64 (2000) 151}172

165

The strongly non-linear dependence on a in Fig. 3 arises from the accumulation of many orders of multiple scattering at larger a. The similarity in the ratios of the intervals between the curves, in each panel and at each r, also calls for an explanation. We may crudely reason that the higher orders of scattering form a geometrical series with the ratio g , the "rst eigenvalue for successive  scattering in a sphere, which equals the second eigenvalue for a slab of optical thickness 2R (Ref. [6, Eq. (16)]). Hence any of the plotted quantities may at any r be approximated by f (a)"f #f /(1!g a), (38)    where f and f are constants. Inserting g "0.503 ([2, Table 4],[15]), we "nd that Eq. (38) leads to    the di!erences f (1)!f (0.9)"f (0.3)!f (0)"0.18[ f (1)!f (0)], f (0.9)!f (0.3)"0.64[ f (1)!f (0)]. The di!erences seen in Fig. 3 indeed follow these ratios pretty closely. 4.4. Functions not dependent on r We shall now discuss problems in which the function sought does not depend on r, but on R and a only. This corresponds to seeking the functions in the 16 squares of Fig. 1, which remain after disregarding line 2 and column b. Quantities which are linked by reciprocity or by complementarity will be said to form a subset. The overview read from Fig. 1 is as follows. Subset A: Subset B: Subset C: Subset D: Not treated:

4 4 4 2 2

squares, squares, squares, squares, squares.

Total count:

16 squares,

3 2 2 1

functions: A, 1P2, 1S2. functions: P(R), S(R). functions: P(0), S(0). function: Z(R).

8 functions.

We shall discuss these subsets consecutively. Subset A: The three functions A, 1P2, and 1S2, are linked by (21) and (23), and in the limit a"1 by (24) and (26). Curves of A over the full range of R are in Ref. [6]. Unfortunately, these curves are a little too high around R"5}10 because by an oversight the tangents at R"R were drawn with twice the correct slope. Table 1, column 2, shows the values of A for R"1 computed by Eq. (55) of Ref. [6] using Loskutov's [16] tables of the moments. The tables of Sobouti [17], supplemented for i"3 by Heaslet and Warming [11], have one decimal less. Column 3 of Table 1 shows the A of Siewert and Grandjean [13], who also cite matching results from two other authors published in less accessible places. The agreement is excellent. Columns (4) and (5) of Table 1 show the values of the two other quantities of the subset, derived from A by Eqs. (21) and (23). The divisors 1!a in this procedure lead to a loss of accuracy towards the higher values of a, particularly for 1S2. As an example of what can be done with relatively primitive means, we found the limit of 1S2 for a"1 by extrapolation. In Fig. 4 we plot 1S2!0.47a against a/(2!a). Neither the abscissa chosen in Fig. 4 nor the ordinate are based on exact reasoning, but the factor 2!a in

166

H.C. van de Hulst / Journal of Quantitative Spectroscopy & Radiative Transfer 64 (2000) 151}172 Table 1 The functions in subset A for R"1 a

A vdH [6]

A SG [13]

1P2

1S2

0 0.3 0.5 0.7 0.8 0.9 0.95 0.99 1

0.29700

0.29700 0.42884 0.54226 0.68657

0.5272 0.6120 0.6866 0.7836 0.8438 0.9149 0.9555 0.9907 1

0.4728 0.5544 0.6268 0.7214 0.7811 0.8508 0.890 0.926 0.935

0.54227 0.68670 0.77499 0.87802 0.93630 0.98677 1

0.87801

1

Values adjusted by interpolation in Fig. 4.

Fig. 4. Extrapolation of 1S2 to a " 1. See explanation in text.

the denominator of the abscissa is suggested by Eq. (38) with the (rounded) value g "0.50. Error  bars arising from the fact that A was known in &only' 5 decimals are shown at the larger a-values. The point at a"0.99 was completely useless because the error is multiplied by 10. Deviations from linearity, if present in the curve of Fig. 4, are below the limit of detection. It seems safe to take the open circle as the limit, with the curved brackets as an ample error estimate. This gives 1S2"0.935$0.002 at a"1. This extrapolation procedure is equivalent to making the double numerical di!erentiation in Eq. (26). The values of 1P2 and 1S2 at a"0.99 in Table 1 have been found by interpolation on the straight line of Fig. 4. For a direct check, we also integrated the S(r) at a"1 (see references in Section 4.3) by Eq. (14). This gave 0.939$0.005 (estimated error), within half a per cent of the value derived above. Clearly, such an extrapolation method cannot be recommended for general use. At values of a close to 1, the reverse order is much better. Start with a computation of 1S2, and then convert to 1P2 and to A, whereby several decimals of accuracy are gained.

H.C. van de Hulst / Journal of Quantitative Spectroscopy & Radiative Transfer 64 (2000) 151}172

167

Table 2 Values of 1#aS(R) a"0.8

a"0.9

a"0.95

a"0.99

a"1

Reference

R"0.25

1.0975

1.1117

1.1190

1.1249

R"0.5

1.1893 1.1888

1.2211 1.219

1.2379 1.240

1.2518

1.3538 1.352

1.4297 1.432

1.4726 1.478

1.5097

1.6044 1.602

1.7922 1.795

1.9142

2.0314

1.1264 1.125 1.2554 1.252 1.258 1.5195 1.528 1.521 1.520 1.530 2.0642 2.061 2.065

Exact Ref. [11, Fig. Exact Ref. [11, Fig. Ref. [11, Fig. Exact Ref. [11, Fig. Ref. [11, Fig. Ref. [12, Fig. > integral Exact Ref. [11, Fig. Ref. [12, Fig.

R"1

R"2

3] 1] 3] 1] 3] 1]

1] 1]

Subset B: This subset links the functions AH, S(R), and P(R), by (20) and (22), and in the limit a"1 by (25). We have converted all functions computed by di!erent authors to the function 1#aS(R)"j(R)/j (39)  and present the values in Table 2. This function is the factor by which the emission by the embedded sources in any volume element has to be multiplied in order to "nd the emission#scattering in a volume element just below the surface of the cloud. The choice of function and the choice of the (a, R) grid in Table 2 have been dictated by the available material. The values marked &exact' for a(1 were computed from






a 1 1 aS(R)" 1! d ! d , 1!a 2  b 

(40)

where d "a !b is the di!erence between the ith moments of the X and > function for a slab of G G G thickness b"2R. A convenient derivation of (40) is to take the AH found in Eq. (54) of Ref. [6], and to convert it to S(R) by (20) and (22). An alternative route to (40) is to use the mapping theorem, which leads to embedded sources proportional to q in the equivalent slab problem, for which Sobolev [5, p. 157, 158] has given the complete solution. In the numerical evaluation of (40) I used a 7-decimal table of the moments prepared by Dr. K. Grossman of the NASA Goddard Institute in New York (around 1970) along with the tables used for Ref. [2] but not published. The high accuracy ensures that in spite of the divisor 1!a all four decimals in Table 2 are probably correct. Checks were made against Loskutov's [16] 5-decimal tables (which contain some misprints) and against Ref. [2] (Table 12 with clues on pp. 194}195). The integral relation (15) with r "R gives upon substitution of (35) the simple and apparently  new equation aS(R)"' (b)!(2/b)' (b),  

(41)

168

H.C. van de Hulst / Journal of Quantitative Spectroscopy & Radiative Transfer 64 (2000) 151}172

where b"2R and ' (b) are functions of a and b de"ned in terms of the '(q) of Eq. (36) by G @ ' (b)" '(q)qG dq (42) G  and tabulated in 3}5 decimals by Nagirner [8]. They are expressible in terms of the moments of the X and > functions in various equivalent ways. Converting (41) by means of these expressions, we recover Eq. (40) for a(1. In the limit a"1 Eq. (40) breaks down. Rather than relying on expressions with derivatives with respect to a at a"1, as may in various ways be obtained from (25) of (40), we opted for the use of the equation



3 3 1 S(R)"!1# bd # d # d , 4  2  b 

(43)

from which the exact values in Table 2, column a"1, have been computed. The correctness of Eq. (43) was established in two ways. Firstly, by means of (20) and (22) from an expression for AH at a"1 given (without derivation) by Heaslet and Warming [11]. Secondly, by inserting into the derivation sketched above the Eq. (167) of Ref. [5], valid for a"1. After some reductions [the denominator of (167) equals 2/3], the same Eq. (43) is recovered. The further entries in Table 2 are valued which I read from graphs published in Refs. [11,12]. They provide an excellent consistency check. The uncertainty, estimated at 0.005 throughout, re#ects the plotting and reading error and is no measure of the accuracy of the computation method. A "nal check, based on a crude integration by Eq. (15) of the data on >(r, R) from Section 4.3, and entered as &> integral' in Table 2, also agreed surprisingly well. Series expansions may be derived as readily as numbers. At small R an expansion in orders of scattering is appropriate. We found in Ref. [6] at R"0.25: AH"0.7869#0.1799a#0.0276a#0.0046a#0.0008a#2 .

(44)

fading into a geometrical series with ratio 0.175a. At larger R, particularly in the vicinity of a"1, an expansion in powers of 1!a is more interesting. By numerical di!erentiation I found at R"1 AH"1!1.0389(1!a)#0.923(1!a)!0.90(1!a)#2 ,

(45)

1#aS(R)"1.51946!0.9810(1!a)#0.914(1!a)#2 ,

(46)

in which the last written digits of the coe$cients are uncertain. Subset C: This subset links the functions P(0), S(0), and an integral over Z(r) by (16), (22), and (25). In sampling the literature, I found numbers from at least six papers that must be related. These data, all converted to P(0), and at a"1 to dP/da, by the equations in Section 3, are collected in Table 3. The physical problem solved in the original paper is identi"ed by the code of Fig. 1. Error estimates re#ect the accuracy of my graph reading or numerical integration and not the accuracy of the original numbers. The consistency of the data is clear. The numbers of Sandell and Mattila [18] actually came from Monte Carlo computations on slightly anisotropic scattering (g"0.05), which I did not attempt to correct to g"0. Flannery et al. [19] derive by singular eigenvalue expansion a form that should be valid for all R and all phase functions, but they do not appear in Table 3, since their graphical

H.C. van de Hulst / Journal of Quantitative Spectroscopy & Radiative Transfer 64 (2000) 151}172

169

Table 3 Escape probability P(0) from the center of a cloud with R"1 Reference

Problem solved

Heaslet and Warming [11] Gruschinske and Ueno [12] Sandell and Mattila
[18] Siewert and Maiorino [14] Siewert and Thomas [10] Stark unpub. Kolesov and Sobolev [20,21] van de Hulst unpub. Code [22]

e4 e4 e1 b5 e1#e4 e1 a5 b5 a5

P(0) at a"0.3

P(0) at a"0.9

dP(0)/da at a"1 1.387$0.003 1.385$0.002

0.47$0.01 0.4672$0.0002

0.45

0.865$0.015 0.8748$0.0001 0.87480 0.87485 0.875 0.85

1.38$0.02 1.396 1.39 1.396$0.010 1.7

Identi"ed by the code of Fig. 1.
For g"0.05, i.e., slightly anisotropic scattering.

results are limited to R'5. Their criticism on the accuracy of the Sandell Mattila data is based on a misreading of that paper and should be dismissed. The entries in the next lines of Table 3 were found by integrating numbers from Siewert and co-workers by means of (16) and by a free-hand extrapolation of the ratios [1!P(0)]/(1!a) to a"1. The number on line 6 of Table 3 was computed by Dr. R. Stark using the method of Ref. [6]. Line 7 was read from a manuscript paper by Kolesov [20,21]. Line 8 lists the integral over the Z(r) values shown in Fig. 3. Finally, line 9 is based on an early approximation by Code (Eq. (12) of Ref. [22]). Curiously, his Eq. (9), which he considers a cruder approximation, is closer to the correct value. Subset D: The equality of Z(1) to the limit of >(r, R)/r for r"0, by Eq. (17), has already been illustrated in Fig. 3. For a more precise test, I newly treated the example R"1, a"0.9, and found Z(1)"1.0266$0.0020. This equals 2J,(0) in the notation of Ref. [6]. Siewert and co-workers [13,14] tabulate for the same parameters Z(1)"1.0258 (error estimate not given). The agreement is quite satisfactory.

5. Extensions One set of assumptions, speci"ed in Section 2, underlies all results we have presented thus far. This set was chosen to include a su$ciently wide variety of problems solved by other authors, while still remaining simple. We shall now sketch how these assumptions can be relaxed. 5.1. Obvious extensions E Extension to homogeneous spherical shells (e.g., [23,12,24]). E Extension to inhomogeneous spherical clouds. The theoretical study of stellar and planetary atmospheres has led to a vast literature of radiative transfer in inhomogeneous spheres. Some

170

E E

E E

H.C. van de Hulst / Journal of Quantitative Spectroscopy & Radiative Transfer 64 (2000) 151}172

useful references are: Schmidt Burgk [25], Sobolev [7], Lenoble [26]. The reciprocity principle can be applied but it would seem that the use of complementarity does not lead to simple equations. Extension to clouds of arbitrary shape. No problem, in principle. A classical example was mentioned in our comments to Eq. (21), a new one in the theorem of Section 5.2. Extension to radiation "elds without spherical symmetry. Such "elds arise, for instance, from a source placed excentrically inside or outside the sphere. Presumably, the reciprocity principle can be readily applied to such "elds or to its Fourier components. A much discussed problem in this class is that of a homogeneous spherical cloud illuminated from one side by a star or the sun. Monte Carlo computations by Mattila [27] gave for a variety of phase functions and radii the perpendicularly emerging intensity into various directions. These data, and unpublished ones made available by Mattila, have been analyzed, and a brief account has been published in Ref. [9]. As expected, the numbers for large radii show a smooth transition to those valid for a semi-in"nite atmosphere bounded by a #at surface. The di!erences with the limiting curves are proportional to R\ and the coe$cients of R\ may be determined with fair precision. Exclusively for the asymptotic-behavior of spheres with isotropic scattering, see also [28,29]. The analytical treatment of the sphere illuminated from one side by Shia and Yung [30] yields for isotropic scattering and an arbitrary radius a solution in the form of a series of Legendre functions resembling the Mie solution for scattering by a sphere with constant refractive index. This solution permits a simple transition to the problem of illumination by a homogeneously bright sky (line 1 of Fig. 1 of the present paper). The integration over all directions of incidence makes all Legendre terms vanish, except the term of order 0. As a check, I veri"ed that in that case Eq. (52) of Ref. [30] is, after the proper conversions, identical to Eq. (11) of Ref. [6], which underlies some of the results reported in Section 4. Extension to polarized light. Applying reciprocity requires care. The intricacies may be appreciated from Hovenier [31] and Hovenier and De Haan [32], where this problem is discussed for slab geometry. Extension to intensity, with its full dependence on position and direction. A classical example is the symmetry of the re#ection function in k and k , which holds not only for a slab but also for  the spherical re#ection function. More surprising is the theorem in Section 5.2.

5.2. A theorem about a cloud of arbitrary shape Take a "nite cloud of arbitrary shape, or an arbitrary con"guration of clouds. Let the cloud or con"guration be homogeneous in the sense that the scattering and/or absorbing material at each interior point x is characterized by the same extinction coe$cient i and the same albedo a. The scattering diagram, which for simplicity we assume without polarization, is arbitrary and is permitted to vary from point to point. Let P(x, n) be the intensity at x in the direction n if the cloud or con"guration is exposed to uniform intensity from outside. Likewise, let S(x, n) be the intensity found at x in direction n if primary sources of radiation are spread uniformly over all points interior to the cloud(s). The normalization is so chosen that the intensity produced by the primary sources over an optical path dq is dq. We then have 1!P(x, n)"(1!a)S(x, n).

(47)

H.C. van de Hulst / Journal of Quantitative Spectroscopy & Radiative Transfer 64 (2000) 151}172

171

Proof. Imagine a photon gun at point x sending radiation into a narrow solid angle around direction !n. A fraction p(x, n) of these photons escapes to in"nity. The remaining fraction 1!p(x, n) is necessarily absorbed somewhere in the con"guration. We now require the photon gun to obey Kirchho! 's law so that it can also serve as a directional detector. Then we perform both "ctitious experiments described in the theorem and postulate detailed balancing, exactly as in the simpler examples of Section 3.1. This makes p(x, n) equal to P(x, n) and 1!p(x, n) equal to (1!a)S(x, n), from which Eq. (47) follows. Comments. We can return from (47) in many ways to problems with less general assumptions. For instance, in a homogeneous spherical cloud Eq. (47) becomes 1!P(r, k)"(1!a)S(r, k),

(48)

where r is the optical distance from the center and k the cosine of the angle between n and the outgoing radial direction. Narrowing the assumptions further down to isotropic scattering does not change the form of (48). In this form the essential identity between the solution of the two di!erent problems was noticed earlier by Siewert and Grandjean [13]. Trying to verify their statement by a reasoning based on reciprocity led me to the wider theorem presented above. As a "nal step, we may, of course, average both sides of (48) over all directions and recover Eq. (22). One open question remained. A colleague at Saint Petersburg hinted at the possibility that Kadomtsev [33] might, in fact, have derived this theorem long before. After studying this paper in order to verify if this suggestion was correct, I conclude that Ref. [33] indeed derives certain relations valid for a cloud of arbitrary shape, but does not contain this theorem.

6. Final remarks In response to sceptical remarks about the method used in this paper, I have sometimes made the claims that the method is rigorous not heuristic, accurate not approximate, and general not limited to isotropic scattering. The "rst claim may be a matter of taste but receives some support from the analytical checks. The second claim is well supported by the numerical checks. No data were found as yet in the literature to support the third claim by numerical examples. The &tools' of reciprocity and complementarity cannot be used to solve a problem in radiative transfer. But, once the solution to one such problem has been found, they can be used to "nd simply the solutions to 1}3 other problems. This o!ers opportunities for checks, or for avoiding lengthy additional computations, or for ruling out empirical relations that these tools prove wrong. The question has puzzled me why so few authors (Ref. [10] forms a favorable exception) seize these opportunities. A possible answer is related to the time-honored way of dealing with radiative transfer problems by solving the equation of radiative transfer with a given set of boundary conditions. Sources and detectors play very di!erent roles in that &classical' approach. The sources are de"ned by the boundary conditions. The detectors are added as an afterthought when the author has "nished his calculation and only then asks: &what value (or integral) of the intensity has my particular interest'? From this traditional viewpoint, the diagram of Fig. 1 is very asymmetric. The squares on one line are regarded as &the same' problem; the squares in one column are regarded

172

H.C. van de Hulst / Journal of Quantitative Spectroscopy & Radiative Transfer 64 (2000) 151}172

as &di!erent' problems. Clearly, relaxing the emphasis on this traditional viewpoint is a "rst necessity to bring the possible use of the reciprocity principle into full view.

Acknowledgements This paper has gone through many versions over more than 10 years. Conversations with V.V. Ivanov and V.V. Sobolev, and many others, and the critical remarks by J.W. Hovenier, have been instrumental in giving it a readable form. Dr. E. Nezhinsky helped me with the translation of Ref. [33].

References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] [30] [31] [32] [33]

Case KM, Zweifel PF. Linear transport theory. Reading, MA: Addison Wesley, 1967. Van de Hulst HC. Multiple light scattering; tables, formulas and applications. New York: Academic Press, 1980. Chandrasekhar S. Radiative transfer. London: Oxford University Press, 1950. de Groot SR, Mazur P. Non-equilibrium thermodynamics. Amsterdam: North-Holland Publ., 1969. Sobolev VV. A treatise on radiative transfer. Princeton, NJ: Van Nostrand Reinhold, 1963 (Russian original 1956). Van de Hulst HC. Astron Astrophys 1987;173:115}23. Sobolev VV. Light scattering in planetary atmospheres. Oxford: Pergamon, 1975. Nagirner DI. Astrophysics 1975;9:196}204 (Russian original 1973). Van de Hulst HC. Trudy Astronomicheskoi Observatorii Saint Petersburg University 1993;44:131}53. Siewert CE, Thomas JR. Radiative transfer calculations in spheres and cylinders. JQSRT 1985;34:59}64. Heaslet MA, Warming RF. JQSRT 1965;5:669}82. Gruschinske J, Ueno S. Publ Astron Soc Japan 1970;22:365}71. Siewert CE, Grandjean P. Nucl Sci Engng 1979;70:96}8. Siewert CE, Maiorino JR. JQSRT 1979;22:435}9. Van de Hulst HC, Stark R. Astron Astrophys 1990;235:511}20. Loskutov VM. Astrophysics 1975;9:205}10 (Russian original 1973). Sobonti Y. Astrophys J 1963;7(Suppl.):411}560. Sandell G, Mattila K. Astron Astrophys 1975;42:357}64. Flannery BP, Roberge W, Rybicki GB. Astrophys J 1980;236:598}608. Kolesov AK, Sobolev VV. The light of a spherical nebula under in#uence of a central star (in Russian). Kolesov AK. Private communication, 1990. Code AD. Radiative transfer in circumstellar dust clouds. In: Greenberg JM, van de Hulst HC, editors, Interstellar dust and related topics. Dordrecht: Reidel, 1973, p. 505}15. Bellman RE, Kagiwada HH, Kalaba RE. J. Computat Phys 1966;1:245}56. Wilson SJ, Sen KK. JQSRT 1973;13:255}66. Schmidt Burgk J. Astrophys J 1973;181:865}74. Lenoble J. Radiative transfer in scattering and absorbing atmospheres: standard computational procedures. Hampton VA: A. Deepak Publ. Cy, 1985. Mattila K. Astron Astrophys 1970;9:53}63. Van de Hulst HC. Astron Astrophys 1988;207:182}93. Hovenier JW, van der Mee CVM, de Heer D. Astron Astrophys 1988;207:194}203. Shia RL, Yung YL. Astrophys J 1986;301:554}70. Hovenier JW. J Atmos Sci 1969;26:488}99. Hovenier JW, de Haan JF. Astron Astrophys 1985;146:185}91. Kadomtsev BB. Dokl Akad Nauk SSSR 1957;112:831}4 (in Russian).