A A New An Analytical S N So Solut ution i n in Sl n Slab Ge b - - PowerPoint PPT Presentation

A A New An Analytical SN So Solut ution i n in Sl n Slab Ge b Geome metry y

Dean Wang, Tseelmaa Byambaakhuu

University of Massachusetts Lowell November 1, 2017

2017 ANS Winter Meeting, Washington DC

Why another solution?

• Previous work:
• Chandrasekhar 1960; Vargas 1997; Warsa 2002; Ganapol

2008; Goncalez 2011, …

• Solution methods: Separation of variables, Green’s function,

Laplace transfer, and decomposition method.

• Our solution techniques:
• Eigen decomposition: a system of coupled SN PDEs is

decoupled into a system of separate ODEs.

• Boundary treatment: the left and right incoming angular flux

vectors are combined into one single vector.

• Derivation: the whole derivation process is based on linear

algebra.

• Solution: a truly closed-form analytical expression.

2

Problem Statement

Find the solution of the monoenergetic SN equation in slab geometry: 𝛎 𝑒 𝑒𝑦 𝛀 + Σ'𝛀 = Σ) 2 𝐗𝛀 + 𝑅 2 1

𝛀 = 𝜔/ 𝜔0 … 𝜔2 𝑼, angular flux vector; 𝛎 = 𝝂 −𝝂 , 𝑂×𝑂 matrix consisting of Gauss-Legendre quadrature direction cosine values, and 𝝂 = diag(𝜈>) > 0, 𝑜 = 1, … ,

2

𝐗 = 𝒙 𝒙 𝒙 𝒙 , 𝑂×𝑂 matrix consisting of Gauss-Legendre quadrature weights, and in which 𝒙 = 𝑥/ 𝑥0 … 𝑥F

G

𝑥/ 𝑥0 … 𝑥F

G

⋮ ⋮ ⋱ 𝑥F

G

𝑥/ 𝑥0 … 𝑥F

G

,

2 0 × 2 0 matrix, and ∑>K/

F G

𝑥> = 1; 𝟐 = 1 1 … 1 𝑼; Σ', total macroscopic cross section; ΣM, macroscopic scattering cross section; 𝑅, constant neutron source.

where L

3

Solution

𝑒 𝑒𝑦 𝛀 + Σ'𝛎N/ 𝐉 − c 2 𝐗 𝛀 = 𝐫

𝑑 = ST

SU , scattering ratio

𝐫 = V

0 𝛎N𝟐𝟐

where Σ'𝛎N/ 𝐉 − c 2 𝐗 = 𝐒𝚳𝐒N/ Matrix eigen decomposition:

𝚳 = 𝚳Y 𝚳N , and in which 𝚳Y = diag(𝜇>), 𝑜 = 1, …

2 0; and

𝚳N = diag(𝜇>), 𝑜 = 2

0 , … 𝑂

where

4

Solution

𝑒 𝑒𝑦 𝐒N/𝛀 + 𝚳𝐒N/𝛀 = 𝐒N/𝐫 Let 𝕑 = 𝑧/ 𝑧0 ⋮ 𝑧2 = 𝐒N𝟐𝛀, and 𝐜 = 𝐒N𝟐𝐫, we have 𝑒 𝑒𝑦 𝕑 + 𝚳𝕑 = 𝐜 Integrating gives the analytical solution: 𝕑 = 𝚳N/𝐜 − eN_𝚳𝒃 where 𝒃 = 𝑏/ 𝑏0 … 𝑏2 𝑼

5

Solution

where 𝒃Y 𝒃N can be determined by the boundary conditions at 𝑦 = 0 and 𝑀: 𝕑Y 𝕑N = 𝚳Y

N/𝐜Y − eN_𝚳d𝒃Y

𝚳N

N/𝐜N − eN_𝚳e𝒃N

𝒃Y = 𝚳Y

N/𝐜Y − 𝕑Y f ,

𝑦 = 0 𝒃N = eg𝚳e𝚳N

N/𝐜N − eg𝚳e𝕑N h ,

𝑦 = 𝑀,

𝛀Y

𝟏

𝛀N

𝐌 = 𝐉

𝟏 𝐒 𝕑Y

f

𝕑N

f + 𝟏

𝐉 𝐒 𝕑Y

h

𝕑N

h

After some algebra: where 𝕑Y

f

𝕑N

h can be determined by the following equation:

𝕑Y

f

𝕑N

h

= 𝐒𝟐𝟐 𝐒𝟐𝟑eg𝚳e 𝐒𝟑𝟐eNg𝚳d 𝐒𝟑𝟑

N/ 𝛀Y 𝟏

𝛀N

𝐌 −

𝐒𝟐𝟐 𝐒𝟐𝟑eg𝚳e 𝐒𝟑𝟐eNg𝚳d 𝐒𝟑𝟑

N/

× 𝐒𝟐𝟑 𝐉 − eg𝚳e 𝐒𝟑𝟐 𝐉 − eNg𝚳d 𝚳Y

N/𝐜Y

𝚳N

N/𝐜N

6

Solution

𝛀 = 𝛀Y 𝛀N = 𝐒 𝕑Y 𝕑N = 𝐒 l 𝐉 − eN_𝚳d e gN_ 𝚳e 𝚳N𝟐𝐒N𝟐

V 0 𝛎N𝟐𝟐

m + eN_𝚳d e gN_ 𝚳e 𝕑Y

f

𝕑N

h

𝛎N𝟐𝟐 = 𝐒𝚳𝐒N𝟐 𝟐

SU 𝐉 − n 0 𝐗 N/

𝟐 = 𝐒𝚳𝐒N𝟐

/ SU /No 𝟐

7

Final Solution

Φ = 𝕏𝑼𝛀 =

V SU /No − 𝕏r𝐒 eN_𝚳d

e gN_ 𝚳e 𝐒N𝟐×

V / SU /No 𝟐 − 𝐒 𝕑Y f

𝕑N

h

Remark:

• Diffusion limit: 𝛀 ≈

V / SU /No 𝟐 = t 0 𝟐,

as Σ' → ∞

• Thin limit:

𝛀 ≈ 𝐒 𝕑Y

f

𝕑N

h

= 𝛀Y

𝟏

𝛀N

𝑴 ,

as Σ' → 0 𝛀 =

V / SU /No 𝟐 − 𝐒 eN_𝚳d

e gN_ 𝚳e 𝐒N𝟐×

V / SU /No 𝟐 − 𝐒 𝕑Y f

𝕑N

h

Particular Solution Homogenous Solution

8

Eigen Decomposition

Conditioning of Eigenvalues: Cond 𝜇 = 𝑣 0 𝑥 0 𝑣, 𝑥 ~1 where 𝑣 and 𝑥 are the right and left eigenvectors associated with 𝜇. 𝐵 ≡ Σ'𝛎N/ 𝐉 − c 2 𝐗 = 𝐒𝚳𝐒N/ Conditioning of Eigenvectors: Cond 𝑣 = 𝑇 𝜇 𝐽 − 𝑄 where 𝑇 𝜇 is the reduced resolvent of 𝐵 at 𝜇, and 𝑄 is the spectral projector associated with 𝜇. for Matlab “eig” function

9

Saad 2011

SN Angular Convergence

1.0E-15 1.0E-13 1.0E-11 1.0E-09 1.0E-07 1.0E-05 1.0E-03 1.0E-01 1 10 100 1000 10000 100000 L1 Error N c = 0 c = 0.4 c = 0.8 c = 0.99

Σ• = 1 cmN/ and L = 1 cm

1.0E-11 1.0E-10 1.0E-09 1.0E-08 1.0E-07 1.0E-06 1.0E-05 1.0E-04 1.0E-03 1.0E-02 1.0E-01 1 10 100 1000 10000 100000 L1 Error N sigma_t = 1 cm^-1 sigma_t = 5 cm^-1 sigma_t = 10 cm^-1

𝑑 = 0.8 and L = 1 cm

1.0E-11 1.0E-10 1.0E-09 1.0E-08 1.0E-07 1.0E-06 1.0E-05 1.0E-04 1.0E-03 1.0E-02 1.0E-01 1 10 100 1000 10000 100000 L1 Error N L = 1 cm L = 5 cm L = 10 cm

Σ• = 1 cmN/ and c = 0.8

𝜁 = 0.98𝑂N0./ˆ

Gauss-Legendre Quadrature

10

Reference: 𝑂 = 2/‰ = 16394

SC Spatial Error – S10 1D Slab

1.0E-16 1.0E-14 1.0E-12 1.0E-10 1.0E-08 1.0E-06 1.0E-04 1.0E-02 0.000001 0.00001 0.0001 0.001 0.01 0.1 L1 Error Mesh Size (cm) c = 0.8 c = 0.2 c = 0.001 c = 0

~𝑃 ℎN/ ~𝑃 ℎ0 10 cm Σ' = 2 cmN/, 𝑅 = 1 cmN•sN/

11

Inhomogeneous Case

• 0.4

0.4 0.8 1.2 10 20 30 40 50 60 70 80 Scalar Flux Mesh Points Analytical DD SC

Σ' = 50 cmN/ 𝑑 = 0.6 𝑀2 = 4 cm ℎ = 0.1 cm 𝑅 = 1 cmN•sN/ Σ' = 2 cmN/ 𝑑 = 0.6 𝑀1 = 2 cm ℎ = 0.1 cm 𝑅 = 1 cmN•sN/ Σ' = 2 cmN/ 𝑑 = 0.6 𝑀3 = 2 cm ℎ = 0.1 cm 𝑅 = 1 cmN•sN/ 12