CUDA/C++ Delta Eddington Solver

[mattldawson]

• Create stub class for Delta Eddington solver in C++ #65
• Set up testing for C++ Delta Eddington solver #66
• Create grid class in C++ #67
• Create profile class in C++ #68
• Variable Initialization Function for Delta Eddington Solver #102
• Function for Assembling the Tridiagonal Linear System for the Delta Eddington Solver  #103
• Function for Computing the Radiation Field  #104
• Translate current solver algorithm into C++ solver
• Replace LAPACK functions with MAGMA versions
• Performance analysis of C++/CUDA Delta Eddington Solver

Solving Radiative Transfer Equation using Delta Eddington Approximation

1. The methodology is described in this paper.
2. Original Delta-Eddington Paper

INPUT and INITIALIZATION

Input Variables	Description
$N$	Number of Layers
Wavelength Grid	Contains wavelengths $\lambda$
Radiator state	Contains optical depths $\tau$, surface albedos $\omega$, assymetry parameters $g$
Spherical geometry	Contains the azimuth ($\phi$) and zenith ($\theta$) angles

[STEP 1 Compute the Delta Eddington coeffcients ($\gamma$ values from [1])] #102

Given:

1. Assymetry parameters ${g_0 \cdots g_n}$
2. Single scattering albedo $\omega_0$
3. Zenith cosine $\mu = cos(\theta)$.

For each layer $n$, do

$$ \begin{align} \gamma_{1n} &= \frac{7 - \omega_0(4+3g_n)}{4} \\ \gamma_{2n} &= -\frac{1 - \omega_0(4-3g_n)}{4} \\ \gamma_{3n} &= \frac{2-3g_n\mu_0}{4} \\ \gamma_{4n} &= 1 - \gamma_{3n} \\ \mu_n &= 0.5 \\ \lambda_n &= \sqrt{\gamma_1^2 - \gamma_2^2} \\ \Gamma_n &= \frac{\gamma_{1n} - \lambda_n}{\gamma_{2n}} \end{align} $$

STEP 2 Define Solar radiation functions : upwards represented as +; downwards -) (#103)

$$ \begin{align} C^+(\tau) &= \frac{\omega_0\pi F_s e^{\frac{-\tau_c+\tau}{\mu_0}}(\frac{\gamma_1 - 1} {\mu_0}\gamma_3+\gamma_4\gamma_2)}{\frac{\lambda^2-1}{\mu_0^2}} \\ C^+(\tau) &= \frac{\omega_0\pi F_s e^{\frac{-\tau_c+\tau}{\mu_0}}(\frac{\gamma_1 + 1} {\mu_0}\gamma_4+\gamma_2\gamma_3)}{\frac{\lambda^2-1}{\mu_0^2}} \end{align} $$

STEP 4 Compute tridiagonal linear system coeffcients. (#103)
For each layer $n$, do

$$ \begin{align} e_{1n} &= 1 + \Gamma_n e^{-\lambda_n \tau_n} \\ e_{2n} &= 1 - \Gamma_n e^{-\lambda_n \tau_n} \\ e_{3n} &= \Gamma_n + e^{-\lambda_n \tau_n} \\ e_{4n} &= \Gamma_n - e^{-\lambda_n \tau_n}. \end{align} $$

[STEP 3 Assemble tridiagonal linear system left hand side.] (#103)

lower diagonal

$$ \begin{equation} A_n = \begin{cases} 0 & |first , element \\ e_{2n}e_{3n}-e_{4n}e_{1n} & |n=odd \\ e_{2n+1}e_{1n}-e_{3n}e_{4n+1} & |n=even \\ e_{1n} - R_{\text{sfc}}e_{3n} & |last , element \end{cases} \end{equation} $$

main diagonal

$$ \begin{align} B_n = \begin{cases} e_{11} & |first , element \\ e_{1n}e_{1n+1}-e_{3n}e_{3n+1} & | n=odd \\ e_{2n}e_{2n+1}-e_{4n}e_{4n+1} & | n=even \\ e_{2n} - R_{\text{sfc}}e_{4n} & |last , element \end{cases} \end{align} $$

upper diagonal

$$ \begin{equation} D_n = \begin{cases} -e_{21} & |first , element \\ e_{3n}e_{4n+1}-e_{1n}e_{2n+1} & |n=odd \\ e_{1n+1}e_{4n+1}-e_{1n}e_{2n+1} & |n=even \\ 0 & |last , element \end{cases} \end{equation} $$

NOTE: $e_{2n+1}$ = e[2][n+1].

STEP 4 Assemble tridiagonal system right hand side #103

$$ \begin{equation} E_n = \begin{cases} F_0^-(0) - C_1^-(0) & |first , element \\ e_{3n}[C_{n+1}^+(0) - C_n^+(\tau_n^+)] + e_{1n}[C_n^-(\tau_n) - C_n^-(0)] & |n=odd \\ e_{2n+1}[C_{n+1}^+(0) - C_n^+(\tau_n^+)] + e_{4n+1}[C_{n+1}^-(0) - C_n^-(\tau_n)] & |n=even \\ S_{\text{sfc}} - C_n^+(\tau_n) + R_{\text{sfc}}C_n^-(\tau_n) & |last , element \end{cases} \end{equation} $$

STEP 5 Solve tridiagonal linear system

$$ \begin{equation} \begin{pmatrix} B_1 & D_1 & 0 & 0 & \dots & 0 \\ A_2 & B_2 & D_2 & 0 & \dots & 0 \\ 0 & A_3 & B_3 & D_3 & \dots & 0 \\ \vdots & \vdots & \vdots & \ddots & \ddots & \vdots \\ 0 & 0 & 0 & \dots & A_n & B_n \end{pmatrix} \begin{pmatrix} Y_1 \\ Y_2 \\ \vdots \\ Y_{2N-1} \\ Y_{2N} \end{pmatrix} = \begin{pmatrix} E_1 \\ E_2 \\ \vdots \\ E_{2N-1} \\ E_2N \end{pmatrix} \end{equation} $$

Solve for $[Y_{1n}, Y_{2n} \cdots Y_{2N-1}, Y_{2N}]^t$

STEP 6 Plug $Y$ back in general solution to compute Flux and Intensity #104

$$ \begin{align} F_n^+ &= Y_{1n} e_{1n} + Y_{2n} e_{2n} + C_n^+( \tau_n ) \\ F_n^- &= Y_{1n} e_{3n} + Y_{2n} e_{4n} + C_n^+( \tau_n ) \\ J_n^+ &= \frac{F^+(\tau)}{ \mu_n} \\ J_n^- &= \frac{F^-(\tau)}{ \mu_n} \end{align} $$