JOURNAL OF THERMOPHYSICS AND HEAT TRANSFER Vol. 6, No. 2, April-June 1992 Two-Dimensional Radiative Transfer in a Cylindrical Layered Medium with Reflecting Interfaces N. M. Reguigui* and R. L. Doughertyt Oklahoma State University, Stillwater, Oklahoma 74078 A system of exact linear integral equations for the source function, intensity, and flux is presented for a twodimensional cylindrical medium consisting of up to four layers with reflecting interfaces between the layers. Properties that may change from layer to layer are the single scattering albedo, optical thickness, and refractive index. The incident radiation is coUimated and has a Bessel function distribution. The Bessel function boundary condition reduces the two-dimensional problem to a one-dimensional problem. Superposition is then used to derive the solution for any other boundary condition that is Hankel transformable. A special case of a Gaussian distribution that models a laser beam is presented. Some one-dimensional numerical results are presented for the source function and intensity within one- and two-layer media. Two-dimensional results are presented for back-scattered intensity due to the laser-beam boundary condition. Only the conservative case, optical thicknesses of 2.0 and 5.0, and refractive indices of 1.00 and 1.33 are considered. Zt z j3 f$e AT, 8 8fj O function defined in the Appendix physical depth coordinate Hankel transform parameter extinction coefficient optical thickness of layer i, rzi - rzi_l Dirac delta function Kronecker delta function angle between the incoming and scattered radiation 6 = polar angle Aiy = kernel function Mext = cosine of the polar angle outside the medium fjii = cosine of the polar angle, 6-, within layer i Miy.cr = cosine of critical polar angle between layers / and j IJLO = cosine of incident polar angle of intensity &J = function defined in the Appendix Pij = Fresnel's coefficient of reflectivity for energy crossing the interface between layer i and layer j r0 = overall optical thickness of all layers of a given medium, rzN - TZO Tr = optical radial coordinate, fier Tro = laser beam optical radius, fijr0 rs = optical coordinate along the general threedimensional s direction TZ — optical depth coordinate, pez rzi = optical depth at the bottom of layer i Tzi_l = optical depth at the top of layer / <&ij - function defined in the Appendix </> = azimuthal angle \lft = function defined in the Appendix ft = solid angle (Oi = single scattering albedo in layer / Nomenclature A, B{J DN E Giy g /, /+ It~ I0 IT /, J0 Kfj W HIJ P Qij qzi Ri r r0 5, Scio T Ti tfj x = multiple reflection coefficient, defined in the Appendix = function defined in the Appendix — function defined in the Appendix = component of kernel function for source function integral equation = function defined in the Appendix = Hankel transformed function = intensity of radiation in layer i = intensity in layer i in + TZ direction = intensity in layer i in — rz direction = magnitude of incident intensity outside the medium = intensity transmitted across an interface = magnitude of incident intensity = zeroth order Bessel function of first kind = function defined in the Appendix = total number of layers = ratio of layer i refractive index to that of layer ; = scattering phase function = function defined in the Appendix = z direction flux in layer i = function defined in the Appendix = radial distance from center of medium = laser beam radius = source function in layer i = source function lead term, defined in the Appendix = transmission function, defined by Eq. (32) = function defined in the Appendix = transmission through the interface between layer i and / = dummy integration variable = = = = = = = = Subscripts cr ext i j r 5 Presented as Paper 90-1779 at the AIAA/ASME 5th Joint Thermophysics and Heat Transfer Conference* Seattle, WA, June 18-20, 1990; received Nov. 20, 1990; revision received March 25, 1991; accepted for publication April 2,1991. Copyright © 1990 by the American Institute of Aeronautics and Astronautics, Inc. All rights reserved. * Graduate Student, School of Mechanical and Aerospace Engineering, EN 218. tAssociate Professor, School of Mechanical and Aerospace Engineering, EN 218. Member AIAA. related to critical angle external to all layers of the medium layer i layer; related to radial direction related to a general three-dimensional direction related to depth direction Hankel transformed variable Superscripts 1,2 = either function 1 or 2 (not both at once) + , - = either positive ( + ) or negative ( - ) z direction (not both at once) 232 233 REGUIGUI AND DOUGHERTY: TWO-DIMENSIONAL RADIATIVE TRANSFER Introduction I N radiative transfer problems, many realistic situations, such as solar pond modeling, laser-material interaction, and underwater radiation back scattering, require a multidimensional analysis that accounts for the effects of the refractive index. However, most of the results found in the literature employ many simplifying assumptions, such as one-dimensional geometries or multidimensional geometries that neglect the effects of the index of refraction and consider only one or two aspects of the radiative transfer problem of interest: 1 reflecting boundaries,1"4 diffusely incident radiation and col4 5 limated incident radiation, isotropic scattering and anisotropic scattering,6 variation in absorption and scattering with depth,7-8 and two or three dimensionality.5-6-8 A more realistic analysis has to include the effects of multidimensionality as well as the effects of possible variations in the pertinent physical properties, such as index of refraction, and absorption and scattering coefficients of the medium. The index of refraction change across an interface, coupled with the directional redistribution of radiation by anisotropic scattering, has been shown to alter considerably the radiative field within the medium.1 Cengel and Ozisik2 presented onedimensional data for a one-layer slab which indicate that the constant reflectivity approximation often leads to severalfold under- or overestimated results as compared to the exact handling of the reflection coefficient; thus disproving the assumption of using average values. Most of the work that includes the effect of index of refraction has been limited to one dimension. Buckius and Tseng1 considered a plane-parallel isothermally emitting medium that scatters anisotropically. Armaly and Lam3 and Dougherty4 have investigated the influence of refractive index on the reflectance from a semi-infinite one-dimensional absorbingscattering medium. Although a lot of work has been done for the two-dimensional geometry, index of refraction effects have been neglected for much of that work. A two-dimensional problem that has received attention in the past 10 years is the back scattering of a laser beam by a multiple scattering medium when the incident beam is normal to the surface of the medium.5'6 Most of these studies were based on the assumption that the index of refraction is unity. Crosbie and Dougherty presented graphical and tabular results for isotropic5 and anisotropic6 scattering in a finite cylindrical medium exposed to a laser beam. This present work builds on these studies to include the effects of the index of refraction and the effects of layering the medium. Some investigators9-10 have considered superposing the unit refractive index results or introducing a factor in the results based on the unit index of refraction to account for index of refraction effects. Crosbie and Dougherty9 have proposed superposing the unit refractive index results to obtain directional and hemispherical reflectances from a two-dimensional cylindrical medium. No numerical results were presented, but some asymptotic solutions were examined. In a comparison between theoretical and experimental results for back scattering from a two-dimensional optically thick medium exposed to a laser beam, Nelson et al.10 have found that agreement between theory and experiment is improved when the theoretical back-scattered intensity for unit refractive index is reduced by (1 - pN)2ln2. However, this approximation was found to be valid only in the thin limit. Studies of radiative transfer in composite (layered) media including the effects of scattering are very limited and mostly deal with one-dimensional situations. Shouman and Ozisik11 considered the problem of radiative transfer in an absorbing, emitting, and isotropically scattering two-layer slab with diffusely and specularly reflecting boundaries. The slab was irradiated externally by diffuse radiation at the boundary surface. The FN method was used to compute results for the transmissivity and reflectivity of the slab. The results showed that by increasing the single scattering albedo of the second slab, the reflectivity and the transmissivity of the composite slab increase. The transmissivity was found also to be influenced by the relative optical thicknesses of the slabs. A few other studies also exist for media with layered properties; but usually the index of refraction is assumed to be unity, which considerably simplifies the analysis. Stamnes and Conklin7 have developed a matrix method to solve the discrete ordinate approximation to the radiative equation in a homogeneous, plane parallel atmosphere consisting of N adjacent layers. However, because of the one-dimensionality of the analysis and the unit index of refraction, this study is applicable to the field of atmospheric physics and meteorologyA general multilayer study was done by Sutton and Kamath.8 They considered a three-dimensional rectangular medium with layered properties and unit index of refraction. The boundaries between the layers were considered transparent, hence neglecting the effects of index of refraction changes. The volume-averaged discrete ordinates method was employed to solve for the transmitted and reflected radiative flux. In the present investigation, the exact solution for the radiative transport equation within an absorbing and isotropically scattering two-dimensional cylindrical medium is presented. The medium is finite in height and infinite in the radial direction, has internally reflecting top and bottom surfaces, and consists of layers whose single scattering albedo, optical thickness, and index of refraction may differ (see Fig. 1). The interfaces are assumed to be smooth so that Fresnel's equation and Snell's law can be used. The incident radiation is collimated and normal to the top surface of the medium and it is assumed to be azimuthally symmetric. It has a radial variation in the form of a Bessel function. The Bessel function boundary condition is used to separate variables, and therefore reduces the two-dimensional equations to one-dimensional form. Superposition is then used to solve the problem in which the incident radiation is in the form of a Gaussian function (i.e., a laser beam). No radiation is incident on the bottom surface of the medium. Development The transfer equation in a given coordinate system for a certain layer i in a medium consisting of many layers and in which scattering and absorption take place is5 = /,(T.)P(cps 6) dft (1) Io Tr-XX> CDI Fig. 1 Geometry of a layered medium, infinite in the rr direction and finite in the TZ direction. REGUIGUI AND DOUGHERTY: TWO-DIMENSIONAL RADIATIVE TRANSFER 234 In the source function 'S-h a>, is assumed to be constant for layer -i, but could differ among the layers. For the remainder of the analysis, P will be assumed equal to one for isotropic scattering. Decomposing Eq. (1) into its positive (7+) and negative (/") directions, and then solving both equations using an integrating factor approach, yields5 x exp[(T2l_! - T2)//xf] + exp[(r; Jrzi-l the intensities at the interior interfaces between the different layers are unknown, and they need to be eliminated from the set of Eq. (6) in order to solve for the source functions for all layers. Intensities and Source Functions We will consider collimated incident radiation at the topmost boundary. Thus, the radiation crossing the interface of an arbitrary layer i has, in general, a collimated part, due to the collimated incident radiation, and a diffuse part due to scattering from the adjacent layers. The collimated radiation can be represented using the delta function, i.e. (2) (7) and (3) /; is the magnitude of the intensity and (/t/,0) determines the direction of the collimated radiation. From the law of conservation of energy, the relationship between the incident and the refracted intensities is3 (8) where r'r is given by (T;)2 = (4) with «2 = Tl = (r'z - r2)2(l - M ; 2 (5) and r'nf . = r, for t2 = r2l (see Fig. 2). All of the quantities with subscript i are related to layer i, and rzi and Tzi_l refer to the bottom arid top of the layer i, respectively (see Fig. 1). 5 is a function of rr and rz in Eqs. (2) and (3) only because of the isotropy of the problem and the axisymmetric loading assumption. ^(r^r^.^p.^) and 7f (T^T,,,/^) are the boundary conditions at the top and bottom of layer /, repectively, which is finite in the depth (rz) direction and infinite in the radial (TV) direction. Substituting Eqs. (2) and (3) into the right-hand side of Eq. (1) and rearranging the result gives the integral equation for the source function in layer / My = [1 - (1 - (9) The imaginary part of the refractive index is assumed to be small relative to the real part, which is true for most dielectrics. When n^ > 1, Eq. (9) cannot be satisfied for ^ < ju,/7cr where /4cr (10) Therefore, when niy > 1, all intensity that is coming at angles Of greater than 0£> cr (or ^ < ^ cr) will be totally reflected (# = 1.0). When the intensity crossing the interface between layers i and j is diffuse, the amount of the refracted intensity in medium j is3 (11) x exp[(Tz|._, - Tz)/p.'t] + /r(T;.,TZI.,/i>') x exp[-(r2/ - T2)/M/]} dM; d<£' + ^ P* P P* 477 JO where 7y is the magnitude of the intensity refracted into layer ;'. The reflectivity function p,7 is given by Fresnel's equation1 and /Xy is given by (SnelPs Law1) JO JTZ/-I (6) Equation (6) actually represents a set of equations to be solved for all layers in the medium. Besides the top and bottom boundary conditions of the medium which should be known, In Eq. (2) /^(T^_i,T z/i _ !,/*,-,<£) is the sum of the incident radiation that crosses the interface and the radiation originating in the layer that is reflected at the top interface in the positive direction. Therefore (12) Similarly, the boundary condition /"(T^T.,,,/^,^) in Eq. (3) is the sum of the interface transmittance of radiation from layer / + 1 plus the interface back reflection of scattered radiation at the bottom of layer /. Therefore (13) Fig. 2 Cylindrical coordinate system. In Eq. (12), /r( T ri-u T z/-i?M'/^) is. replaced by Eq. (3), where TZ is set equal to r2l_!. The transmitted intensity Jr( T r/-i> T z/-i>A t /><£) could be collimated and/or diffuse. Similarly, in Eq. (13), /^(r^r^,^,^) is written from Eq. (2) where TZ is set equal to r2l and /f (r^,T2/,^/,<^) is determined from Eq. (8) and/or Eq. (11). Both resulting equations can be solved simultaneously for /^(r^.^T^-.!,^,^) and Irfr'rhTziMd) [using Eqs. (2) and (3) a second time in the process]. This will result in a system of linear simultaneous REGUIGUI AND DOUGHERTY: TWO-DIMENSIONAL RADIATIVE TRANSFER equations in terms of the interface intensities. When the solution to this sytem is put back into Eqs. (2) and (3), we get 235 Equation (16) represents a system of simultaneous integral equations in terms of the source functions of each layer for a cylindrical medium consisting of N layers. 5f0(/^l,r2) is an initial term due only to the collimated boundary condition, and it is being evaluated at //,,- = 1 because of the normal incidence of the radiation [see Eq. (Al)]. Again, //,, should be greater than the critical angle /x,/y cr when considering radiation from layers other than layer i. The auxiliary function T/y is another collection of terms due to the transmission, reflection, and attenuation of the radiation as it crosses the different layers, and it depends on the relation between / and ; as given in the Appendix [Eqs. (A16)-(A18)]. Note that Eqs. (14), (15), and (16) reduce to those of Buckius and Tseng1 when there is only one layer (except that their boundary condition is diffuse while this boundary condition is collimated). Radiative Flux Once the intensity of radiation is determined, the radiative flux in the z direction may be computed from substituting Eqs. (14) and (15) into the following equation: P" Stt.TiWO^T;) dr' JTzi-1 (14) and ~i (i?) Having the basic two-dimensional equations derived for the collimated boundary condition, it is easy to write those equations for other forms of boundary conditions. This is explained in more detail in the next section. Special Boundary Conditions /=! Bessel Function Boundary Condition exp[(r, S,(T;,T;)*KM»,T;) (15) where because of the axisymmetric loading, the intensity does not depend on </> directly, except through rr. The equations for AX/0, DM, Z+GO, Zr(/0, *,W Rft*,), 77fo), T?(/ij), BJO^Tv), Bl{ni9Tr), *ft(Hi,Tz)9 and ttf(Hi,Tz) are given in the Appendix. The first lead term in each of Eqs. (14) and (15) has an infinite value when IJLO = 1 due to the Dirac delta. The remaining terms in the equations include the radiation contributions from all the layers other than layer /, as well as the contribution from within the layer. Notice that when the number of layers is greater than one, r'r will be given by a more general formula than Eq. (4), i.e., Eq. (A28). In the process of obtaining Eqs. (14) and (15), the following equation for the source function was obtained: Even for the isotropic case, the two-dimensional equations are still difficult to solve. However, if the incident radiation is assumed to be in the form of a Bessel function, then the source function, the intensity, and the flux can be shown5 to be separable functions of rr and r2, and the two-dimensional governing equations can be reduced to simpler one-dimensional forms. The Bessel solution can then be superposed to obtain the solution for any Hankel transformable boundary condition.5 A Bessel function boundary condition has the following form: (18) where /, is a constant, the magnitude of the intensity. Source Function The Bessel-varying boundary condition suggests that, for this case, the source function can be separated,5 i.e. (19) and S0,XTz) represents the j3-varying one-dimensional source function for layer /. Substituting Eqs. (18) and (19) into Eq. (16) and then, using a technique similar to Crosbie and Dougherty,6 the /3-varying source function can be found12'13 (20) REGUIGUI AND DOUGHERTY: TWO-DIMENSIONAL RADIATIVE TRANSFER 236 where a), - r'z) + f /.(flr •rpftcr Jrzi-l (21) (26) -(T; - Tz)/M,]5p,(r;) dr'z where raij is defined by Eq. (A29), and where it was necessary to split the integral over /*/, because when/ is different from i, /Lty will be restricted to the interval between /*l>fCr [given by Eq. (10)] and 1.0. E in Eq. (21) is a Bessel-varying exponential integral, and it is defined as Flux Similarly, the flux can be written as a separable function of Tr and TZ, i.e. (22) (27) which is related to the generalized exponential integral of Ref . 5. This integral reduces to the exponential integral when ft is zero (the one-dimensional case). I\ can be numerically computed for any given values of /i,', TZ, and r'z from Eqs. (A16)(A18). Then Eq. (21) can be integrated numerically to compute A/y, and this can be used in Eq. (20) to solve for the /3varying source function. Substituting this equation and Eqs. (23), (25), and (26) into Eq. (17), and following a similar development as for the source function gives wto = [exp[(T,,_, - T,)]Z,+( Intensity Equation The Bessel-varying boundary condition suggests that the intensity can be separated,5 i.e. (23) ,,.,, - T, - Any/*,'] - where /^,/(TZ,^,<^) represents the /3-varying one-dimensionall intensity for layer i. Similarly, the effective intensity terms B and B2 [of Eqs. (14) and (15), defined in Eq. (A19)] can be separated as ] i-l + 2 /=! (24) } where B^j(fit) is now defined in terms of the integral of Spj(rz), as in Eq.' (A20). Substituting Eqs. (18), (23), and (24) into Eqs. (14) and (15) gives the j8-varying intensity equation, viz. (x) f 1 [Tzi d/Lt/ + Y Jo JTZ i sgn(rz - T'Z) .n x exp(- An/ft) 2 ttflj/M . AM) drz + where Bqft is a modification of the function (A20), replacing S^r'z) by S,,(r'z)J0(ftra^ - T^J^KM.) f" S,XT;)W(MI,T ft9 given in Eq. Other Boundary Conditions JTZ/-I (25) and (28) After the basic solution for the Bessel function boundary condition has been found, this can be used in superposition to find the solution for any other boundary 5condition that is expressible in terms of a Hankel transform. Representation of the rr-varying boundary condition, /0(rr), by a Hankel transform yields5 exp[ - (rr/ I0(rr) (29) - 1) + The solution of any other problem with a boundary condition that is Hankel transformable can be found from an inversion process similar to Eq. (29)—substituting the source function or flux for intensity as required. 237 REGUIGUI AND DOUGHERTY: TWO-DIMENSIONAL RADIATIVE TRANSFER Gaussian Distribution In this work, a Gaussian functional form is used to model a laser beam,5 such that I0(rr) and is /.(TV) = (30) where rro is the effective optical radius of the beam (at the e~* point). The source function inversion equation becomes12 (31) where the integration over /3 has been replaced by integration over x ( = /3Tr), and rjrr has been replaced by rjr (assuming uniform properties within each layer). The corresponding intensity and flux equations are similar to Eq. (31).12-13 sults. Tables 1 and 2 present a comparison of results generated by this multilayer program with those of Dougherty.4 In those tables, "Direct" and "Ambar" refer to solution methods employed by Dougherty, solving the integral equation for the source function by discretizing rz (Direct), and solving only for the source function at the boundary by a variation of Ambartzumian's method (Ambar). Table 3 gives a comparison between the current results and those of Dougherty14 for a one-layer medium with both boundaries having reflective interfaces. Current 1L and Current 2L refer to employing the modeling of this paper, and using either one uniform layer (1L) of a given number of rz quadrature or two uniform layers (2L) having the same properties, but twice the number of rz quadrature. In addition, Table 3 uses T(p), which is given by12 (32) Results Numerical Procedure The method of solution consists basically of using Gaussian quadrature to do the different integrals involved, and applying successive approximation to determine the source function of Eq. (20). First, the different A/y(j8,rz,rO values given by Eq. (21) are computed at the quadrature points TZ and T'Z for every P. These values necessitate a large amount of storage, about one-third megabyte per )3 for a one-layer problem. In addition to storage, the computations require approximately 45 CPU minutes on an IBM 3090 per /3 calculation for a one-layer problem. Increasing the number of layers obviously increases this computational time. But once the Al7 values are computed, the iterations to compute 5^ become quite fast. Also, since the A/y values do not depend on albedo (w), the source function can be determined for a variety of a) values without having to evaluate the Aiy values for every (o. To start the iteration on Eq. (20), S is first set equal to 5f0(l,rz) of Eq. (Al). Test runs have shown that the function riy(//>,Tz,T^)vhas several irregularities in the curve as a function of the cosine of the angle p. This is due to the sharp change in the reflectivity coefficient function1 whenever the argument (a function of ju,) becomes larger than the appropriate critical angle. To alleviate this problem, the integral was divided into several subintervals determined by all of the critical angles that can exist when going from layer / to layer;; and a dense quadrature was placed around each of these breakpoints. This procedure increased the execution time considerably, but on the other hand, helped to assure the accuracy of the results. Convergence in computing the source function was achievable for any accuracy desired and for any set of variables (w,n,AT, . . .). However, for albedo equal to unity, the convergence was relatively slower than for other values. To test the accuracy of the computer program, some published results were first reproduced. Then, some particular two-dimensional multilayer problems were considered. Most of the cases from the literature that have been considered are in good agreement with the generated results. Unfortunately, all of those cases deal with one-dimensional geometries, as there is no other work in the open literature that considers index of refraction effects for two-dimensional geometries. One-Dimensional Results Two main problems were considered for comparison. The first one considers a single-layer one-dimensional medium with a top reflecting surface.4 The second case has, in addition to the first case, a bottom reflecting surface.14 Both cases have collimated intensity incident on the top surface and no incident intensity on the bottom surface of the medium. By setting )8 equal to zero and the number of layer (N) to one, the current results were comparable to one-dimensional single-layer re- In all of these tables, the maximum error between any two entries is less than 0.07%. If the layer is divided into two sublayers that have the same properties, then the two-layer model (Current 2L) predicts the same values as compared to those from the one-layer model (Current 1L) as shown in the tables. The slight differences in results between the 1L cases and the 2L cases are due to the increased number of quadrature for rz. Some additional one-dimensional results will be presented first in order to show the effects of changing the index of refraction between two layers. The one-dimensional results are a special case (/3 = 0) for the /3-varying solution. Figure 3 presents a plot of the source function as a function of the medium depth (rz). The medium consists of two layers. The top layer models a thin film of paraffin oil15 (n^ = 1.48) with an optical thickness of 0.1. The second layer has an optical thickness of 5.0 and an index of refraction of 1.33 (water). The albedo is assumed to be unity for the water and 0.9 for the oil. These albedos can be experimentally simulated by adding small latex spheres to the fluids.15 The graph shows a Table 1 Comparison of current results with those of Dougherty4 4077^/0 = o,i (0,/i)/ s~ _fnv Direct4 Ambar4 Current 1L Current 2L 'm(l) 1.486318 1.486321 1.485868 1.486116 4)i(l) mf-i n fji = 0.5 /Lt = 1.0 ~ <?z0 = 0,l(0)] 2.72818 2.7282 2.72752 2.72772 4.25353 4.25366 4.25246 4.25710 0.60330 0.60331 0.60315 0.60323 «10 = 1.33, «12 = 1.00, TO = 6.0, to, = 0.5, 0 = 0.0. Table 2 Comparison of current results with those of Dougherty4 Direct4 Ambar4 Current 1L Current 2L 'oi(l) /LI = 0.5 /Lt = 1.0 2.001094 2.001031 1.999625 2.00037 5.44812 5.44790 5.44448 5.44628 8.74655 8.74616 8.74038 8.74066 - ^=0,i(0)] 0.92609 0.92604 0.92544 0.92577 n10 = 1.50, «12 = 1.00, T0 = 8.0, o>! = 0.7, /3 = 0.0. Table 3 Comparison of current results with those of Dougherty14 'oi(l) Direct14 Current 1L 3.857896 3.80523 , 'm(l) 2.557394 2.556909 T(fJL 2.101162 2.100681 /i,o = 1.50, nu = 1.20, TO = 1.0, o>, = 1.0, 0 = 0.0. 1.965850 1.965385 REGUIGUI AND DOUGHERTY: TWO-DIMENSIONAL RADIATIVE TRANSFER 238 0.40 6.0 1 Oil/Water Layers 2 Water Layer, To = 5.0 3 Water Layer, ro = 2.0 5.0 (0.0, 1.0), TO= 2.0 (To, 1.0), TO • 2.0 (0.0,1.0), r0 : ' 5.0 (To, 1.0), TO = 5.0 0-0 4.0 CO . 3.0 2.0 0.0 1.0 2.0 3.0 4.0 5.0 0.00 10-3 1Q-2 1Q-1 1.0 101 102 1Q3 ft OPTICAL DEPTH (Tz) Fig. 3 One-dimensional source function inside the medium (for oil: at = 0.9 and n = 1.48; for water: o> = 1.0 and n = 1.33). Fig. 6 Effect of optical thickness on the 0-varying intensity for o> = 1.0 and n = 1.33. dramatic drop in the source function just at the interface between the two layers. This drop is due to the change in 0.225 index of refraction from 1.48 to 1.33, which results in more radiation being reflected back into the thin layer. The thin layer exerts less influence on the source function in the medium as we move away from the interface. This is shown by comparing curve 2—the source function for the same layer of water (r0 = 5.0) but without the second thin layer—to the i? 0.150 curve resulting from the addition of the thin layer. The third curve (TO = 2.0, water only) on Fig. 3 is included to show the effects of varying the optical thickness of the layer on the source function. Small optical thickness means less scattering, and hence a smaller source function. 1 I , water-oil — r. 0.075 Another special case of interest is the optical thickness and q 2 I*, water-oil the index of refraction effects on the angular distribution of o 3 IT TO = 2.0, water _ the back-scattered intensity (7~) and the transmitted intensity 4 I* TO = 2.0, water " i'5 (/+). The one-dimensional solution, which is still the special 5 n TO = 5.0, water " case of j8 = 0, gives satisfactory insight into these effects. 6 I* TO = 5..0, water 0.000 Figure 4 gives the back-scattered intensity (7~) and the trans0.00 0.20 0.40 0.60 0.80 1.00 mitted intensity (7+) outside the medium for the cases of oneand two-layer media. For the one-layer medium (n = 1.33, CD = 1.0), comparing r0 = 2 and r0 = 5, the back-scattered Fig. 4 The intensity outside a one-dimensional medium: one- and intensity is increased for the larger TO. This trend is reversed two-layer effects (for water only: w = 1.0; for water-oil: o> = 1.0 for the transmitted intensity. When a thin layer (r0 = 0.1 and and TO = 5.0 for water, while o> = 0.9 and r0 = 0.1 for oil). (o = 0.9) of paraffin oil is introduced on top of the water layer (for the case where TO = 5), both transmitted and backscattered intensities are reduced drastically. n = 1.33, TO = 2.0 2 80,1(0.0), n - 1.33, T0 = 5.0 3 S^j (T0), n = 1.33, TO = 2.0 4 S£,J(TO), n = 1.33, TO - 5.0 5 S^j(O.O), n = 1.00, TO = 5.0 6 S^j(T0), n - 1.00, r0 - 5.0. O.OI___i_ 10-3 1Q-2 io~1 _____________ 1 2 1.0 10 10 103 104 Fig. 5 Effect of refractive index and optical thickness on the /3-varying source function for o> = 1.0. (Curves 5 and 6 are from Ref. 16). Two-Dimensional Results Before presenting some two-dimensional results, we will examine the /3-varying solution. Figure 5 presents the influence of optical thickness TO, spatial frequency /3, and index of refraction on the source function. The effects of these parameters on back-scattered intensity and transmitted intensity are illustrated in Figs. 6 and 7. The source function at the top and bottom of the medium exhibits similar trends for both optical thicknesses. For small /3 values, the source function at the top is larger for the smaller r0. This trend reverses for the bottom boundary source function. These functions are compared to the results by Crosbie and Koewing16 for unit refractive index. It is found from the comparison that the index of refraction increases the /3-varying results for the source function, but reduces the intensity results (Figs. 5 and 7). Note that, regardless of optical thickness and refractive index, the /3-varying source function results approach nearly the same asymptote as /3 increases. As ft becomes very large, the results represent small Tr behavior (see Refs. 5 and 6). After obtaining the IB-varying solution, numerical integration [e.g., Eqs. (31)] gives the rr-varying solution. Figures 8 239 REGUIGUI AND DOUGHERTY: TWO-DIMENSIONAL RADIATIVE TRANSFER 3.00 2.0, 2.0, 5.0, 5.0, 5.0, 2.0, n n n n n n = 1.33 = 1.33 = 1.33 = 1.33 = 1.0d = 1.00 2 0.60 0.00 . 10-3 1Q-2 1Q-1 1.0 101 10)2 103 Fig. 7 Effects of refractive index and optical thickness on the ftvarying intensity for o> = 1.0. (Curves 5 and 6 are from Ref. 16.) 10-3 1 2 3 4 10-4 n n n n = 1.00, T0 = 2.0 = 1.33, T0 = 2.0 = 1.00, r0 = $.0 = 1.33, T0 = 5.0 X 10-5 o b o 1(T ^ and 9 are plots of the reflected intensity outside the medium. The first plot (Fig. 8) provides a close look at15the trends over the normal experimental measurement range. Figure 9 gives the intensity over a wider range, specifically within the original incident laser beam. The trends clearly show a dramatic increase in the reflected intensity within the original laser beam—which is to be expected. Curves 2 and 4 are for the reflected intensity from a one-layer medium with top and bottom reflecting boundaries (ri = 1.33). These curves show that the back reflected intensity increases in value as the optical thickness of the medium increases. This is due to the increased amount of scattering [events]* Some results for unit refractive index16 are also included in Fig. 9 for comparison. It is apparent from the plots that assuming unit refractive index in the case of a water medium would lead to severe overestimatioh of the reflected intensity. However, the general trends of the two sets of solutions are quite similar. Conclusions The theoretical work presented herein will provide a more flexible and general tool to study a wider range in indices of refraction, and to investigate their effect on the back-scattered intensity, than is available in the existing models. The computerized model has the capability of varying the refractive index, along with handling internal reflection at the boundaries and allowing the albedo to vary through layers of variable optical thickness of the medium being analyzed. With such a model, validated by experimental work (now being conducted at Oklahoma State University15), radiative transfer within simulated solar ponds or any other multilayered medium can be studied for a much wider variety of parameters than those assessed by the hardware. Appendix: Auxiliary Functions In this appendix, several functions that were referenced in the body of this paper are written explicitly. In the following equations, N represents the number of layers, and i represents a particular layer: 7 & 10~ 10,-8 ^ 10-9 100.0 10.0 1000.0 (Al) 10000.0 (A2) r/r0 Fig. 8 Effect of refractive index and optical thickness on reflected intensity far from the incident beam for <D = 1.0. O-Varying results of Ref. 16 were used to obtain curves 1 and 3.) ) = [1 - p,,..^^ (A3) (A4 is defined in the following equations: 1 n = 1.00, T0 = 2.0 2 n = 1.33, T0 = 2.0 3 n = 1.00, TO = 5.0 i) = n[exp[-AV/i,]], for abs(/ - i) > 1 (A5) 4 n = 1.33, TO = 5.0 - where v - minimum (/,;') + 1, and w = maximum (/,;) 1, and £,(/*./) = 1, for abs(/ — 1)^1. The <l>'s are defined as follows: ], i >J (A6) *=7 where M, = [1 - (1 - M?>4]"2. I f / = A then OJG*,) = 1 10 10-10 0.01 0.1 1.0 10.0 ...._______ 100.0 1000.0 10000.0 r/r0 Fig. 9 Effect of refractive index and optical thickness on reflected intensity for <w = 1.0. (0-Varying results of Ref. 16 were used to obtain curves 1 and 3.) *S(M») = fl I4JM*K.*-ifaj],i</ (A7) *-!+! If i = /, the Q^GI/) = 1- The definitions for $ are (A8) REGUIGUI AND DOUGHERTY: TWO-DIMENSIONAL RADIATIVE TRANSFER 240 where A: = (i 4- ;)/2, m = (3; - 0/2, / = (3i - y)/2. And for abs(/ - /) = 3, Q and # are given as (A9) Z is defined as X tkj(fJLk)Gjm(lJLj) (A10) X x exp(-2Ar//x/)[p/n(/A/) + (All) x C (A25) + , and T, are given by x (A12) x (A26) (A13) for (i + 2) < f (A14) where A: = (2; + i)/3, / - (/ + 20/3, m = (4; - 0/3, and r = (4/ - y)/3. For all of the above definitions of Qnj and Kv, the IJL values are related as follows: JJL£ = [1 - (1 jit?)n?e]1/2, where e = y, k, /, m, n, or r. In addition, the G,y values are If i > N - 2, then 77(/O = 1 *?(/**) = 1 - t2 Ci*(fO, ^r i > 2 =2 Gfa) = ^Mi)p«OA/X^-XM/)exp(-2AT//il) (A15) If i < 2, then 77(ju,,) = 1. For F, when j < / if M/ < ^,2/-y,cr, then Gifa) = 0. Lastly, for cases where the number of layers is greater than one, the equation for f'r [Eq. (4)] should be modified to include passage through all intermediate layers as follows: = (rr)2 + (raiy)2 - (A16) When j > ; (A28) where max(/,y) (A29) k = min(/,y) and when / = ;', then = [!-(!- /Lt/2)n?J1/2. Note that when k = min (/,/), zk-i — Tz> an<^ when k = max (/,/), r2A. = T^. In the special case where i = ;, Eq. (A29) becomes Eq. (5). r,,(M,',T2,T2) = RHjp T x exp[-(T2,_, - (A18) x exp[-(Tz( - rz For the B values r (A27) S/T;,T2)*;-2(M/,T2) dr2//ty. (A19) ./TZ/-1 Acknowledgments This work was supported in part by the University Center for Energy Research at Oklahoma State University, Grants 1150714 and 1150718, and by National Science Foundation Grant CTS 8907149. In addition, the authors are grateful to IBM and its Palo Alto Scientific Center for a grant of supercomputing time through the Research Support Program. References 2 f" sw(T;)v;- Ovr;) dry/*, JT,/-I (A20) The two functions Q and ^ are given only for four layers. For abs(/ - 1) = 1, Q and K are given as (A21) (A22) For abs(/ - /) = 2, 0 and /C are given as (A23) (A24) ^uckius, R. O., and Tseng, M. M., "Radiation Heat Transfer in a Planar Medium with Anisotropic Scattering and Directional Boundaries," Journal of Quantitative Spectroscopy & Radiative Transfer, Vol. 20, 1978, pp. 385-402. 2 Cengel, Y. A., and Ozisik, M. N., "Radiation Transfer in an Anisotropically Scattering Slab with Directional Dependent Reflectivities," American Society of Mechanical Engineers Paper 86-HT28, 1986. 3 Armaly, B. F., and Lam, T. T., "Influence of Refractive Index on Reflectance from a Semi-Infinite Absorbing-Scattering Medium with Collimated Incident Radiation," International Journal of Heat and Mass Transfer, Vol. 18, 1975, pp. 893-900. 4 Dougherty, R. L., "Numerical Results for Radiative Transfer in a Semi-Infinite Absorbing/Scattering Slab Exhibiting Fresnel Reflection," AIAA Paper 88-0076, AIAA 26th Aerospace Sciences Meeting, Reno, NV, 1988. 5 Crosbie, A. L., and Dougherty, R. L., "Two-Dimensional Isotropic Scattering in a Finite Thick Cylindrical Medium Exposed to a Laser Beam," Journal of Quantitative Spectroscopy & Radiative Transfer, Vol. 27, No. 2, 1982, pp. 149-183. 6 Crosbie, A. L., and Dougherty, R, L., "Two-Dimensional Ra- REGUIGUI AND DOUGHERTY; TWO-DIMENSIONAL RADIATIVE TRANSFER diative Transfer in a Cylindrical Geometry with Anisotrppic Scattering," Journal of Quantitative Spectroscopy & Radiative Transfer, Vol. 25, 1980, pp. 551-569. 7 Stamnes, K., and Conklin, P., "A New Multi Layer Discrete Ordinate Approach to Radiative Transfer in Vertically Inhomogeneous Atmospheres," Journal of Quantitative Spectroscopy & Radiative Transfer, Vol. 31, 1984, pp. 273-282. "Suttori, W. H., and Kamath, R., "Participating Radiative Heat Transfer in a Three-Dimensional Rectangular Medium with Layered Properties," American Society of Mechanical Engineers Paper 86HT-25, 1986. 9 Crosbie, A. L., and Dougherty, R. L., "Influence of Refractive Index on the Two-Dimensional Back-Scattering of a Laser Beam: Asymptotic Solutions," Journal of Quantitative Spectroscopy & Radiative Transfer, Vol. 40, No. 2, 1988, pp. 123-129. >°Nelson, H. F., Look, D. C., Jr., and Crosbie, A. L., "TwoDimensional Radiative Back-Scattering from Optically Thick Media," Journal of Heat Transfer, Vol. 108, 1986, pp. 619-625. "Shouman, S. M., and Ozisik, M, N., "Radiation Transfer in an Isotropically Scattering Two-Region Slab with Reflecting Boundaries," Journal of Quantitative Spectroscopy & Radiative Transfer, Vol. 26, 1981, pp. 1-9. 12 Reguigui, N. M., "Radiative Transfer in a Four-Layer Cylindrical Medium with Reflecting Boundaries," Master of Science Thesis, Oklahoma State Univ., 1990. 13 Reguigui, N. M., and Dougherty, R. L., "Two-Dimensional Radiative Transfer in a Cylindrical Layered Medium with Reflective Boundaries," AIAA Paper 90-1779, AIAA/ASME 5th Joint Thermophysics and Heat Transfer Conf., Seattle, WA, 1990. 14 Dougherty, R. L., private communication, 1990. 15 Dorri-Nowkoorani, F., Reguigui, N. M., and Dougherty, R. L., "Back-Scattering of a Laser Beam from a Layered Cylindrical Medium: Refractive Index Effects," AIAA Paper 90-1764, AIAA/ASME 5th Joint Thermophysics and Heat Transfer Conf., Seattle, WA, 1990. 16 Crosbie, A. L., and Koewing, J. W., "Two-Dimensional Radiative Transfer in a Finite Scattering Planar Medium," Journal of Quantitative Spectroscopy & Radiative Transfer, Vol. 21, 1979, pp. 573595. Recommended Reading from the AIAA **MMM M Progress in Astronautics and Aeronautics Series . . . Thermal Design of Aeroassisted Orbital Transfer Vehicles H. F. Nelson, editor Underscoring the importance of sound thermophysical knowledge in spacecraft design, this volume emphasizes effective use of numerical analysis and presents recent advances and current thinking about the design of aeroassisted orbital transfer vehicles (AOTVs). Its 22 chapters cover flow field analysis, trajectories (including impact of atmospheric uncertainties and viscous interaction effects), thermal protection, and surface effects such as temperature-dependent reaction rate expressions for oxygen recombination; surface-ship equations for low-Reynolds-number multicomponent air flow, rate chemistry in flight regimes, and noncatalytic surfaces for metallic heat shields. TO ORDER: Write, Phone, or FAX: American Institute of Aeronautics and Astronautics c/o Publications Customer Service, 9 Jay Qould Ct., P.O. Box 753, Waldorf, MD 20604 Phone: 301/645-5643 or 1-800/682-AIAA, Dept. 415 • FAX: 301/843-0159 Sales Tax: CA residents, 8.25%; DC, 6%. For shipping and handling add $4.75 for 1-4 books (call for rates for higher quanties). Orders under $50.00 must be prepaid. Foreign orders must be prepaid. Please allow 4 weeks for delivery. Prices are subject to change without notice. Returns will be accepted within 15 days. View publication stats 241 1985 566 pp., illus. Hardback ISBN 0-915928-94-9 AIAA Members $54.95 Nonmembers $81.95 Order Number V-96