The nodal transport method for general triangular meshes in (X,Y)-geometry

Progress in Nuclear Energy. Vol. 18, No. 1/2, pp. 153-160, 1986. Printed in Great Britain. All rights reserve&

0149-1970/86 $0.00+.50 Copyright ~ 1986 Pergamon Journals Ltd.

THE NODAL TRANSPORT METHOD FOR GENERAL TRIANGULAR MESHES IN (X, Y)-GEOMETRY* RICHARD R. PATERNOSTERand WALLACE F. WALTERS Theoretical Applications Division, Los Alamos National Laboratory, Los Alamos, New Mexico 87545, U.S.A.

Abstract--A nodal differencing scheme is derived for two-dimensional (x,y)-geometry neutron transport on general triangular meshes. Two cases must be solved for triangular mesh cells; inflow through one face and inflow through two faces of the triangle. The one inflow side problem is solved by a linear characteristic method. The two inflow side problem is treated using a nodal solution. The linear characteristic-nodal (LCN) solution is compared with other two-dimensional differencing methods for spatial convergence and positivity. The LCN method was found to be third-order asymptotically convergent for the test cases examined.

1. INTRODUCTION The triangle is the most versatile geometric figure available for construction of two-dimensional computational meshes. Using simple subdivision it can form the basic mesh unit for complex polygonal grids. Design problems with curved boundaries or requiring precise geometric resolution can be modeled better with triangular mesh cells than with orthogonal meshes. The discrete ordinates spatial differencing methods currently available for solving general triangular mesh problems are diamond-difference-like schemes 1'2 and the linear discontinuous (LD) method based on finite element techniques. 3'4 The linear characteristic (LC) method 5 and the linear nodal (LN) transport method 6 are well developed spatial differencing schemes for solving two-dimensional rectangular mesh transport problems. Walters 7 has combined the LC and LN methods to solve equilateral triangle mesh problems. This technique has been successfully extended to the case of general triangles with the same methodology.8

(x3,%)

(0,01

,.~8,88

(x2,o~ x

/am Fig. la. The reference one inflow side triangle cell problem to be solved by the linear characleristic method.

Y

2. CHARACTERISTIC SOLUTION FOR ONE INFLOW SIDE TRIANGLES "- X

The general triangle problem is solved by considering two different cases, inflow through one side of the triangle (Fig. la), and inflow through two faces of the triangle (Fig. lb). The first case uses the analytic

* This work performed under the auspices of the U.S. Department of Energy.

Fig. lb. The reference two inflow sides triangle cell problem to be solved by the transport nodal method. 153

154

R.R. PATERNOSTERand W. F. WALTERS

solution of the discrete ordinates equation with the source approximated by a linear moment expansion about the centroid of the triangle

a~,

/*

1

a¢,

+Y3(~G-P,) G,

" a x + ' ~y + ~O(x, y) = sA + - ( x - xc)Sx

Xc

1

+ -(y-yc)Sy, Yc

(1)

where (xc, Yc) are the coordinates of the centroid. The source terms Sa, S~ and Sy are the zeroth source moment, first x-moment and first y-moment, respectively. The remaining symbols follow the usual nomenclature. Here the angular direction and energy group subscripts are omitted for clarity. The angular flux on the left (L) and right (R) outgoing faces is approximated by linear moment expansions

O,=O(x=y/m,, Y)= ~bL+ O L ( ~ - - 1)

O
(2) and

t,kR=@(x=y/mR, y)=~R+OR(2yy3 -- 1)

S

(5)

where PL = 1 --m L ~ q

and

PR = 1 - m R ~,

and

Po =

1

(I-e-9,

1

P , = - (gI - P o ) ,

and 1

P2 = - ( 1 - 2 P O. 8

Expressions (4) and (5) are the precise analytic expressions for the angular flux moments on the outgoing faces of the triangle given the linear source approximation within the triangle. The 0 coefficients are obtained from finite difference expressions of the angular flux at the vertices of the triangle.

0=
where ~ and 0 on the right side are the zeroth and first moment averages over the respective faces, and mL and mR are the slopes of the left and right faces of the triangle. The terms ~kL and tkR are obtained in terms of the knowns ~kB and 0B directly from the characteristic solution of (1) as

20nX3 Y3 tltL =(~n--OB)Po + PL(Po- PO + PaSA

3. NODAL SOLUTION FOR TWO INFLOW SIDE TRIANGLES The case with inflow through two sides of the triangle (Fig. lb) is solved by the linear nodal method. The coordinate system is first rotated to have the outflow face of the triangle parallel to the abscissa. The angular flux on the right (R) face is represented in the new (s, t) coordinate frame by a linear moment expansion

x2

(6)

-I- 2r/(x2 + X3) +Y3 /2

S

+ Y3 (3P2 _ P1)Sr,

(4)

where lR=s2--s 3 and ht=t2=t 3. The OR and OR coefficients are obtained from the zeroth and first spatial moments of the transformed discrete ordinates equation with a linear source approximation about the transformed centroid

u, 8¢9~+ ,80

and

2

,7 at +'1'(~' t ) = & + F~(~-~Os~

//X 3 - - X 2 \

3

+ ~ (t- tOs,.

(7)

+ Y~ P1 SA r/

The ~bR coefficient in terms of the known incoming terms ~ks, 0B, ~kL and 0L is given by

+

~bR=

3Y3

[p'(~O s -- 08) + (1 -- p')(~OL. --

+ [p'O 8 + (1 -- p')OL]2P ,

OL)]Po

General triangular meshes

+

--

h, +h'2 [ is 2+s3\ q, SA soP1

.',._

-

-

1870

q

t )"q

3

(8)

(t
!55

9 I

~etgenvoLue 865

where

NTWONTRAN (DD)

186C

P'--\m'L--rn'lt]\m' n q'l]'

[ 3 8 cm i,

a n d m~ a n d m~ are the slopes o f the b o t t o m a n d left faces in the t r a n s f o r m e d c o o r d i n a t e system. T h e 0 R t e r m is o b t a i n e d in a n a n a l a g o u s m a n n e r f r o m the first

1855

spatial m o m e n t of (7). C o m p l e t e details are given in

0

i

i

i

8cm I i I i [ i I ! I i I I L , I , I , I , 1 ~ I , 1J 02040608 I 0 12 14 16 1 8 2 0 2 2 2 4 2 6 2 . 8 3 0

reference 8.

Mesh

I J

spocing (cm)

Fig. 2. kef f eigenvalue for problem l (pure absorber). 4. N U M E R I C A L RESULTS T h e first test case is a c o m p a r i s o n w i t h the d i a m o n d difference spatial differencing s c h e m e as i m p l e m e n t e d in the T W O T R A N code. P r o b l e m 1 is a single material, high leakage, p u r e a b s o r b e r d i m e n s i o n e d 8 c m o n each side. T h e r e c t a n g u l a r T W O T R A N m e s h is split o n the d i a g o n a l s to f o r m triangles as s h o w n in the Fig. 2 insert. T h e L C N code, n a m e d C T R A N , used the s a m e q u a d r a t u r e set as T W O T R A N . T h e scalar fluxes a n d the eigenvalue are each c o n v e r g e d to 10 - 6

in b o t h codes. T h e results of the m e s h r e f i n e m e n t s t u d y for P r o b l e m 1 are s h o w n in T a b l e 1. T h e m e s h n o r m listed in T a b l e 1 is a m e a s u r e of the width. I n the case of the r e a c t a n g u l a r n o r m is the cell width. T h e m e s h n o r m cell was t a k e n to be length of the

m a x i m u m cell cells the m e s h for the triangle bisector. T h e

eigenvalues v e r s u s m e s h n o r m for S 4 is p l o t t e d in Fig. 2.

Table 1. Comparison of two-dimensional eigenvalues for problem 1 (pure absorber) CTRAN (LCN)

T W O T R A N (DD) Angular quadrature

Mesh Mesh"

norm

keff

Relative error

Mesh a

Mesh norm

kcff

Relative error

S2 S2 S2 S2 S2 S2 S2

4x4 2.0 8×8 1.0 16 x 16 0.5 32 x 32 0.25 64x64 0.125 128 x 128 0.0625 Infinite mesh

1.851104 1.858159 1.859846 1.860263 1.860368 1.860394 1.860404

9.30x 10 -3 2.25×10 -3 5.58 × 10 -4 1.41 × 10 -4 3.6 x 10 5 1.0 x 10- 5

4×2 8x4 16 x 8 32 × 16 64×32 128 × 64

4.472 2.236 1.118 0.559 0.280 0.140

1.862616 1.859164 1.860097 1.860365 1.860399 1.860403 1.860404

2.21X 10 - 3 1.24×10 -3 3.07 x 10 -4 3.90 x 10 -5 5.0×10 -0 1.0 x 10 6

S4 S4 S4 S4 S4 S4 S4

4x4 2.0 8x8 1.0 16×16 0.5 32 x 32 0.25 64x64 0.125 128 × 128 0.0625 Infinite mesh

1.859246 1.865236 1.866693 1.867054 1.867144 1.867166 1.867176

7.93×10 -3 1.94x10 -3 4.83 x 10 -4 1.22 x 10 -4 3.2 × 10 - s 1.0 x 10- s

4x2 8x4 16×8 32 × 16 64x32 128 × 64

4.472 2.236 1.118 0.559 0.280 0.140

1.872463 1.868327 1.867273 1.867175 1.867175 1.867176 1.867176

5.29 x 10 -3 1.15×10 a 9.7 x 10 s 1.0 × 10 -6 1.0× 10 6 0.0

S6 S6 S6 S6 S6 S6 S6

4x4 2.0 8×8 1.0 16x16 0.5 32 x 32 0.25 64×64 0.125 128 x 128 0.0625 Infinite mesh

1.858675 1.864607 1.866084 1.866450 1.866541 1.866565 1.866574

7 . 9 0 x 10 - 3

4x2 8 ×4 16×8 32 × 16 64×32 128 x 64

4.472 2.236 1.118 0.559 0.280 0.140

1 . 8 7 2 4 8 1 5.91 x 10 - 3 1.867997 1.42 × 10 - 3 1.866764 1.90 x 10 -'~ 1.866590 1.6 × 10 - 5 1.866575 1.0x 10 6 1.866574 0.0 1.867176

1.97 × 4.90 × 1.24 × 3.3 × 9.0 x

10 3 10 -4 10 -4 10 -5 10 6

A T W O T R A N mesh size shows number of rectangles, whereas CTRAN gives number of triangles (also in Tables 3 and 5).

156

R.R. PATERNOSTERand W. F. WALTERS

A c o m p a r i s o n of the eigenvalue versus mesh cell size for the d i a m o n d difference and linear characteristicnodal m e t h o d s shows the superior accuracy of the L C N solution. The eigenvalue relative error is also c o m p u t e d and given on Table 1. This error is plotted versus the mesh n o r m on Fig. 3 for P r o b l e m 1 with S 4 quadrature. It has been shown that the d i a m o n d difference m e t h o d is uniformly second-order convergent in two-dimensional (x,y)-geometry. 9 This is verified from the slope of the T W O T R A N convergence curve on Fig. 3. It can be seen that the L C N method, as implemented in this research, is asymptotically convergent. The linear characteristic-nodal results are used to calculate the asymptotic order of convergence, n<. The theoretical order of convergence, n,, for the d i a m o n d difference scheme is from Larsen, 9 and for the linear characteristic scheme on rectangular meshes from Lecot. lo These results are shown in Table 2. A second test problem was run to test the order of convergence of the linear characteristic-nodal m e t h o d with scattering included. The same square geometry of Problem 1 was used. The results for S 2 and S 4 quadratures are shown in Table 3. The convergence plot for these results is shown in Fig. 4 for S 4

Table 2. Convergence results for problem 1 (pure absorber TWOTRAN(DD)

Sz S4 S6

CTRAN(LCN)

nt

nc

nt

nc

2 2 2

1.972 1.926 1.956

3 3 3

2.974 3.304 3.790

102

8cm

P

I-f' r,,l,,ll

11,1o /I '1

//

: ,o-4 I-~:<,=to,,.'~:,:2o

J

f

/

p

f

/

keN

l

~> ~>

I0 ':

I0 "1

I00

I0 ;

Mesh norm {cm)

Fig. 3. Convergence plot for problem 1 (pure absorber).

quadrature. The asymptotic order of convergence for the L C N m e t h o d was again calculated from the slope of the convergence curve at the finest mesh spacing. The order of convergence results for Problem 2 are s h o w n in Table 4. The characteristic and transport nodal m e t h o d form a family of semi-analytic solutions depending on the order of the analytic functions chosen to represent the source t h r o u g h o u t the cell interior and the angular flux on the cell boundary. Lecot 1° has studied the theoretical spatial convergence properties of characteristic m e t h o d s on rectangular meshes in two-dimensional (x,y)-geometry. The family of characteristic solutions are denoted C k l , where k denotes the order of the polynomial a p p r o x i m a t i o n used for the interior source, and I denotes the order of the polynomial

Table 3. Comparison of two-dimensional eigenvalues for problem 2 (with scattering) TWOTRAN (DD) Angular quadrature S2 S2 S2

Mesh

Mesh norm

kcrr

CTRAN (LCN) Relative error

Mesh

Mesh norm

k©ff

Relative error

S2 S2 S2

4x4 2.0 8x 8 1.0 16 x 16 0.5 32 x 32 0.25 64 x 64 0.125 128 x 128 0.0625 Infinite mesh

0. 907893 0.919309 0.922069 0.922753 0.922924 0.922967 0.922985

1.51 x 10 -2 3.68 x 10- 3 9.16x10 4 2.32 x 10 -4 6.1 x 10 -5 1.8 x 10 -5

4x2 8x4 16x8 32×16 64 x 32 128 x 64

4.472 2.236 1.118 0.559 0.280 0.140

0.926621 0.920952 0.922478 0.922919 0.922975 0.922982 0.922985

3.64x 10 -3 2.03 x 10 - 3 5.07x10 -4 6.6x 10 -s 1.0x 10 -5 3.0x 10 6

S4 S4 S4 S4 S4 S4 S4

4x4 2.0 8x8 1.0 16 x 16 0.5 32 x 32 0.25 64 x 64 0.125 128 x 128 0.0625 Infinite mesh

0.921082 0.930961 0.933387 0.933989 0.934139 0.934177 0.934191

1.31 x 10 -z 3.23 x 10 -3 8.04x10 -4 2.02 x 10 -4 5.2x 10 -5 1.4x 10 -5

4x2 8x4 16x8 32 x 16 64x32 128 x64

4.472 2.236 1.118 0.559 0.280 0.140

0.943071 0.936113 0.934352 0.934188 0.934189 0.934190 0.934191

8.88 x 10 - 3

S2

1.92x 10 -3 1.61x10 4 3.0 x 10 -6 2.0x 10 -6 1.0x 10 -6

157

General triangular meshes

Hexagon test problem

I0'

(x,y)-geometry has been calculated by Lecot. 1° It is given by n, = Min (k + 2, 21 + 1). The theoretical order of convergence and the asymptotic order of convergence from the P r o b l e m 2 results are given in Table 6. The numerical order of convergence is consistently lower than the theoretical prediction. This is found to be consistent with results of L e c o d ° for the linear characteristic m e t h o d on two-dimensional rectangular meshes.

0 - THREETRAN(he~,z)-Se(90°)~ 102 ~° - TWOHEX-S4(60°) ~ / ? - CTRAN-S6(90°/ / 9

g # g

//

y

2/~mhe*ogon

I0 4

Table 6. Convergence results for problem 2 for COO,C01 and C10 characteristic approximations

B _

~o= 001 vT,, =005

[

_

I0

2,1=o~5,'Z~=ol4

I0 ~

IO

CO0

COl

CIO

IO"!

Hetght of mesh triangte {cm) S2 &

Fig. 4. Convergence plot for problem 2 (scattering problem).

g/t

nc

~lt

/lc

nt

t/c

1 1

0.932 0.937

2 2

1.895 1.835

1 1

0.831 0.822

Table 4. Convergence results for problem 2 (with scattering) TWOTRAN(DD)

S: &

CTRAN(LCN)

F/t

t/c

r/t

nc

2 2

1.942 1.974

3 3

2.836 3.170

a p p r o x i m a t i o n along the cell edges. The COO (step source, step edge flux), C10 (linear source, step edge flux), CI0 (linear source, step edge flux), and C I I (linear source, linear edge flux) a p p r o x i m a t i o n s were c o m p u t e d for Problem 2 for the S 2 and & angular quadratures. The results are presented in Table 5. The theoretical order of convergence of the various characteristic a p p r o x i m a t i o n s in two-dimensional

A series of test problems was also run to examine the properties of the linear characteristic-nodal m e t h o d for triangular meshes and to c o m p a r e results against other triangular mesh transport codes. Problem 3 is a hexagon subdivided into triangles. It was run on C T R A N and three other two-dimensional discrete ordinates transport codes; T H R E E T R A N (hex, z), T R I D E N T - C T R , and T W O H E X . T H R E E T R A N (hex, z) is a three-dimensional ( x , y , z ) - g e o m e t r y transport code using the DITR1 ( D i a m o n d TRiangle) spatial differencing scheme. 2 The T R I D E N T - C T R code is a two-dimensional banded triangle code employing the linear discontinuous (LD) differencing scheme. 4 T W O H E X is a new two-dimensional hexagon code using the TLC (Triangle Linear Characteristic) differencing scheme. 7 The TLC differencing

Table 5. Comparison of two-dimensional eigenvalues for various characteristic approximations problem 2 (with scattering) COO approximation Angular quadrature S2

Mesh Mesh

norm

C01 approximation

Relative keff

error

C10 approximation

Relative keff

error

keff

Relative error

4.472 2.236 1.118 0.559 0.280 0.0625

0.770290 0.826642 0.870110 0.895517 0.908993

1.53 × 10- 1 9.63 x 10- z 5.29×10 -z 2.75 × 10- z 1.39x10 2

0.776970 0.850790 0.898185 0.916038 0.921187

1.46 × 10- l 7.22 x 10 2 2.48×10 2 6.95 × 10- 3 1.80×10 3 . . .

0.879677 0.875785 0.887408 0.900779 0.910457 .

4.33 × 10- z 4.72 × 10 - 2 3.56×102 2.22 x 10 2 1.25×10 z

S2

4×2 8 ×4 16x8 32 × 16 64×32 128 × 64

& S4 & & S,~ S,,

4×2 8 ×4 16×8 32×16 64 × 32 128 x 64

4.472 2.236 1.118 0.559 0.280 0.0625

0.797764 0.846561 0.885988 0.909357 0.921667 . .

1.36x 10 ~ 8.76 x 10 -2 4.82×10 -2 2.48 x 10 2 1.25 × 10 2 . . . .

0.802356 0.866715 0.909676 0.926959 0.932261 . . .

1.32× 10 -~ 6.75 x 10 -2 2.45 × 10 -z 7.23×10 -3 1.93 × 10-3

0.901572 0.895572 0.903272 0.914619 0.923132

3.26× 10 z 3.86 × I 0 z 3.09×10 2 1.96×10 -1 1.11 × 10- z

Sz

S2 Sz Sz

158

R. R. PATERNOSTERand W. F. WALTERS Table 7. Comparison of two-dimensional eigenvalues for problem 3 (hexagon) THREETRAN(DD) S 6 (90°) Height of equilateral triangle

kef f

12.0 6.0 4.5 3.0 2.4 1.5 0.75 Infinite mesh

0.53983 0.58919 0.59796 0.60115 0.60265 0.60424 0.60501 0.605270

TRIDENT(LD)

Relative

Relative keff

error

6.54x10 -2 1.61x10 2 7.31 × 1 0 - 3 4.12x10 -3 2.62 x 10 -3 1.03x10 -3 2.60x 10 -4

0.593802 0.598321

5.98x10 3 1.46x10 -3

--

0.599510 -0.599694 0.599756 0.599780

, ,,

--

2.70x10 -4 -8.60x10 -5 2.4x 10 5

,

,o-+r tl [i l"'iJ"i',lgcm

/r,,l,,14J /

2

,o~

+6

,0

5

I-Zo=O25,,z,=o3o/ I

io-6 /

tO 2

EITWOTRAN (DO) /

I

I0 I Mesh norm (cm)

/ / LCN

I

IO0

POI

Fig. 5. Convergence plot for problem 3 (hexagon problem).

CTRAN(LCN) S 6 (90°)

Relative keff

error

scheme is a linear characteristic-nodal m e t h o d specificaily for equilateral triangles. The keff eigenvalue results are presented on Table 7. C T R A N a n d T H R E E T R A N (hex, z) used a n S 6 angular quadrature. The q u a d r a t u r e set used in b o t h codes is a 90 ° rotationally i n v a r i a n t from the T W O T R A N code. T R I D E N T - C T R was r u n With an $4 angular q u a d r a t u r e using a n EQNq u a d r a t u r e set. The T W O H E X code was r u n with a S 4 angular q u a d r a t u r e using a 60 ° rotationally i n v a r i a n t q u a d r a t u r e set. The use of different q u a d r a t u r e sets causes the eigenvalues from each code to converge to different values as the mesh size is reduced. However, as the n u m b e r of discrete angular directions increases these results would converge to the same value. The eigenvalue e r r o r results for the P r o b l e m 3 is s h o w n in Figs. 5 a n d 6. The L C N a n d the L D convergence results are clearly superior to the diam o n d difference-based T H R E E T R A N (hex, z) results. The C T R A N a n d the T W O H E X are b o t h linear characteristic-nodal schemes a n d produce parallel convergence results. The T R I D E N T - C T R linear dis-

/

TWOHEX(LCN) S, (60°)

S, (EQN)

0.62363 0.60111 0.59950 0.59912 0.59900 0.59891 0.59890 0.59889

I0 2

~~"

~o-3

0.631378 0.607884 -0.605561 0.605407 0.605297 0.605275 0.605270

2.61x10 -2 2 . 6 1 x l 0 -3 -2.91x10 -4 1.37 x 10 -4 2.7x10 -5 5.0x 1 0 - 6

Hexagontest problem

to-,

w=

2,47x10 -2 2,22x10 -3 6,10 x 10 -4 2.30x10 4 1.10 x 10 -4 2.00x10 -5 1.0x 10 -5

Relative error

keff

error

+ - TRIDENT-CTR- S4(EON) / -~ CTRAN- S6 (90°) /j+

f/ i

24 cm Hexagon

_~ s ~o IO 5

/

\/VV 2 a = O.OI#Y~ = 0 0 3

~o-6IO I

I ~t I= 0"151 E~ =0 14 i0o i0 j 102 Height of mesh tnangte (cm)

Fig. 6. Convergence plot for problem 3 (hexagon problem).

c o n t i n u o u s scheme appears to give more accurate results o n courser meshes, but becomes less accurate o n finer meshes. It is also observed t h a t the linear characteristic-nodal m e t h o d is almost uniformly convergent o n equilateral triangle meshes. The calculated order of convergence for the hexagon test p r o b l e m is: 1.993 for T H R E E T R A N (hex, z), 1.990 for TRID E N T - C T R , 3.397 for T W O H E X , a n d 3.296 for CTRAN. The final p r o b l e m is a c o m p a r i s o n of the positivity of the linear characteristic-nodal scheme a n d the diamond-difference method. The p r o b l e m is a severe test of the positivity of a spatial difference scheme. Region 1 is a highly scattering material with a uniform volume fission source. Region 1 s u r r o u n d s a highly a b s o r b i n g central Region 2. N o negative flux fixups were used in either calculation. The scalar flux across the p r o b l e m is plotted in Fig. 7. The results indicate that the linear characteristic-nodal m e t h o d is m u c h more positive t h a n the d i a m o n d difference scheme.

159

General triangular meshes Positivity test probtem

30 25

20

I

>,

tO o 8 co

-x

z~- LCN Tnangutar mesh

05

[3- T W O T R A N

0 -05

I

~,-~r-

~- . . . .

J',

=-~-~'

i ,,

--

i

Reg,~)r, 2

F~egio~ ,

to Or

02

I I 0 3 04.

05

06

i ~ecji~13

I 07

I 08

I 09

I0

II

12

X - coordinate (cml

Fig. 7. Scalar flux across problem 4 (positivity testl.

5. SUMMARY AND CONCLUSIONS This study has shown that the linear characteristicnodal semianalytic numerical solution can be extended to general triangular meshes. The one inflow side arbitrary triangle has been solved by the linear characteristic method. This method uses a linear polynomial to represent the angular fluxes on the triangle faces. The zeroth and first spatial moments of the angular flux are calculated from the analytic solution of the discrete ordinates equation. Vertex source values are calculated and stored for each triangle. The two inflow side arbitrary triangle has been solved by the transport nodal method. In this method the zeroth and first spatial moment equations of the discrete ordinates equations are solved to obtain face-centered values of the zeroth and first moments of the angular flux on the outflow side of the triangle. This study demonstrates the superior accuracy of the linear characteristic nodal spatial differencing scheme over rectangular and triangular diamond difference (DD) schemes. It has been shown that the linear characteristic-nodal method for arbitrary triangle meshes is at least asymptotically third-order convergent for eigenvalue calculations. For rectangular meshes the order of convergence appeared to approach fourth-order for very fine meshes. The COO, C01, and C10 characteristic approximations were examined for rectangular mesh cells. It was found that the calculated order of convergence of the arbitrary triangle solution asymptotically approached the theoretical order of convergence predicted by Lecot. ~o However, these solutions are of such poor accuracy that they can be of little value in any two-dimensional transport calculation. It would appear from the COO, COl, and the C10 results that the representation of the

cell edge fluxes is of greater importance than the representation cell interior source. The C01 (step source, linear edge flux) demonstrated a higher order of convergence and a more uniform convergence behavior than the CO0 or C10 approximations. The LCN solution for the arbitrary triangle was examined for positivity on a severe cross section gradient problem. The cell average and cell edge fluxes did become negative in the 'black' interior region, but their magnitude was four orders less than the diamond difference fluxes. The comparison with two-dimensional triangular mesh transport codes proved interesting, but frustrating due to the variety of angular quadratures used in these codes. The mesh refinement technique proved useful here, as the convergence properties of each code could be compared irrespective of infinite mesh eigenvalue result. The arbitrary triangle and the equilateral triangle linear characteristic-nodal solutions generated almost identical results. The calculated order of convergence for eigenvalues in both cases was considerably greater than three. It would appear that there is some enhanced spatial convergence from the equilateral triangle cell symmetry. The asymptotic order of convergence of the linear discontinuous (LD) spatial differencing scheme was found to be approximately the same as the linear characteristic-nodal differencing scheme. The accuracy of LD was much greater than DD, and for course meshes, even greater than the linear characteristicnodal scheme. The linear characteristic eigenvalues became more accurate at mesh triangle heights less than 3 cm. This effect is possibly due to the first-order accurate finite difference schemes used to calculate the first moment term in the one inflow side linear characteristic solution and the zeroth moment contributions in the rotated coordinate frame of the first moment term in the two inflow side nodal solution. A calculation of these terms using the balance equation would result in greater accuracy on the courser meshes.

REFERENCES

1. Reed W. H. (1971)Triangular mesh differenceschemes for the transport equation. LA-4769, Los Alamos Scientific Laboratory. 2. Waiters W. F. and O'Dell R. D. (1979) A threedimensional hexagonal-z difference scheme for discrete ordinates codes, Proceedings q[the ANS Topical Meetin9 on Computational Methods in Nuclear Engineering, Williamsburg, Virginia, April 23 25. 3. Seed T. J., Miller W. F. and Brinkley F. W. (1977) TRIDENT: A two-dimensional, multigroup, triangular mesh, discrete ordinates, explicit neutron transport code. LA-6735-M, Los Alamos ScientificLaboratory.

160

R.R. PATERNOSTERand W. F. WALTERS

4. Seed T. J. (1979) TRIDENT-CTR user's manual. LA7835-M, Los Alamos Scientific Laboratory. 5. Larsen E. W. and Alcouffe R. E. (1981) The linear characteristic method for spatially discretizing the discrete ordinates equations in (x,y)-geometry, Proceedings of the ANS-ENS International Topical Meetin9 on Advances in Mathematical Methods for the Solution of Nuclear Engineering Problems, Munich, Germany, April 27 29. 6. Walters W. F. and O'Dell R. D. (1981) Nodal methods for the dicrete-ordinates transport problems in (x,y)-geometry, Proceedings of the ANS-ENS International Topical Meeting on Advances in Mathematical Methods for the Solution of Nuclear Engineering Problems, Munich, Germany, April 27--29. 7. Walters W. F. (1983) The TLC scheme for numerical solution of the transport equation on equilateral triangle meshes, Proceedings of the ANS Topical Meeting on

Advances in Reactor Computations, Salt Lake City, Utah, March 28-31. 8. Paternoster R. R. (1984) A linear characteristic-nodal transport method for the two-dimensional (x,y)-geometry multigroup discrete ordinates equations over an arbitrary triangular mesh. LA-10092-T, Los Alamos National Laboratory. 9. Larsen E. W. and Nelson P. (1982) Finite-difference approximation and superconvergence for the discreteordinate equations in slab geometry, S l A M J. Numer. Anal. 19, 334. 10. Lecot C. (1983) Spatial convergence properties of some characteristic methods for the SN equations in (x,y)geometry, Proceedings of the ANS Topical Meeting on Advances in Reactor Computations, Salt Lake City, Utah, March 28-31.