A version of the moment method of calculating the transfer of selective radiation

222

Yu. D. Shmyglevskii

where fi and gi are the approximate respectively.

Equations

values at the point ti of the free terms of Eqs. (2) and (3)

(5) and (6) are solved by Gauss’s method. The approximate

values of the

u(t) and J/(t) in the whole interval [a, b] are found by means of a five-point parabolic

functions

interpolation,

and all the integrals encountered

along the curves I and I’ in the determination

the derivatives I,,‘, I,,’ were evaluated by a four-point

or five-point Gaussian formula.

For the case where the values of the gradient of the functional the boundary

of each of the domains 52, and Q2,8 iterations

of the functional

I on the curve I, was 1=0.103.10-‘,

of

I were taken at 20 points on

were performed.

The initial value

and on the curve Z8 the value was 1=O.2O2.1O-3

(the contours lo and I, are shown in Fig. 1). The values of the function By (x,y) on the segment I’ were found in the strip -2.739+,<-2.052

for the contour Z. and in the strip -2.152<&,<-2.146

for the contour I,. Translated by J. Berry REFERENCES

1.

TOZONI, 0. V. Calculation of electromagnetic fields on computers (Raschet electromagnetitnykh vychislitel’nykh mashinakh), “Tekhnika”, Kiev, 1967.

2.

KANTOROVICH, L. V. and AKILOV, G. P. Functional analysis in normed spaces (Funktsional’nyi analiz v normirovannykh prostranstvakh), Fizmatgiz, Moscow, 1959.

3.

DZERGACH, A. I. and RADZIN’SH, G. A. Calculation of the two-dimensional field of electromagnets with unsaturated iron by means of integral equations. Tr. Radiotekhn. in-ta Akad. Nauk SSSR,

polei na

No. 14,70-75,1973. U.S.S.R. Comput. MathsMath. Phys. Vol. 17, pp. 222-228 0 Pergamon Press Ltd. 1978. Printed in Great Britain.

A VERSION

OF THE MOMENT

THE TRANSFER

0041-5553/77/0601-0222$07.50/O

METHOD

OF SELECTIVE

OF CALCULATING RADIATION*

Yu. D. SHMYGLEVSKII Moscow (Received 18 March 1976) A MOMENT method using Laguerre frequency transfer for arbitrary optical cell dimensions calculation numerical

polynomials

is used to calculate the radiative

and large temperature

drops. The results of the

of radiative transfer in an isothermal half-space illustrate the accuracy of the scheme.

In [l]

the moment method is used for the transfer equation in a differential

form and

enables calculations on cells with small optimal dimensions to be performed at all frequencies. An improvement of the method was described in [2] . The results of [3,4] permitted the method to be extended [5] to the case of arbitrary optical dimensions of the cells. The version of [S] permits calculations to be performed for comparatively small temperature differences at adjacent computing points. This is due, for example, to the fact that when radiation from an absolutely black body at the temperature T2 falls on a particle at a temperature T, , it is necessary within a semiinfinite interval of measurement of the frequency v to expand in special polynomials the ratio of Planck functions B (Y, T,) /B (Y, Ti), which for T+T, increases exponentially with frequency. The stability of the method is proved in [6] , and the method itself was used to calculate the flow of air in a circular tube with transparent walls [7]. The method is further improved here. *Zh. vychisl. Mat. mat. Fir., 17, 3, 785-790,

1977.

223

Short communications

The equation of radiative energy transfer has the form E

B=_

= -x(l--B),

2hv3

1

c2

exp (hvlkT) - 1 ’

Here s is the linear dimension along the ray considered in space, Z is the intensity of radiation in the direction of increase of the linear dimension, x is the volume coefficient of absorption, depending, for example, on the pressure p and temperature T, h is Elan&s constant, c is the velocity of light, and k is Boltzmann’s constant. In the case of arbitrary optical dimensions of the computing mesh [4], just as in [5], we introduce the function l=Z-B. The transfer equation is transformed to the form dJ

b=-,

bT+3

-XL

--&-=

2kS

(1)

l&!Z,

h2c2

(ex-1)2 hv z=---, kT.

where T* is some fixed temperature. We consider one step of integration of the transfer equation from s to s+g. To obtain a scheme applicable to any optical steps s+E AZ =

s

x as,

in accordance with [4] the quantities X, 3T3/&, and also the quantity T occurring in x, are assumed to be constants. In this case Eq. (1) has the integral

J=Ye-En--

I~T bT2xb-,

l-e-cx x

(2)

where I=I(z, s+E), J +(z, s). Of the quantities occurring in (2), the variables J, J, x, p depend on z, and the remainder are constant. If the absorption coefficient x(z) contains lines or multiplets, then the background of the absorption coefficient, denoted by K(z), must initially be considered separately. The background is formed by cutting the lines and multiplets (see Fig. 1) by lines x=const so that K(z) will not be a sharply changing function. For example, in the frequency interval shown in Fig. 1, the graph of K (z) is the line 12345678. Equation (2) for the background has the form J=Te-EK-

8T I-e-En bT2-~_._.___-----_. a.S

K

This equation is subjected to instantaneous processing. The variables occurring in it are represented approximately in the form

224

Yu. D. Shmyglevskii M

J =

e-‘i”

M

c

y =

J”lLn(Z),

e-z/2

i,,=ll

Here the 1, (z) are orthonormed

Iaguerre

are known coefficients

The coefficients

JnJIn (2))

Tit=0

(5z)‘exp(+z+$)

the H,

c

-1]_1I(-l(r)

[l-e-~K(z)][exp(~z)

polynomials

[8] , the J,

are the unknown

coefficients,

1 exp (z/2) at the initial point of the interval.

of the expansion

Am and uqrn are calculated by the formulas

X [

exp

(+z)

-I]

~ZK-i(z)E,(z)dz,

m

psm =

s

e-z-EK(z)Zp (z)1, (2) dz.

0

FIG. I If T* and [ have fixed values, then the matrices X, and pQrn are calculated beforehand for those values ofp and T which permit interpolation

for performing

calculations

with

arbitrary p and T. Equation (3) is multiplied term by term by ~~(2)e-z/2 for q=O, 1, . . . , M, the expressions (4) are substituted in it and integration is performed with respect to z from 0 to 00. This leads to the equations

Equations (5) give the unknown quantities J4 as a result of the integration step. The quantities Jm on the boundary of the range in which the solution is sought, are determined from

(4)

225

Short communications

the known value of Jo(z),

corresponding

the formulas for the coefficients

to the radiation incident

of an expansion

from outside, by means of

in Laguerre polynomials

OD

J7, =

s

e-zJo(z)lm(z)dz.

0

To complete the succeeding steps it must be borne in mind that Jm at the new step equals Jm on the previous interval. The calculation

in the frequency

which the lines and multiplets coefficient,

is performed

zi,GzGzz,,.The

quantity

intervals of lines and multiplets,

were cut for the formation

that is, in intervals in

of the background

as in [5] . Each such interval numbered

of the absorption

n is defined by the inequalities

N varies from 1 to N, where N is the total number of intervals. The

intervals are assumed to be so small that with permissible accuracy the (aBIaT4 regarded as constant.

In accordance with [5] the quantity

(J), is represented

L

(J) n

=

z

r;=

Jndm (U,

22-z

may be

in the form

in-Zan

(6)

I

zzn-zin

77L=O

where the Pm(b) are Legendre polynomials.

The coefficients J,,, n are given by ’

2m+l J mn = 2

s

(J)nPmK)df.

--1

At an integration

step along the ray the quantities Jmn are determined

from the system

of equations

m-o

l=O, I,...,

The coefficients

L,

plrnn and L,

n=l,

2 ,...,

N.

are calculated by the formulas

2m+l p’mn = 2

The quantities

(7)

na=o

?W=O

i s -1

r~,,,*satisfy the recurrence

m-‘“I

n Pm(t) a,

relation

m+l rtm = m rl--l,m--i + rf-I.m+l, 2m+3 2m-1. rl,-t=O,

USSR

17-3-P

ram =

i, 1 0,

m==O, m*O,

E=l,2 ,...,

m=O,1,...

.

226

Yu. D. Shmyglevskii

In the same frequency intervals we calculate the quantity J formed because of the background of the absorption coefficient with constant values of the coefficient K,. This quantity is denoted by Jn* and for one step is determined by the equation

The rate of supply of heat due to radiation to unit volume in unit time is denoted by Q and is given by Q = -divq,

q=jj(~~dv)odw, &r 0

where q is the vector flux of radiative energy, o is the unit vector of solid angle, do is the differential of solid angle, and the sign 4n denotes integration over the whole solid angle. The quantity J is calculated by addition of the quantity J found at the abso~tion background K(z), to the quantities (J)n found for x(z) in all the intervals z~,,
J=~-z,2~l,h(.)+j;,~,J~~~~~~~-~Jn*. n-1 m-o m-o

n-1

On inte~ating the letter expression with respect to frequency it is necessary to remember the connections of z and 5 with v from (1) and (6). A remarkable equation exists [9] for the polynomials determined in [8] : 02

e-Z4,(z)

cl2 =

2,

m=O, 1,.

. .

.

J 0

Calculations give

(9)

The stability of the calculation by formulas (51, (7), (8) can be proved by almost literal repetition of 161. When the background of the abso~tion coefficient is used the problem is simplified in the limit. Indeed, it follows from the formula for plqrn, that this matrix is symmetric, and all its eigenvalues are positive. As an example a calculation was made of the radiative transfer in an isothermal half-space with temperature 10 OOO’K,on which falls isotropic radiation from an absolutely black body at a temperature of 100 000°K. The dependence of the absorption coefficient on frequency is given by the formula it= (IO5 k/hv) 2. Here ail the quantities are taken in the CGS system. The purpose of the calculations is to check the moment method for the background of the absorption coefficient; therefore the lines depending on x(z) are not included. Calculations for the absorption coefficient with the lines were carried out, for example, in [2,7]. The value of the integration step g was taken in two versions equal to 10 cm and 100 cm.

227

Short communications

TABLE 1

T$=io

I

T*

-7

M=l8

M=6

100

25000 12500 50000 100000 200000 400000

Y:

0.8 0.02 0.03 0.03 0.19

3 11 17

1600000 8°oooo

i-z 17’ 8

44;

The value of T, may be chosen depending

M=6

-

6250

s-i00

-_

-

100 100

733.1 1.0 5: 428 182

I I M=iO

M=iS

100

12

91

0.15

07 5 1:8

8’:: 0:35 3.0

: 242 385

578 20

on the problem considered.

In this case various

values of T* were chosen, given in the table, which contain the quantity &=I (Am--A) /A I iOO%, where A,

is the integral, calculated by Eq. (9) for calculating

the radiative transfer along the

normal to the surface of the halfspace, the quantity A is the exact value of the integral. The calculations

were performed

for various approximations

M.

E was calculated for T* = 100 000°K with steps .$= 10 cm for the interval

The same quantity

0
The calculations

were performed

The rapid monotonic

values of s were:

and kindly made available for this paper by I. N. Naimova.

decrease of X as z increases creates adverse conditions

for the use

of the moment method. It follows from the results given in the table that in practical calculations

situations

are

possible in which at different stages it is useful to use different values of T*. On passing from one value T*, to another T*, it is necessary to carry out a re-expansion

of the function J in

accordance with (4) by the formula

J 1-f m

(Jk)2

=

exp

M

(l+t)zz

0

where the subscript

1 applies to quantities

lc

(1,)

1 &(k?)lk(Z2)

a%,

t=-,

T l2 T -1

PI=0

obtained with T el, and the subscript 2 to quantities

obtained with Tc2. Translated by J. Berry

S. 0. Belotserkovskiiand V. A. Gushchin

228

REFERENCES 1.

SHMYGLEVSKII,Yu.D.~Icu~tion of radiative transfer by Galerkin’s method. Zh. v&h&Z.MN. mat. Fiz., 13,2,389-407,1973.

2. KRIVTSOV,V.M.,NAUMOVA,I.N.,SHULISHNINA,N.P. and SHMYGLEVSKII,Yu.D.Checkof two methods of calculating radiative transfer. Zh. vj%hisl.Mat. mat. Fiz., 15,1, 163-171,

1975.

3.

CHARAKHCH’YAN, A. A. A numerical scheme for the transfer equation on an optically coarse net. Zh. vychisl.Mat. mat. Fiz., 15,4,999-1005,1975.

4.

CHARAKHCH’YAN, A. A. An approach to the calculation of the transfer equation for problems of the dynamics of a radiating gas. In: Dynamics ofa radiatinggas (Dinamika izluchayushchego gaza), No. 2,16-35. VTs Akad. Nauk SSSR, Moscow, 1976.

5.

SHMYGLEVSKII, Yu. D. The moment method of calculating the transfer of selective radiation. In: Dynamics of a radiatinggas (Dinamika ~luchayushchego gaza), No. 2,42-60, VTs Akad. Nauk SSSR, Moscow, 1976.

6.

KRIVTSOV, V. M. Stability of the moment method of calculating the transfer of selective radiation. In: Dynamics ofa radiatinggas(Dinamika izluchayushchego caza). No. 2,61-62, VTs Akad. Nauk SSSR, Moscow, 1976.

7.

NAUMOVA, I. N. and SHMYGLEVSKII, Yu. D. Calculation of flows of radiating air in a tube. In: Dy~mi~s of (i radiat~nggas(Dinamika ~luc~yu~chego gaza), No. 2,99-108, VTs Akad. Nauk SSSR, Moscow, 1976.

8.

JACKSON, D. Fourier series and orthogonal polynomials (Ryady Fur’e i ortogonal’nye polinomy), Izd-vo in. lit., Moscow, 1948.

9.

DITKIN, V. A. and PRUDNIKOV, A. P. The operational calculus (Operatsionnoe shkola”, Moscow, 1975.

ischislenie), ‘“Vysshaya

U.S.S.R. Comput. Maths Math. Phys. Vol. 17, pp. 228-234 0 Pergamon Press Ltd. 1978. Printed in Great Britain.

NUMERICAL SIMULATION OF THE PLANE FLOW OF A VISCOUS FLUID UNDER THE ACTION OF AN EXTERNAL FORCE PERIODIC IN SPACE* S. 0. BELOTSERKOVSKII

and V. A. GUSIICHIN

Moscow (Received 2 1 June 1976) USING the splitting method, previously used to solve external flow problems, the motion of a viscous incompressible space is simulated. methodical

fluid in a rectangular

domain under the action of a force periodic in

The effect of the flow parameters

calculations

performed

on its structure is investigated.

The

enable the accuracy of the results to be judged.

1. The use of modern computers makes it possible to stimulate complex flows of a viscous incompressible fluid by solving the partial differential equations describing these phenomena. The essential non-linearity of these equations hampers the derivation of an analytic solution, except for special cases. At the same time it is precisely the non-linearity which frequently explains many interesting physical properties of this class of flows. The presence of a small parameter as the coefficient of the highest-order derivatives in the initial equations requires the development of high-precision numerical approaches, which are especially important in the field of boundary layer phenomena, for flows with large gradients etc. Therefore, the construction of highly accurate stable numerical algorithms for the non-stationary Navier-Stokes equations is an extremeIy urgent problem at the present time. The numerical approach complements, on the one hand, the results of physical experiments, and on the other, the data of analytic investigations. *Zh. vychisl.Mat. mat. Fiz., 17,3,791-795,

1977.