# Extension of the Reduced Integration Scheme to Calculate the Direct Exchange Areas in 3D Rectangular Enclosures with Nonscattering Media

Authors: Guo-Jun Li; Ben-Wen Li; Ya-Song Sun
Journal: Mathematical Problems in Engineering
 (2015)

---

## Abstract

The evaluation of direct exchange areas (DEA) in zonal method is the most important task due to the heavy computer cost of multi-integrals together with the existing of singularities. A technique of variable transformation to reduce the fold of integrals, which was developed originally by Erkku (1959) to calculate the DEAs of uniformly zonal dividing cylindrical system, was extended by Tian and Chiu (2003) for nonuniformly zonal dividing cylindrical system with large thermal gradients. In this paper, we further extend the reduced integration scheme (RIS) to calculate the DEAs in three-dimensional rectangular system. The detail deductions of six-, five-, and fourfold integrals to threefold ones are presented; the DEAs in a rectangular system with assumption of gray medium are computed by the Gaussian quadrature integration (GQI) and the RIS comparatively. The comparisons reveal that the RIS can provide remarkable higher accuracy and efficiency than GQI. More interestingly and practicably, the singularities of DEAs can be decomposed and weakened obviously by RIS.


---

## Body

## 1. Introduction

In recent decades, the zonal method [1] has been being an outstanding tool for thermal radiation calculations, whether for the developments of the numerical methods themselves [2–10] or the applications in thermal engineering [11–15]. It is well known that the most heavy computation cost of zonal method comes from the direct exchange areas (DEAs), which include the radiative properties of the participating medium and boundaries, together with the geometrical relations. After the DEAs, the total exchange areas (TEAs) can be solved, then the solution of zonal method. Usually, the calculations of DEAs involve six-, five-, and fourfold integrals for three-dimensional systems, as can be seen in the following contents. Further, the analytical or exact solutions of DEAs are difficult or impossible for most cases. Therefore, considerable efforts have been made to reduce or simplify the computation cost of DEAs till now.Larsen and Howell [3] defined the exchange factor method to reduce the computer cost using the reciprocity and symmetry. Sika [4] defined the DEAs of all kinds, say, surface to surface, surface to volume, volume to volume in axisymmetric cylindrical geometries with gray gas, and reduced multifold (four, five, and six) integrals to single- or double-fold integrals, and simplified formulae with high efficiency. Ma [5] illustrated that the GZM (generalized zoning method) can treat the anisotropic scattering using scattering-exchange areas; the GZM was used in one-dimensional case and compared with existing solutions. Goyhénèche and Sacadura [6] presented the explicit matrix relation for total exchange areas; the integration order can be reduced; then the TEAs in anisotropically scattering medium and anisotropically reflecting walls were obtained, and the weighted sum of gray gases model (WSGG) could be adopted in the zonal method. Tian and Chiu [7] deduced the formulas to reduce the multifold to less fold integrations by variable transformation for nonuniform zones in cylindrical systems. Two years later, Tian and Chiu [8] proposed a new approach to compute the DEAs by hybridization of the finite volume method (FVM) with the midpoint integration scheme; the FVM is for adjacent and overlapping zones; the midpoint integration scheme is for nonadjacent zones. El Hitti et al. [9] used zone method plus the thin layer approximation to simulate the glass sheets under high temperature condition; they evaluate the DEAs by flux planes approximation. At the same time, El Hitti et al. [10] evaluated the TEAs from DEAs of nongray system by replating algorithm.The applications of zonal method can be found in many thermal processes, whether for offline analysis of parameter design or for online control of operating parameter optimization. Chapman et al. [11] performed a parametric study for a batch reheating furnace; the DEAs were evaluated by Monte Carlo method and the WSGG (weighted sum of gray gas) was adopted to consider different species of gases. Jeong and Ha [12] simulated the continuously annealing furnace using zonal method for radiation with the measured gas temperature and imaginary surface and discussed the correction of the soot effects on radiative properties. Barr [13] reviewed the thermal model of steady-state pusher-type reheating furnace; more details like the treatment of combination of radiation and forced convection were mentioned also. Chen and Jaluria [14] used the parallelization for DEA computation in their research. Steinboeck et al. [15] presented the thermal model of a slab reheating furnace, in which the heat conduction in the slab was solved by the weighted residual approach and the radiation was solved by the net radiation method based on zone method.Recently, we have made the mathematical modeling of thermal processes both in conventional reheating furnace and regenerative reheating furnace [16] using zonal method to consider the radiative heat transfer. After cautious analyzing and comparing of above-mentioned approaches, the reduced integration scheme (RIS), which was described in detail in [7], is extended to simplify the computation of DEAs in rectangular systems, for uniform or nonuniform zone arrangement. Although Tian and Chiu [7] have already pointed out that this technique can be easily extended to three-dimensional rectangular enclosures, we still would like to give the detail deduction and results. More interestingly, the believable evidences of accuracy and efficiency will be given, especially for the treatment of singularities. This paper is organized as follows. In Section 2, the definition of DEAs for different zone pairs and detailed formulae are given. In Section 3, the detailed description of the RIS for DEAs is described. In Section 4, verifications and comparisons for specified cases are made. Finally, the last section gives the conclusions.


## 2. Definitions of DEAs in Rectangular Enclosures

In fact, the definitions of DEAs of all kinds of zone pair, say, surface-surface, gas-gas (or volume-volume), and gas-surface (or volume-surface), can be found in [17, 18]. Here we present them briefly for the subsequent deductions. Figure 1 shows the radiative exchange between surface zones si and sj, whose areas are Fi and Fj, respectively. Then the formula of DEA between them can be obtained with the guarantee of energy conservation for semitransparent medium of extinction coefficient β. While, as indicated in the title, the media considered here are nonscattering, hereafter, only the absorption coefficient κa is used. Consider(1)sisj¯=∫Fi∫Fjcos⁡ηicos⁡ηjπr2e-κardfjdfi.Figure 1
Radiative exchange between surface zonesi with area Fi and surface zone sj with area Fj.Similarly, as shown in Figure2, the formula of DEA between two gas zones gi and gj, whose volumes are Vi and Vj individually, is as follows:(2)gigj¯=∫Vi∫Vje-κarκa2πr2dvjdvi.Figure 2
Radiative exchange between gas zonegi with volume Vi and gas zone gj with volume Vj.The formula of DEA between one gas zonegi and one surface zone sj, whose volume is Vi and whose area is Fj individually, is as follows (see Figure 3):(3)gisj¯=∫Vi∫Fje-κarcos⁡ηπr2κadfjdvi.Figure 3
Radiative exchange between gas zonegi with volume Vi and surface zone sj with area Fj.Due to the feature of radiative heat transfer, there exists the reciprocity for DEAs and it implies(4)sisj¯=sjsi¯,gigj¯=gjgi¯,gisj¯=sjgi¯.The energy conservation of radiative heat transfer leads to the following relations for surface zone of areaFi and gas zone of volume Vi, respectively:(5)∑jsisj¯+∑jsigj¯=Fi,(6)∑jgisj¯+∑jgigj¯=4κaVi.According to the reciprocity, only half of the DEAs need to be evaluated. The energy conservation can be used to check the accuracy of computation.


## 3. Simplified Calculations of DEAs by Reduced Integration Scheme

### 3.1. Double-Fold Integral into Three Single-Fold Integrals

As described in [7], the original technique to transform the double-fold integral into three single-fold integrals was proposed by Erkku [19]. Here we briefly give its principle.There is a twofold integral; for example,(7)I=∫0a∫0bfx-ydxdy,where the assumptions are made that a and b are constants, b≥a>0, and f is a function of x-y.Setz=a+x-y; we get(8)I=∫0a∫0bfx-ydxdy=∫0a∫a-ya+b-yfz-adzdy.After the variable substitution, the integral can be executed in three regions as shown in Figure4:(9)I=∫0a∫a-yafz-adzdy+∫0a∫abfz-adzdy+∫0a∫ba+b-yfz-adzdy.Figure 4
Regions of integration after the transformation.Changing the order of integrationdzdy, we get(10)I=∫0a∫a-zafz-adydz+∫ab∫0afz-adydz+∫ba+b∫0b+a-zfz-adydz.Sincef is just a function of z, we have(11)I=∫0azfz-adz+a∫abfz-adz+∫ba+ba+b-zfz-adz.Please note that this variable transformation is available for all computation of DEAs in rectangular systems. Now the twofold integral (7) is transformed to three onefold integrals. As pointed out by Tian and Chiu [7], the computation time can be saved remarkably using (11) instead of (7).


### 3.2. Simplification of DEAs in Rectangular System

Considering a rectangular enclosure as shown in Figure5, where the numbers mean the zone number of surface zone and the numbers with circle mean the zone number of gas zone, the simplification of multifold integral of DEAs is given for gas-gas, gas-surface, and surface-surface, sequentially. Please note that the lengths of each zone (gas or surface) in three directions are not necessary to be equal; in another words, the zone dividing can be nonuniform.Figure 5
Zone arrangement in a rectangular enclosure. (a) Gas zonesg1~g8, (b) surface zones s1~s12 and s21~s24 which are parallel to the x-axes, and (c) surface zones s13~s20 which are perpendicular to the x-axes.
(a)
(b)
(c)(i) Gas-Gas Zones
g
i
g
j
¯. In Cartesian coordinates system, formula (2) can be rewritten as(12)gigj¯=∫xi1xi1+Δxi∫xj1xj1+Δxj∫yi1yi1+Δyi∫yj1yj1+Δyj∫zi1zi1+Δzi∫zj1zj1+Δzjfrdzjdzidyjdyidxjdxi,where fr=κa2e-κar/πr2, r=xi-xj2+yi-yj2+zi-zj2,xi,xj,yi,yj,zi,zj are coordinates corresponding to gas zone gi and gas zone gj, xi1,xj1,yi1,yj1,zi1,zj1 are lower limits of integral, and Δxi,Δxj,Δyi,Δyj,Δzi,Δzj are lengths of gi and gj in three directions, respectively.Set(13)xi=xi1+x1,xj=xj1+x2yi=yi1+y1,yj=yj1+y2zi=zi1+z1,zj=zj1+z2and assume Δxj≥Δxi; then one can get the following after integral transformation:(14)gigj¯=∫0Δxi∫0Δxj∫0Δyi∫0Δyj∫0Δzi∫0Δzjfrdz2dz1dy2dy1dx2dx1.This time,(15)r=x1-x2+xi1-xj12+y1-y2+yi1-yj12+z1-z2+zi1-zj121/2.Setz=z1-z2, x=x1-x2, y=y1-y2 and use the RIS which is described in Section 3.1; formula (14) can be expressed as(16)gigj¯=∑n=02∑l=02∑k=02∫anan+1∫blbl+1∫ckck+1φ(x,n)ψ(y,l)η(z,k)×f(r)dzdydx,where (they are used later again)(17)an=-Δxin=00n=1Δxj-Δxin=2Δxjn=3,bl=-Δyil=00l=1Δyj-Δyil=2Δyjl=3,ck=-Δzik=00k=1Δzj-Δzik=2Δzjk=3,φx,n=Δxi+xn=0Δxin=1Δxj-xn=2,ψy,l=Δyi+yl=0Δyil=1Δyj-yl=2,ηz,k=Δzi+zk=0Δzik=1Δzj-zk=2.If the zone dividing is uniform in all directions, the above integrals can be further simplified as(18)gigj¯=∫-Δxi0∫-Δyi0∫-Δzi0Δxi+xΔyi+yΔzi+z×frdzdydx+∫0Δxj∫-Δyi0∫-Δzi0Δxj-xΔyi+yΔzi+z×frdzdydx+∫-Δxi0∫0Δyj∫-Δzi0Δxi+xΔyj-yΔzi+z×frdzdydx+∫0Δxj∫0Δyj∫-Δzi0Δxj-xΔyj-yΔzi+z×frdzdydx+∫-Δxi0∫-Δyi0∫0ΔzjΔxi+xΔyi+yΔzj-z×frdzdydx+∫0Δxj∫-Δyi0∫0ΔzjΔxj-xΔyi+yΔzj-z×frdzdydx+∫-Δxi0∫0Δyj∫0ΔzjΔxi+xΔyj-yΔzj-z×frdzdydx+∫0Δxj∫0Δyj∫0ΔzjΔxj-xΔyj-yΔzj-z×frdzdydx,where(19)r=x+xi1-xj12+y+yi1-yj12+z+zi1-zj12.(ii) Gas-Surface Zones
g
i
s
j
¯. Taking the same meanings of symbols after formula (12), formula (3) becomes(20)gisj¯=∫xi1xi1+Δxi∫yi1yi1+Δyi∫yj1yj1+Δyj∫zi1zi1+Δzi∫zj1zj1+Δzjfrdzjdzidyjdyidxi.There exist six cases for the surface zonesj according to its different locations when the gas zone gi is fixed; they are as follows:(A)
the surface zones are perpendicular to thex-axis but away from the origin, as the zones numbered 13~16 in Figure 5(c), for example;(B)
the surface zones are perpendicular to thex-axis but in the same plane with the origin, as the zones numbered 17~20, which are opposite to the zones of case (A), in Figure 5(c), for example;(C)
the surface zones are perpendicular to they-axis but away from the origin, as the zones numbered 5~8 in Figure 5(b), for example;(D)
the surface zones are perpendicular to they-axis but in the same plane with the origin, as the zones numbered 1~4, which are opposite to the zones of case (C), in Figure 5(b), for example;(E)
the surface zones are perpendicular to thez-axis but away from the origin, as the zones numbered 21~24 in Figure 5(b), for example;(F)
the surface zones are perpendicular to thez-axis but in the same plane with the origin, as the zones numbered 9~12, which are opposite to the zones of case (E), in Figure 5(b), for example.Consider the case (A); now the cosine function in (3) is(21)cos⁡η=X-xiX-xi2+yj-yi2+zj-zi2.Then(22)fr=e-κarκaX-xiπr3,r=X-xi1-x12+yj1+y2-yi1-y12+zj1+z2-zi1-z121/2.Set(23)xi=xi1+x1,yi=yi1+y1yj=yj1+y2,zi=zi1+z1zj=zj1+z2,x=x2-x1,y=y2-y1,z=z2-z1.After changing the order of integration, one can get(24)gisj¯=∫0Δxi∫0Δyi∫0Δyj∫0Δzi∫0Δzjfrdz2dz1dy2dy1dx1,where(25)fr=e-κarκaX-xi1-x1πr3,r=X-xi1-x12+y+yj1-yi12+z+zj1-zi12,and the capital X means the distance along the x-axis.Using the RIS, formula (24) becomes(26)gisj¯=∑l=02∑k=02∫0Δxi∫blbl+1∫ckck+1ψ(y,l)η(z,k)f(r)dzdydx1.If the gas zone can be projected onto the surface zone completely along one direction, or its projection can overlay the surface zone completely after parallel moving along one direction, formula (26) can be further simplified as(27)gisj¯=∫0Δxi∫-Δyi0∫-Δzi0Δyi+yΔzi+zfrdzdy+∫0Δyj∫-Δzi0Δyj-yΔzi+zfrdzdy+∫-Δyi0∫0ΔzjΔyi+yΔzi-zfrdzdy+∫0Δyj∫0ΔzjΔyj-yΔzi-z×∫0Δyj∫0ΔzjΔyj-yΔzi-zfrdzdydx1.As to cases (B)~(F), the cosine functions andfr are different, but the similar formulae as (26) and (27) can be obtained.(iii) Surface-Surface Zones
s
i
s
j
¯. As shown in Figure 5, there are two cases for the pairs of surface-surface zones as follows:(A)
one surface zone is perpendicular to the other;(B)
one surface zone is parallel to the other.For case (A), taking one surface zone on thez-y plane and another one on the x-y plane, for example, the cosine functions in formula (1) are(28)cos⁡ηi=xjxj2+yj-yi2+zi2,cos⁡ηj=zixj2+yj-yi2+zi2.Then formula (1) can be rewritten as(29)sisj¯=∫zi1zi1+Δzi∫yi1yi1+Δyi∫yj1yj1+Δyj∫xj1xj1+Δxjfrdxjdyjdyidzi,where(30)fr=e-κarcos⁡ηicos⁡ηjπr2,r=xj2+yj-yi2+zi2.Set(31)xj=xj1+x2,yi=yi1+y1,yj=yj1+y2,zi=zi1+z1,zj=zj1+z2.After changing the order of integration, one gets(32)sisj¯=∫0Δzi∫0Δyi∫0Δyj∫0Δxjfrdx2dy2dy1dz1.Now(33)cos⁡ηi=xj1+x2r,cos⁡ηj=zi1+z1r,r=xj1+x22+yj1-yi1+y2-y12+zi1+z12.Using the RIS, formula (32) becomes(34)sisj¯=∫0Δzi∫-Δyi0∫0ΔxjΔyi+yfrdx2dydz1+∫0Δzi∫0Δyj∫0ΔxjΔyj-yfrdx2dydz1.For other situations, one surface zone on thez-y plane or parallel to z-y plane, and another one on the x-y plane or parallel to x-y plane, for example, the deductions can be carried out similarly.For case (B), taking two surface zones, which are perpendicular to thex-axis, for example, the cosine functions in formula (1) are(35)cos⁡ηi=XX2+yj-yi2+zj-zi2,cos⁡ηj=cos⁡ηi.Then formula (1) becomes(36)sisj¯=∫zi1zi1+Δzi∫zj1zj1+Δzj∫yi1yi1+Δyi∫yj1yj1+Δyjfrdyjdyidzjdzi,where(37)fr=e-κarcos⁡ηi2πr2,r=X2+yj-yi2+zj-zi2.Set(38)yi=yi1+y1,yj=yj1+y2,zi=zi1+z1,zj=zj1+z2.After changing the order of integration, one gets(39)sisj¯=∫0Δzi∫0Δzj∫0Δyi∫0Δyjfrdyjdyidzjdzi,where(40)fr=e-κarcos⁡ηi2πr2,r=X2+yj1+y2-yi1-y12+zj1+z2-zi1-z12.Using the RIS, formula (39) becomes(41)sisj¯=∑k=02∑l=02∫ckck+1∫blbl+1ψ(y,l)η(z,k)f(r)dydz.If the lengths of both surface zones are equal iny- and z-directions, separately, formula (41) can be further simplified as(42)sisj¯=∫-Δzi0∫-Δyi0Δyi+yΔzi+zfrdydz+∫0Δzj∫-Δyi0Δyi+yΔzj-zfrdydz+∫-Δzi0∫0ΔyjΔyj-yΔzi+zfrdydz+∫0Δzj∫0ΔyjΔyj-yΔzj-zfrdydz.The formulae of other situations for case (B) can be deduced in the same way. We state that the situations with shielding effect do not belong to the cases studied here.Up to now we know that, if the zone grids along three directions are uniform, the integrals of gas-gas zones can be reduced from one sixfold integral to eight threefold integrals and those of gas-surface zones from one fivefold to four threefold integrals. For the case of parallel surface zones, one fourfold integral is reduced to four twofold integrals; for the case of perpendicular surface zones, one fourfold integral is reduced to two threefold integrals. If the quadrature point numberM is chosen for the single-fold integral, the M6 and M3 points are needed for a sixfold integral and a threefold integral, separately. Therefore, the quadrature points needed for one DEA of gas-gas zones can be reduced from M6 to 8M3, from M5 to 4M3 for one DEA of gas-surface zones, and from M4 to 2M2 (perpendicular) or 4M2 (parallel) for DEA of surface-surface zones. Once the integral regions are determined, the greater M is, the higher the accuracy is. On the other hand, more computational time can be saved with the increasing of M, and more higher accuracy can be achieved as well.


## 3.1. Double-Fold Integral into Three Single-Fold Integrals

As described in [7], the original technique to transform the double-fold integral into three single-fold integrals was proposed by Erkku [19]. Here we briefly give its principle.There is a twofold integral; for example,(7)I=∫0a∫0bfx-ydxdy,where the assumptions are made that a and b are constants, b≥a>0, and f is a function of x-y.Setz=a+x-y; we get(8)I=∫0a∫0bfx-ydxdy=∫0a∫a-ya+b-yfz-adzdy.After the variable substitution, the integral can be executed in three regions as shown in Figure4:(9)I=∫0a∫a-yafz-adzdy+∫0a∫abfz-adzdy+∫0a∫ba+b-yfz-adzdy.Figure 4
Regions of integration after the transformation.Changing the order of integrationdzdy, we get(10)I=∫0a∫a-zafz-adydz+∫ab∫0afz-adydz+∫ba+b∫0b+a-zfz-adydz.Sincef is just a function of z, we have(11)I=∫0azfz-adz+a∫abfz-adz+∫ba+ba+b-zfz-adz.Please note that this variable transformation is available for all computation of DEAs in rectangular systems. Now the twofold integral (7) is transformed to three onefold integrals. As pointed out by Tian and Chiu [7], the computation time can be saved remarkably using (11) instead of (7).


## 3.2. Simplification of DEAs in Rectangular System

Considering a rectangular enclosure as shown in Figure5, where the numbers mean the zone number of surface zone and the numbers with circle mean the zone number of gas zone, the simplification of multifold integral of DEAs is given for gas-gas, gas-surface, and surface-surface, sequentially. Please note that the lengths of each zone (gas or surface) in three directions are not necessary to be equal; in another words, the zone dividing can be nonuniform.Figure 5
Zone arrangement in a rectangular enclosure. (a) Gas zonesg1~g8, (b) surface zones s1~s12 and s21~s24 which are parallel to the x-axes, and (c) surface zones s13~s20 which are perpendicular to the x-axes.
(a)
(b)
(c)(i) Gas-Gas Zones
g
i
g
j
¯. In Cartesian coordinates system, formula (2) can be rewritten as(12)gigj¯=∫xi1xi1+Δxi∫xj1xj1+Δxj∫yi1yi1+Δyi∫yj1yj1+Δyj∫zi1zi1+Δzi∫zj1zj1+Δzjfrdzjdzidyjdyidxjdxi,where fr=κa2e-κar/πr2, r=xi-xj2+yi-yj2+zi-zj2,xi,xj,yi,yj,zi,zj are coordinates corresponding to gas zone gi and gas zone gj, xi1,xj1,yi1,yj1,zi1,zj1 are lower limits of integral, and Δxi,Δxj,Δyi,Δyj,Δzi,Δzj are lengths of gi and gj in three directions, respectively.Set(13)xi=xi1+x1,xj=xj1+x2yi=yi1+y1,yj=yj1+y2zi=zi1+z1,zj=zj1+z2and assume Δxj≥Δxi; then one can get the following after integral transformation:(14)gigj¯=∫0Δxi∫0Δxj∫0Δyi∫0Δyj∫0Δzi∫0Δzjfrdz2dz1dy2dy1dx2dx1.This time,(15)r=x1-x2+xi1-xj12+y1-y2+yi1-yj12+z1-z2+zi1-zj121/2.Setz=z1-z2, x=x1-x2, y=y1-y2 and use the RIS which is described in Section 3.1; formula (14) can be expressed as(16)gigj¯=∑n=02∑l=02∑k=02∫anan+1∫blbl+1∫ckck+1φ(x,n)ψ(y,l)η(z,k)×f(r)dzdydx,where (they are used later again)(17)an=-Δxin=00n=1Δxj-Δxin=2Δxjn=3,bl=-Δyil=00l=1Δyj-Δyil=2Δyjl=3,ck=-Δzik=00k=1Δzj-Δzik=2Δzjk=3,φx,n=Δxi+xn=0Δxin=1Δxj-xn=2,ψy,l=Δyi+yl=0Δyil=1Δyj-yl=2,ηz,k=Δzi+zk=0Δzik=1Δzj-zk=2.If the zone dividing is uniform in all directions, the above integrals can be further simplified as(18)gigj¯=∫-Δxi0∫-Δyi0∫-Δzi0Δxi+xΔyi+yΔzi+z×frdzdydx+∫0Δxj∫-Δyi0∫-Δzi0Δxj-xΔyi+yΔzi+z×frdzdydx+∫-Δxi0∫0Δyj∫-Δzi0Δxi+xΔyj-yΔzi+z×frdzdydx+∫0Δxj∫0Δyj∫-Δzi0Δxj-xΔyj-yΔzi+z×frdzdydx+∫-Δxi0∫-Δyi0∫0ΔzjΔxi+xΔyi+yΔzj-z×frdzdydx+∫0Δxj∫-Δyi0∫0ΔzjΔxj-xΔyi+yΔzj-z×frdzdydx+∫-Δxi0∫0Δyj∫0ΔzjΔxi+xΔyj-yΔzj-z×frdzdydx+∫0Δxj∫0Δyj∫0ΔzjΔxj-xΔyj-yΔzj-z×frdzdydx,where(19)r=x+xi1-xj12+y+yi1-yj12+z+zi1-zj12.(ii) Gas-Surface Zones
g
i
s
j
¯. Taking the same meanings of symbols after formula (12), formula (3) becomes(20)gisj¯=∫xi1xi1+Δxi∫yi1yi1+Δyi∫yj1yj1+Δyj∫zi1zi1+Δzi∫zj1zj1+Δzjfrdzjdzidyjdyidxi.There exist six cases for the surface zonesj according to its different locations when the gas zone gi is fixed; they are as follows:(A)
the surface zones are perpendicular to thex-axis but away from the origin, as the zones numbered 13~16 in Figure 5(c), for example;(B)
the surface zones are perpendicular to thex-axis but in the same plane with the origin, as the zones numbered 17~20, which are opposite to the zones of case (A), in Figure 5(c), for example;(C)
the surface zones are perpendicular to they-axis but away from the origin, as the zones numbered 5~8 in Figure 5(b), for example;(D)
the surface zones are perpendicular to they-axis but in the same plane with the origin, as the zones numbered 1~4, which are opposite to the zones of case (C), in Figure 5(b), for example;(E)
the surface zones are perpendicular to thez-axis but away from the origin, as the zones numbered 21~24 in Figure 5(b), for example;(F)
the surface zones are perpendicular to thez-axis but in the same plane with the origin, as the zones numbered 9~12, which are opposite to the zones of case (E), in Figure 5(b), for example.Consider the case (A); now the cosine function in (3) is(21)cos⁡η=X-xiX-xi2+yj-yi2+zj-zi2.Then(22)fr=e-κarκaX-xiπr3,r=X-xi1-x12+yj1+y2-yi1-y12+zj1+z2-zi1-z121/2.Set(23)xi=xi1+x1,yi=yi1+y1yj=yj1+y2,zi=zi1+z1zj=zj1+z2,x=x2-x1,y=y2-y1,z=z2-z1.After changing the order of integration, one can get(24)gisj¯=∫0Δxi∫0Δyi∫0Δyj∫0Δzi∫0Δzjfrdz2dz1dy2dy1dx1,where(25)fr=e-κarκaX-xi1-x1πr3,r=X-xi1-x12+y+yj1-yi12+z+zj1-zi12,and the capital X means the distance along the x-axis.Using the RIS, formula (24) becomes(26)gisj¯=∑l=02∑k=02∫0Δxi∫blbl+1∫ckck+1ψ(y,l)η(z,k)f(r)dzdydx1.If the gas zone can be projected onto the surface zone completely along one direction, or its projection can overlay the surface zone completely after parallel moving along one direction, formula (26) can be further simplified as(27)gisj¯=∫0Δxi∫-Δyi0∫-Δzi0Δyi+yΔzi+zfrdzdy+∫0Δyj∫-Δzi0Δyj-yΔzi+zfrdzdy+∫-Δyi0∫0ΔzjΔyi+yΔzi-zfrdzdy+∫0Δyj∫0ΔzjΔyj-yΔzi-z×∫0Δyj∫0ΔzjΔyj-yΔzi-zfrdzdydx1.As to cases (B)~(F), the cosine functions andfr are different, but the similar formulae as (26) and (27) can be obtained.(iii) Surface-Surface Zones
s
i
s
j
¯. As shown in Figure 5, there are two cases for the pairs of surface-surface zones as follows:(A)
one surface zone is perpendicular to the other;(B)
one surface zone is parallel to the other.For case (A), taking one surface zone on thez-y plane and another one on the x-y plane, for example, the cosine functions in formula (1) are(28)cos⁡ηi=xjxj2+yj-yi2+zi2,cos⁡ηj=zixj2+yj-yi2+zi2.Then formula (1) can be rewritten as(29)sisj¯=∫zi1zi1+Δzi∫yi1yi1+Δyi∫yj1yj1+Δyj∫xj1xj1+Δxjfrdxjdyjdyidzi,where(30)fr=e-κarcos⁡ηicos⁡ηjπr2,r=xj2+yj-yi2+zi2.Set(31)xj=xj1+x2,yi=yi1+y1,yj=yj1+y2,zi=zi1+z1,zj=zj1+z2.After changing the order of integration, one gets(32)sisj¯=∫0Δzi∫0Δyi∫0Δyj∫0Δxjfrdx2dy2dy1dz1.Now(33)cos⁡ηi=xj1+x2r,cos⁡ηj=zi1+z1r,r=xj1+x22+yj1-yi1+y2-y12+zi1+z12.Using the RIS, formula (32) becomes(34)sisj¯=∫0Δzi∫-Δyi0∫0ΔxjΔyi+yfrdx2dydz1+∫0Δzi∫0Δyj∫0ΔxjΔyj-yfrdx2dydz1.For other situations, one surface zone on thez-y plane or parallel to z-y plane, and another one on the x-y plane or parallel to x-y plane, for example, the deductions can be carried out similarly.For case (B), taking two surface zones, which are perpendicular to thex-axis, for example, the cosine functions in formula (1) are(35)cos⁡ηi=XX2+yj-yi2+zj-zi2,cos⁡ηj=cos⁡ηi.Then formula (1) becomes(36)sisj¯=∫zi1zi1+Δzi∫zj1zj1+Δzj∫yi1yi1+Δyi∫yj1yj1+Δyjfrdyjdyidzjdzi,where(37)fr=e-κarcos⁡ηi2πr2,r=X2+yj-yi2+zj-zi2.Set(38)yi=yi1+y1,yj=yj1+y2,zi=zi1+z1,zj=zj1+z2.After changing the order of integration, one gets(39)sisj¯=∫0Δzi∫0Δzj∫0Δyi∫0Δyjfrdyjdyidzjdzi,where(40)fr=e-κarcos⁡ηi2πr2,r=X2+yj1+y2-yi1-y12+zj1+z2-zi1-z12.Using the RIS, formula (39) becomes(41)sisj¯=∑k=02∑l=02∫ckck+1∫blbl+1ψ(y,l)η(z,k)f(r)dydz.If the lengths of both surface zones are equal iny- and z-directions, separately, formula (41) can be further simplified as(42)sisj¯=∫-Δzi0∫-Δyi0Δyi+yΔzi+zfrdydz+∫0Δzj∫-Δyi0Δyi+yΔzj-zfrdydz+∫-Δzi0∫0ΔyjΔyj-yΔzi+zfrdydz+∫0Δzj∫0ΔyjΔyj-yΔzj-zfrdydz.The formulae of other situations for case (B) can be deduced in the same way. We state that the situations with shielding effect do not belong to the cases studied here.Up to now we know that, if the zone grids along three directions are uniform, the integrals of gas-gas zones can be reduced from one sixfold integral to eight threefold integrals and those of gas-surface zones from one fivefold to four threefold integrals. For the case of parallel surface zones, one fourfold integral is reduced to four twofold integrals; for the case of perpendicular surface zones, one fourfold integral is reduced to two threefold integrals. If the quadrature point numberM is chosen for the single-fold integral, the M6 and M3 points are needed for a sixfold integral and a threefold integral, separately. Therefore, the quadrature points needed for one DEA of gas-gas zones can be reduced from M6 to 8M3, from M5 to 4M3 for one DEA of gas-surface zones, and from M4 to 2M2 (perpendicular) or 4M2 (parallel) for DEA of surface-surface zones. Once the integral regions are determined, the greater M is, the higher the accuracy is. On the other hand, more computational time can be saved with the increasing of M, and more higher accuracy can be achieved as well.


## 4. Cases Study and Comparisons

Returning to Figure5 again, we take the cubic enclosure with the sizes of 2m×2m×2m and assume the medium is gray with κa=0.5m-1 to check the accuracy and efficiency of the RIS for DEAs. The cube is divided into 8 gas zones (numbers with circle) and 24 surface zones uniformly. Therefore the side length of each gas zone and surface zone is 1m. All the following comparisons are made on a personal computer with Intel(R) Core (TM)2 Duo CPU P8700 (2.53 GHz) and 3.0 GB (DDR2) memory.Same as in [7], the Gaussian quadrature integration is adopted as the direct method for DEAs themselves or the reduced three- or twofold integrals for comparisons. Although there may be some other simple methods which can be used for DEAs or the reduced three- or twofold integrals, we will not do much more comparisons between the RIS and them. In the following tables and figures, the symbols GQI and RIS mean the results by Gaussian quadrature integration and reduced integration scheme, respectively.


### 4.1. Accuracy Validation and Comparisons

Taking the number of integral stepsN=8 (the corresponding integral step length is 1/N=0.125m) in each dimension, the DEAs of gas zone g1 to all gas zones gi (i=1,…,8) in the enclosure are computed by GQI and RIS separately, and the results are listed in Table 1. It can be seen that both GQI and RIS almost give the same results except for the DEA of g1g1¯. The difference for g1g1¯ between by GQI and by RIS  comes from the so called “self-irradiation” which exits in g1g1¯. The “self-irradiation” is one kind of singularities, which happens in the calculation of DEAs of gas zone to itself. Those results of DEAs of gas zone g1 to all surface zones si (i=1,…,24) in the enclosure are listed in Table 2. It is clear that both GQI and RIS almost give the same results except for g1s1¯, g1s9¯, and g1s17¯. The differences between by GQI and by RIS are due to the singularity when surface zones are adjacent to the gas zone g1. Those results of DEAs of surface zone s1 to all surface zones si (i=1,…,24) in the enclosure are listed in Table 3. Again, both GQI and RIS almost give the same results except for s1s9¯ and s1s17¯. The differences between by GQI and by RIS are caused by the singularity that the surface zones s9 and s17 are just neighbor to s1.Table 1
DEAs of gas zoneg1 to all gas zones gi  (i=1,…,8) in the enclosure for N=8 in each dimension.
DEAs(×10-1)
g
1
g
1
¯
g
1
g
2
¯
g
1
g
3
¯
g
1
g
4
¯
g
1
g
5
¯
g
1
g
6
¯
g
1
g
7
¯
g
1
g
8
¯
GQI
3.55586
0.64400
0.64400
0.24433
0.64400
0.24433
0.24433
0.13122
RIS
3.79281
0.64415
0.64415
0.24433
0.64415
0.24433
0.24433
0.13123Table 2
DEAs of gas zoneg1 to all surface zones si  (i=1,…,24) in the enclosure for N=8 in each dimension.
DEAs(×10-1)
g
1
s
1
¯
g
1
s
2
¯
g
1
s
3
¯
g
1
s
4
¯
g
1
s
5
¯
g
1
s
6
¯
g
1
s
7
¯
g
1
s
8
¯
GQI
3.05222
0.43781
0.43781
0.14722
0.30129
0.16836
0.16836
0.10348
RIS
2.70782
0.43819
0.43819
0.14723
0.30129
0.16836
0.16836
0.10348
DEAs(×10-1)
g
1
s
9
¯
g
1
s
10
¯
g
1
s
11
¯
g
1
s
12
¯
g
1
s
13
¯
g
1
s
14
¯
g
1
s
15
¯
g
1
s
16
¯
GQI
3.05222
0.43781
0.43781
0.14722
0.30129
0.16836
0.16836
0.10348
RIS
2.70782
0.43819
0.43819
0.14723
0.30129
0.16836
0.16836
0.10348
DEAs(×10-1)
g
1
s
17
¯
g
1
s
18
¯
g
1
s
19
¯
g
1
s
20
¯
g
1
s
21
¯
g
1
s
22
¯
g
1
s
23
¯
g
1
s
24
¯
GQI
3.05222
0.43781
0.43781
0.14722
0.30129
0.16836
0.16836
0.10348
RIS
2.70782
0.43819
0.43819
0.14723
0.30129
0.16836
0.16836
0.10348Table 3
DEAs of surface zones1 to all surface zones si  (i=1,…,24) in the enclosure for N=8 in each dimension.
DEAs(×10-1)
s
1
s
1
¯
s
1
s
2
¯
s
1
s
3
¯
s
1
s
4
¯
s
1
s
5
¯
s
1
s
6
¯
s
1
s
7
¯
s
1
s
8
¯
GQI
0
0
0
0
0.24342
0.15657
0.15657
0.10507
RIS
0
0
0
0
0.24342
0.15657
0.15657
0.10507
DEAs(×10-1)
s
1
s
9
¯
s
1
s
10
¯
s
1
s
11
¯
s
1
s
12
¯
s
1
s
13
¯
s
1
s
14
¯
s
1
s
15
¯
s
1
s
16
¯
GQI
1.65273
0.15096
0.25190
0.07788
0.15096
0.11658
0.07788
0.07402
RIS
1.57526
0.15096
0.25234
0.07788
0.15096
0.11658
0.07788
0.07402
DEAs(×10-1)
s
1
s
17
¯
s
1
s
18
¯
s
1
s
19
¯
s
1
s
20
¯
s
1
s
21
¯
s
1
s
22
¯
s
1
s
23
¯
s
1
s
24
¯
GQI
1.65273
0.15096
0.25190
0.07788
0.15096
0.11658
0.07788
0.07402
RIS
1.57526
0.15096
0.25234
0.07788
0.15096
0.11658
0.07788
0.07402It is a common knowledge that the accuracy of integrals will increase with the number of integral stepN (the corresponding integral step length is 1/N). Here we choose g1g1¯, g1s1¯, and s1s9¯ as examples to compare the DEAs by GQI and by RIS against N in each dimension. The results are illustrated in Figure 6. Here, in plot (a), the results by GQI are given only for N being up to 16 because one week or even more days may be needed when N is up to 32 on the present personal computer. Obviously, the results by GQI are more sensitive to N, and the relative errors (that is, with respect to the final result of RIS at N=32) of g1g1¯ (plot (a)) and g1s1¯ (plot (b)) at N=16 are still greater than 20%, while the results of RIS almost go to convergence when N≥8. Note that the singularities are involved in both g1g1¯ and g1s1¯ as mentioned above. Thus we know that the RIS can provide much higher accuracy than the GQI especially for those DEAs in which the singularities are included.Figure 6
DEAs against the number of integral stepsN in each dimension by GQI and RIS: (a) g1g1¯, (b) g1s1¯, and (c) s1s9¯.
(a)
(b)
(c)Figure7 depicts the differences of DEAs by GQI and RIS for various N. Whether for DEAs of g1gi¯ (i=1,…,8) as shown in plot (a) or for DEAs of g1si¯ (i=1,…,24) as shown in plot (b), the differences between GQI and RIS are almost invisible except for those DEAs who include the singularities. However, the differences of both for DEAs with singularities, g1g1¯, g1s1¯, g1s9¯, and g1s17¯, for example, are decreasing remarkably when N increases from 2 to 8. For higher values of N, the compared differences are not given in Figure 7 since almost several weeks computer time is needed by GQI on our personal computer.Figure 7
Differences of DEAs against the number of integral stepsN=2,4,8 in each dimension by GQI and RIS: (a) g1gi¯ (i=1,…,8); (b) g1si¯ (i=1,…,24).
(a)
(b)As mentioned in Section2, the energy conservation of DEAs in one enclosure can be applied as the criterion to check the accuracy of numerical methods. Taking the gas zone g1 as the example, there is the exact value of 2 for formula (6); in other words, ∑i=18g1gi¯+∑i=124g1si¯=2m2. The results by GQI and RIS, respectively, are listed in Table 4 for different N. One can see that the relative errors of RIS with respect to the exact value are even lesser than 0.1% when N is only up to 4, while, for the same conditions, the relative errors of GQI are still greater than 1%. The reason may be mainly because of the eliminating of singularities of DEAs by GQI and RIS differently. In RIS, one integral region is always decomposed into several ones as described in Section 3.1; at the same time, the singularity is also decomposed and weakened.Table 4
Verifications and relative errors of conservation for∑i=18g1gi¯+∑i=124g1si¯ = 2 m2 with the change of N.
N
GQI
Relative errors of GQI (%)
RIS
Relative errors of RIS (%)
2
2.31407
15.70
1.99713
0.14
4
2.02711
1.36
1.99936
0.03
8
2.04049
2.03
2.00097
0.05As to the singularities that appeared in the computation of DEAs, there exit three cases. The first case is the DEA of one gas zone to itself, the second case is the DEA of one gas zone to its adjacent surface, and the third case is the DEA between two surfaces with an intersecting line. In the computations of all DEAs by RIS, the folds are always reduced. That means the possibilities of singularity can be reduced by RIS, just as said in [7] like “reduction of severity of singularities.”


### 4.2. CPU Time Comparisons

The best way to evaluate the efficiency of one approach is the comparison of CPU time. Again we first take the single DEA ofg1g1¯, g1s1¯, and g1s9¯, respectively, to compare the CPU time cost by GQI and RIS as N increases from 2 to 64. The results are listed in Table 5. There, the reason for the absence of some values by GQI is just because almost several weeks are needed for multifold of DEAs. For gas-gas DEAs, the CPU time of GQI is about 104 times that of RIS; for gas-surface and surface-surface DEAs, they are about 103 and 102 times, separately. This phenomenon reveals that the larger folds of the original integrals are, the more CPU time can be saved when using RIS.Table 5
Comparison of single DEA CPU time in seconds for differentN.
N
g
1
g
1
¯
g
1
s
1
¯
s
1
s
9
¯
GQI
RIS
GQI
RIS
GQI
RIS
2
0.0781
0.0000
0.0780
0.0000
0.0000
0.0000
4
1.3122
0.0000
0.0781
0.0000
0.0156
0.0000
8
84.5153
0.0156
2.5781
0.0156
0.1094
0.0000
16
5897.4692
0.1564
63.4372
0.0625
2.3281
0.0781
32
—
1.5621
2569.1563
0.6402
43.1250
0.7500
64
—
10.4374
—
8.6440
633.6250
5.9844Then taking the sum of DEAs ofgigj¯ (i,j=1,…,8) and gisj¯ (i=1,…,8;j=1,…,24) as examples, CPU time costs by both GQI and RIS for different N=2,4,8 are listed in Table 6. As one can see, at least one hundred times of CPU time can be saved using RIS instead of GQI for the same N. With the increasing of N, even more CPU time can be saved when using RIS.Table 6
Comparison of sum DEAs CPU time in seconds for differentN.
N
g
i
g
j
¯
(
i
,
j
=
1
,
…
,
8
)
g
i
s
j
¯
(
i
=
1
,
…
,
8
;
j
=
1
,
…
,
24
)
GQI
RIS
GQI
RIS
2
1.9066
0.0469
1.0165
0.0781
4
421.7854
0.2548
26.7033
0.6092
8
24649.4421
1.1723
9799.5624
12.8447


## 4.1. Accuracy Validation and Comparisons

Taking the number of integral stepsN=8 (the corresponding integral step length is 1/N=0.125m) in each dimension, the DEAs of gas zone g1 to all gas zones gi (i=1,…,8) in the enclosure are computed by GQI and RIS separately, and the results are listed in Table 1. It can be seen that both GQI and RIS almost give the same results except for the DEA of g1g1¯. The difference for g1g1¯ between by GQI and by RIS  comes from the so called “self-irradiation” which exits in g1g1¯. The “self-irradiation” is one kind of singularities, which happens in the calculation of DEAs of gas zone to itself. Those results of DEAs of gas zone g1 to all surface zones si (i=1,…,24) in the enclosure are listed in Table 2. It is clear that both GQI and RIS almost give the same results except for g1s1¯, g1s9¯, and g1s17¯. The differences between by GQI and by RIS are due to the singularity when surface zones are adjacent to the gas zone g1. Those results of DEAs of surface zone s1 to all surface zones si (i=1,…,24) in the enclosure are listed in Table 3. Again, both GQI and RIS almost give the same results except for s1s9¯ and s1s17¯. The differences between by GQI and by RIS are caused by the singularity that the surface zones s9 and s17 are just neighbor to s1.Table 1
DEAs of gas zoneg1 to all gas zones gi  (i=1,…,8) in the enclosure for N=8 in each dimension.
DEAs(×10-1)
g
1
g
1
¯
g
1
g
2
¯
g
1
g
3
¯
g
1
g
4
¯
g
1
g
5
¯
g
1
g
6
¯
g
1
g
7
¯
g
1
g
8
¯
GQI
3.55586
0.64400
0.64400
0.24433
0.64400
0.24433
0.24433
0.13122
RIS
3.79281
0.64415
0.64415
0.24433
0.64415
0.24433
0.24433
0.13123Table 2
DEAs of gas zoneg1 to all surface zones si  (i=1,…,24) in the enclosure for N=8 in each dimension.
DEAs(×10-1)
g
1
s
1
¯
g
1
s
2
¯
g
1
s
3
¯
g
1
s
4
¯
g
1
s
5
¯
g
1
s
6
¯
g
1
s
7
¯
g
1
s
8
¯
GQI
3.05222
0.43781
0.43781
0.14722
0.30129
0.16836
0.16836
0.10348
RIS
2.70782
0.43819
0.43819
0.14723
0.30129
0.16836
0.16836
0.10348
DEAs(×10-1)
g
1
s
9
¯
g
1
s
10
¯
g
1
s
11
¯
g
1
s
12
¯
g
1
s
13
¯
g
1
s
14
¯
g
1
s
15
¯
g
1
s
16
¯
GQI
3.05222
0.43781
0.43781
0.14722
0.30129
0.16836
0.16836
0.10348
RIS
2.70782
0.43819
0.43819
0.14723
0.30129
0.16836
0.16836
0.10348
DEAs(×10-1)
g
1
s
17
¯
g
1
s
18
¯
g
1
s
19
¯
g
1
s
20
¯
g
1
s
21
¯
g
1
s
22
¯
g
1
s
23
¯
g
1
s
24
¯
GQI
3.05222
0.43781
0.43781
0.14722
0.30129
0.16836
0.16836
0.10348
RIS
2.70782
0.43819
0.43819
0.14723
0.30129
0.16836
0.16836
0.10348Table 3
DEAs of surface zones1 to all surface zones si  (i=1,…,24) in the enclosure for N=8 in each dimension.
DEAs(×10-1)
s
1
s
1
¯
s
1
s
2
¯
s
1
s
3
¯
s
1
s
4
¯
s
1
s
5
¯
s
1
s
6
¯
s
1
s
7
¯
s
1
s
8
¯
GQI
0
0
0
0
0.24342
0.15657
0.15657
0.10507
RIS
0
0
0
0
0.24342
0.15657
0.15657
0.10507
DEAs(×10-1)
s
1
s
9
¯
s
1
s
10
¯
s
1
s
11
¯
s
1
s
12
¯
s
1
s
13
¯
s
1
s
14
¯
s
1
s
15
¯
s
1
s
16
¯
GQI
1.65273
0.15096
0.25190
0.07788
0.15096
0.11658
0.07788
0.07402
RIS
1.57526
0.15096
0.25234
0.07788
0.15096
0.11658
0.07788
0.07402
DEAs(×10-1)
s
1
s
17
¯
s
1
s
18
¯
s
1
s
19
¯
s
1
s
20
¯
s
1
s
21
¯
s
1
s
22
¯
s
1
s
23
¯
s
1
s
24
¯
GQI
1.65273
0.15096
0.25190
0.07788
0.15096
0.11658
0.07788
0.07402
RIS
1.57526
0.15096
0.25234
0.07788
0.15096
0.11658
0.07788
0.07402It is a common knowledge that the accuracy of integrals will increase with the number of integral stepN (the corresponding integral step length is 1/N). Here we choose g1g1¯, g1s1¯, and s1s9¯ as examples to compare the DEAs by GQI and by RIS against N in each dimension. The results are illustrated in Figure 6. Here, in plot (a), the results by GQI are given only for N being up to 16 because one week or even more days may be needed when N is up to 32 on the present personal computer. Obviously, the results by GQI are more sensitive to N, and the relative errors (that is, with respect to the final result of RIS at N=32) of g1g1¯ (plot (a)) and g1s1¯ (plot (b)) at N=16 are still greater than 20%, while the results of RIS almost go to convergence when N≥8. Note that the singularities are involved in both g1g1¯ and g1s1¯ as mentioned above. Thus we know that the RIS can provide much higher accuracy than the GQI especially for those DEAs in which the singularities are included.Figure 6
DEAs against the number of integral stepsN in each dimension by GQI and RIS: (a) g1g1¯, (b) g1s1¯, and (c) s1s9¯.
(a)
(b)
(c)Figure7 depicts the differences of DEAs by GQI and RIS for various N. Whether for DEAs of g1gi¯ (i=1,…,8) as shown in plot (a) or for DEAs of g1si¯ (i=1,…,24) as shown in plot (b), the differences between GQI and RIS are almost invisible except for those DEAs who include the singularities. However, the differences of both for DEAs with singularities, g1g1¯, g1s1¯, g1s9¯, and g1s17¯, for example, are decreasing remarkably when N increases from 2 to 8. For higher values of N, the compared differences are not given in Figure 7 since almost several weeks computer time is needed by GQI on our personal computer.Figure 7
Differences of DEAs against the number of integral stepsN=2,4,8 in each dimension by GQI and RIS: (a) g1gi¯ (i=1,…,8); (b) g1si¯ (i=1,…,24).
(a)
(b)As mentioned in Section2, the energy conservation of DEAs in one enclosure can be applied as the criterion to check the accuracy of numerical methods. Taking the gas zone g1 as the example, there is the exact value of 2 for formula (6); in other words, ∑i=18g1gi¯+∑i=124g1si¯=2m2. The results by GQI and RIS, respectively, are listed in Table 4 for different N. One can see that the relative errors of RIS with respect to the exact value are even lesser than 0.1% when N is only up to 4, while, for the same conditions, the relative errors of GQI are still greater than 1%. The reason may be mainly because of the eliminating of singularities of DEAs by GQI and RIS differently. In RIS, one integral region is always decomposed into several ones as described in Section 3.1; at the same time, the singularity is also decomposed and weakened.Table 4
Verifications and relative errors of conservation for∑i=18g1gi¯+∑i=124g1si¯ = 2 m2 with the change of N.
N
GQI
Relative errors of GQI (%)
RIS
Relative errors of RIS (%)
2
2.31407
15.70
1.99713
0.14
4
2.02711
1.36
1.99936
0.03
8
2.04049
2.03
2.00097
0.05As to the singularities that appeared in the computation of DEAs, there exit three cases. The first case is the DEA of one gas zone to itself, the second case is the DEA of one gas zone to its adjacent surface, and the third case is the DEA between two surfaces with an intersecting line. In the computations of all DEAs by RIS, the folds are always reduced. That means the possibilities of singularity can be reduced by RIS, just as said in [7] like “reduction of severity of singularities.”


## 4.2. CPU Time Comparisons

The best way to evaluate the efficiency of one approach is the comparison of CPU time. Again we first take the single DEA ofg1g1¯, g1s1¯, and g1s9¯, respectively, to compare the CPU time cost by GQI and RIS as N increases from 2 to 64. The results are listed in Table 5. There, the reason for the absence of some values by GQI is just because almost several weeks are needed for multifold of DEAs. For gas-gas DEAs, the CPU time of GQI is about 104 times that of RIS; for gas-surface and surface-surface DEAs, they are about 103 and 102 times, separately. This phenomenon reveals that the larger folds of the original integrals are, the more CPU time can be saved when using RIS.Table 5
Comparison of single DEA CPU time in seconds for differentN.
N
g
1
g
1
¯
g
1
s
1
¯
s
1
s
9
¯
GQI
RIS
GQI
RIS
GQI
RIS
2
0.0781
0.0000
0.0780
0.0000
0.0000
0.0000
4
1.3122
0.0000
0.0781
0.0000
0.0156
0.0000
8
84.5153
0.0156
2.5781
0.0156
0.1094
0.0000
16
5897.4692
0.1564
63.4372
0.0625
2.3281
0.0781
32
—
1.5621
2569.1563
0.6402
43.1250
0.7500
64
—
10.4374
—
8.6440
633.6250
5.9844Then taking the sum of DEAs ofgigj¯ (i,j=1,…,8) and gisj¯ (i=1,…,8;j=1,…,24) as examples, CPU time costs by both GQI and RIS for different N=2,4,8 are listed in Table 6. As one can see, at least one hundred times of CPU time can be saved using RIS instead of GQI for the same N. With the increasing of N, even more CPU time can be saved when using RIS.Table 6
Comparison of sum DEAs CPU time in seconds for differentN.
N
g
i
g
j
¯
(
i
,
j
=
1
,
…
,
8
)
g
i
s
j
¯
(
i
=
1
,
…
,
8
;
j
=
1
,
…
,
24
)
GQI
RIS
GQI
RIS
2
1.9066
0.0469
1.0165
0.0781
4
421.7854
0.2548
26.7033
0.6092
8
24649.4421
1.1723
9799.5624
12.8447


## 5. Conclusions

The reduced integration scheme is extended to calculate the direct exchange areas in three-dimensional rectangular system. The detailed reductions of multifold (six-, five-, and fourfold) integrals to threefold integrals are presented for nonuniform zonal dividing. A uniform zonal dividing case with gray medium is adopted to comparatively check the accuracy and efficiency of the reduced integration scheme. The standard Gaussian quadrature integration is also applied for the purpose of comparisons. The comparisons reveal that only10-4 times, 10-3 times, and 10-2 times CPU of Gaussian quadrature integration for DEAs of gas-gas with sixfold integrals, gas-surface with fivefold integrals, and surface-surface with fourfold integrals are needed, respectively. More attractively, the reduced integration scheme can remarkably decompose and weaken the singularities, which exist in the DEAs of adjacent zones or the volume zones themselves, obviously.



---

*Source: 703823-2015-03-24.xml*