Solving the Radiative Transfer Equation using the Radon Transform - - PowerPoint PPT Presentation

1 Motivation The Radon Transform (RT) Hyperbolic PDEs Results

Solving the Radiative Transfer Equation using the Radon Transform

Megan Oeltjenbruns (Wayne State College)

Collaborators: Eappen Nelluvelil (Rice University), Jacob Spainhour (Florida State University), Christine Vaughan (Iowa State University), and Dr. James Rossmanith (Iowa State University)

NCUWM 2020

2 Motivation The Radon Transform (RT) Hyperbolic PDEs Results

Outline

1

Motivation

2

The Radon Transform (RT) RT Derivation in 2D Discretizing the RT Computing the IRT and Backprojection

3

Hyperbolic PDEs The Radiative Transfer Problem RTs of Derivatives Discretization of the PDE

4

Results Results

3 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Goal

Imaging Problem

The imaging problem is to reconstruct a distribution from several profiles taken at different angles The Radon transform is a simple way to produce profiles Certain advection problems are easier to solve profile-by-profile

4 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Goal

Our Goal

To efficiently solve multi-dimensional systems of hyperbolic partial differential equations using the Radon transform hi Three main steps:

1

Forward Radon transform multi-dimensional problem

2

Solve family of 1D advection problems

3

Use inverse Radon transform to original space

5 Motivation The Radon Transform (RT) Hyperbolic PDEs Results RT Derivation in 2D

Outline

1

Motivation

2

The Radon Transform (RT) RT Derivation in 2D Discretizing the RT Computing the IRT and Backprojection

3

Hyperbolic PDEs The Radiative Transfer Problem RTs of Derivatives Discretization of the PDE

4

Results Results

6 Motivation The Radon Transform (RT) Hyperbolic PDEs Results RT Derivation in 2D

Coordinate Rotation

ω s z y x

7 Motivation The Radon Transform (RT) Hyperbolic PDEs Results RT Derivation in 2D

Radon Transform

Suppose f : R2 → R is a function with compact support The Radon transform of f is formally defined as follows R( f) = f(s, ω) := ∞

−∞

f( x(s, z; ω), y(s, z; ω) ) dz = ∞

−∞

f(s cos(ω) − z sin(ω), s sin(ω) + z cos(ω)) dz

8 Motivation The Radon Transform (RT) Hyperbolic PDEs Results RT Derivation in 2D

Radon Transform at (x0, y0)

(x0, y0) ω s z

9 Motivation The Radon Transform (RT) Hyperbolic PDEs Results RT Derivation in 2D

Radon Transform

(s1, ω) (s2, ω) (s3, ω) (s4, ω)

ω s z

10 Motivation The Radon Transform (RT) Hyperbolic PDEs Results RT Derivation in 2D

Radon Transform Example

−1.0 −0.5 0.0 0.5 1.0 x −1.0 −0.5 0.0 0.5 1.0 y

f

0.00 0.12 0.24 0.36 0.48 0.60 0.72 0.84 0.96

−1.0 −0.5 0.0 0.5 1.0 s

π 4 π 2 3π 4

π ω

• f

0.00 0.09 0.18 0.27 0.36 0.45 0.54 0.63 0.72

11 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Discretizing the RT

Outline

1

Motivation

2

The Radon Transform (RT) RT Derivation in 2D Discretizing the RT Computing the IRT and Backprojection

3

Hyperbolic PDEs The Radiative Transfer Problem RTs of Derivatives Discretization of the PDE

4

Results Results

12 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Discretizing the RT

Grid Structure

Need discretized domain for inverse computation Profile-by-profile construction suggests Nω evenly-spaced diameters, with ω held constant

Chose a method that allows for spectrally accurate interpolation, differentiation

13 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Discretizing the RT

Grid Structure Example

−1.0 −0.5 0.5 1.0 x −1.0 −0.5 0.5 1.0 y

One diameter of the mesh

14 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Discretizing the RT

Grid Structure Example

−1.0 −0.5 0.5 1.0 x −1.0 −0.5 0.5 1.0 y

Full mesh

15 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Discretizing the RT

Computing Integrals

Compact support of f means we compute line integrals across chords of the domain Compute line integrals with Clenshaw-Curtis quadrature R( f) = ∞

−∞

f(x, y) dz ≈

Nq

• i=1

wi f(zi)

16 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Discretizing the RT

Quadrature Example

−1 1 x −1 1 y

17 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Discretizing the RT

Quadrature Example

−1 1 x −1 1 y

(s0, ω0)

18 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Discretizing the RT

Quadrature Example

−1 1 x −1 1 y

19 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Discretizing the RT

Quadrature Example

−1 1 x −1 1 y

20 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Discretizing the RT

Quadrature Example

−1 1 x −1 1 y

21 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Discretizing the RT

Interpolation Scheme

Need to sample f at arbitrary quadrature nodes, use interpolation Use a series of one-dimensional interpolation schemes:

Spectrally accurate on each diameter Fourth order accurate across angles

22 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Discretizing the RT

Interpolation Scheme Example

23 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Discretizing the RT

Interpolation Scheme Example (cont.)

24 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Discretizing the RT

Interpolation Scheme Example (cont.)

25 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Discretizing the RT

Interpolation Scheme Example (cont.)

26 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Computing the IRT and Backprojection

Outline

1

Motivation

2

The Radon Transform (RT) RT Derivation in 2D Discretizing the RT Computing the IRT and Backprojection

3

Hyperbolic PDEs The Radiative Transfer Problem RTs of Derivatives Discretization of the PDE

4

Results Results

27 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Computing the IRT and Backprojection

Ill-posedness of the IRT

Suppose f(x, y) : R2 → R is a function with compact support, and R( f) = f In discretized case, we only know f, f at mesh points We can approximate R( f) with a matrix-vector multiplication R f = f

28 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Computing the IRT and Backprojection

Backprojection

The adjoint of the Radon transform, denoted R∗, is known as backprojection If f (s, ω) : R × [0, π] → R, the backprojection of f is R∗

• f
• (x, y) :=

π f(x cos(ω) + y sin(ω), ω) dω

29 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Computing the IRT and Backprojection

Graphical Representation of Backprojection

−1 1

x

−1 1

y

(x, y)

30 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Computing the IRT and Backprojection

Using BiCGSTAB to Solve the Normal Equations

We can approximate the action of R∗ using RT We want to solve RT R f = RT f We do not explicitly need RT R, just the matrix-vector products R f, RT R f

• We use an iterative method (BiCGSTAB) to solve

RT R f = RT f

31 Motivation The Radon Transform (RT) Hyperbolic PDEs Results The Radiative Transfer Problem

Outline

1

Motivation

2

The Radon Transform (RT) RT Derivation in 2D Discretizing the RT Computing the IRT and Backprojection

3

Hyperbolic PDEs The Radiative Transfer Problem RTs of Derivatives Discretization of the PDE

4

Results Results

32 Motivation The Radon Transform (RT) Hyperbolic PDEs Results The Radiative Transfer Problem

Application: Radiative Transfer

We want to solve the Radiative Transfer Equation F,t + Ω · ∇F + σtF = σs 4π

• S2 F dΩ

A kinetic model for subatomic particles propagating through a homogeneous medium

Ω · ∇F is a transport term

σs 4π

• S2 F dΩ − σt F is a collision term

A linear transport equation in 1 + 5 dimensions

33 Motivation The Radon Transform (RT) Hyperbolic PDEs Results The Radiative Transfer Problem

Spherical Harmonics

Use spherical harmonics to create PN equations F(t, x, Ω) ≈

N

• ℓ=0

ℓ

• m=−ℓ

Fm

ℓ (t, x)Ym ℓ (µ, φ)

Already know spherical harmonics, Ym

ℓ (µ, φ); removes

angular dependence from F Need to solve for (most) of the Fm

ℓ (t, x)

34 Motivation The Radon Transform (RT) Hyperbolic PDEs Results RTs of Derivatives

Outline

1

Motivation

2

The Radon Transform (RT) RT Derivation in 2D Discretizing the RT Computing the IRT and Backprojection

3

Hyperbolic PDEs The Radiative Transfer Problem RTs of Derivatives Discretization of the PDE

4

Results Results

35 Motivation The Radon Transform (RT) Hyperbolic PDEs Results RTs of Derivatives

Solving the PN Equations

We now consider the following system of M differential equations q,t + A q,x + B q,y = C q q(t = 0, x, y) = q0(x, y)

36 Motivation The Radon Transform (RT) Hyperbolic PDEs Results RTs of Derivatives

Radon Transforms of Partial Derivatives

The Radon transform is a spatial transformation, thus: R ( f,t) = ∂ ∂tR( f) = f,t Therefore, R( f,t) = f,t R( f,x) = cos(ω) f,s R( f,y) = sin(ω) f,s All spatial derivatives in xy space become derivatives of s Note: R is a linear operator

37 Motivation The Radon Transform (RT) Hyperbolic PDEs Results RTs of Derivatives

Radon Transform of the PDE

We can now complete the transformation and add like terms: R

• q,t + A q,x + B q,y = C q
• =

⇒ R

• q,t
• + A R
• q,x
• + B R
• q,y
• = C R
• q
• =

⇒ q,t + cos(ω) A q,s + sin(ω) B q,s = C q = ⇒ q,t +

• cos(ω) A + sin(ω) B
• q,s = C

q

38 Motivation The Radon Transform (RT) Hyperbolic PDEs Results RTs of Derivatives

Radon Transform of the PDE

Given

• q,t +
• cos(ω) A + sin(ω) B
• q,s = C

q let A(ω) := cos(ω) A + sin(ω) B so that

• q,t +

A(ω) q,s = C q

• q(t = 0, s, ω) =

q0(s, ω) For each ω, we now have a system of one dimensional PDEs which can be solved with traditional methods

39 Motivation The Radon Transform (RT) Hyperbolic PDEs Results RTs of Derivatives

Defining Hyperbolicity

We say the PDE q,t + A q,x + B q,y = C q is hyperbolic if

• A(α) := cos(α) A + sin(α) B

is diagonalizable with only real eigenvalues for all α ∈ R Physically, this means that information in the system travels at a finite speed The wave equation and PN equations are hyperbolic

40 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Discretization of the PDE

Outline

1

Motivation

2

The Radon Transform (RT) RT Derivation in 2D Discretizing the RT Computing the IRT and Backprojection

3

Hyperbolic PDEs The Radiative Transfer Problem RTs of Derivatives Discretization of the PDE

4

Results Results

41 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Discretization of the PDE

Advection Example

2 4 6

Physical Solution

t0 tf

2 4

Characteristic Solution

−1.0 −0.5 0.0 0.5 1.0 −2 2 −1.0 −0.5 0.0 0.5 1.0 −4 −2

42 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Discretization of the PDE

Time-stepping Methods

We discretize in space and have M semi-discrete equations Wp,t + λpDp Wp =

M

• q=1

Fpq Wq More freedom in selecting a time-stepping method We use a third-order method

43 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Results

Outline

1

Motivation

2

The Radon Transform (RT) RT Derivation in 2D Discretizing the RT Computing the IRT and Backprojection

3

Hyperbolic PDEs The Radiative Transfer Problem RTs of Derivatives Discretization of the PDE

4

Results Results

44 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Results

Full Example: Wave Equation

x y

f0

R →

s ω

• f0

↓ solve 1D equations

x y

fT

← R−1

s ω

• fT

45 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Results

Full Example: P1 Approximation

−1.5 −0.75 0.0 0.75 1.5 x −1.5 −0.75 0.0 0.75 1.5 y

P1 at t = 0

10 20 30 40 50 60 70 80 90

−1.5 −0.75 0.0 0.75 1.5 s

π 4 π 2 3π 4

π ω

P1 at t = 0 in Radon space

0.00 0.75 1.50 2.25 3.00 3.75 4.50 5.25 6.00 −1.5 −0.75 0.0 0.75 1.5 x −1.5 −0.75 0.0 0.75 1.5 y

P1 at t = 1.0

−1.50 −0.75 0.00 0.75 1.50 2.25 3.00 3.75 4.50 −1.5 −0.75 0.0 0.75 1.5 s

π 4 π 2 3π 4

π ω

P1 at t = 1 in Radon space

0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8

46 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Results

Acknowledgments

NSF Grant DMS-1457443 Iowa State University

• Dr. Rossmanith and Christine Vaughan

Eappen and Jacob

47 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Results

Thank You!

Future Work

Adding spatially-dependent collision terms to the PN equations Implementing more sophisticated/higher order timestepping schemes Improving efficiency through a parallelization of transport computations hi

47 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Results

Thank You!

Future Work

Adding spatially-dependent collision terms to the PN equations Implementing more sophisticated/higher order timestepping schemes Improving efficiency through a parallelization of transport computations hi

Questions?

48 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Results

References

(1) Brunner, T. & Holloway, J. (2005). Two-dimensional time dependent Riemann solvers for neutron transport, Journal of Computational Physics 210 (2005) 386-399. hi (2) Peterson, L. (2018). An asymptotic-preserving spectral method based on the radon transform for the PN approximation

• f radiative transfer, Master’s thesis, Iowa State University,

2018. hi (3) Pieraccini, S. & Puppo, G. (2007). Implicit-Explicit Schemes for BGK Kinetic Equations, Journal of Scientific Computing, Vol. 32, No. 1, July 2007.

49 Motivation The Radon Transform (RT) Hyperbolic PDEs Results Results

References (cont.)

(4) Press, W. (2006). Discrete Radon transform has an exact, fast inverse and generalizes to operations other than sums along lines, PNAS December 19, 2006 103 (51) 19249-19254 hi (5) Rim, D. (2018). Dimensional splitting of hyperbolic partial differential equations using the Radon transform, SIAM J. Sci. Comput., 40(6) (2018), A4184-A4207 hi (6) Shin, M. (2019). Hybrid discrete (HT

N) approximations to the

equation of radiative transfer, Ph.D. thesis, Iowa State University, 2019. hi (7) Trefethen, L. (2000). Spectral Methods in MATLAB, Society for Industrial and Applied Mathematics, 2000.