Skip to content

2 · Vector RTE & Discretization

For: readers building the equation-level mental model — atmospheric scientists, theorists, retrieval developers who want to understand what the solver actually solves.

Prev: 1 · The Problem & MOM Thesis · Next: 3 · Layer Optical Properties

This page writes the equation that vSmartMOM solves, then discretizes it. By the end, the per-layer problem reduces to four arrays of shape (NquadN, NquadN, nSpec) — the arrays the rest of the RT basics arc manipulates.

The vector RTE

For a plane-parallel atmosphere illuminated by direct solar flux along , the Stokes vector for diffuse light obeys (Sanghavi 2014, Eq. 2):

Reading the right-hand side term by term:

  1. — extinction along the path.

  2.   — thermal emission (zero in solar bands).

  3. The single-scatter direct-beam source, exponentially attenuated from TOA.

  4. The diffuse multiple-scatter integral over all incoming directions .

Sign convention:   is downward,   is upward.   is the single-scattering albedo. is the 4×4 phase matrix (Sanghavi 2014, Eq. 8). The diffuse solution adds to the direct delta-source    to give the total radiation field; we solve only for the diffuse part to avoid the singularity.

Boundary conditions

At the top of atmosphere (Sanghavi 2014, Eq. 6):

At the surface, an emissivity and a polarized BRDF close the system (Sanghavi 2014, Eqs. 33–37). Concepts/05 covers the surface in detail.

Fourier decomposition in azimuth

Because the problem is azimuthally symmetric about the solar direction, can be expanded as a finite Fourier series in   (Sanghavi 2014, Eq. 8):

is block-diagonal (mixing intensity with linear polarization ); is block-off-diagonal (mixing with the rest). Each Fourier moment solves an independent RTE — the -loop in rt_run.jl is therefore embarrassingly parallel.

Quadrature in μ

After Fourier decomposition, each per-moment integral over    is replaced by a finite sum. vSmartMOM offers two quadrature schemes (src/CoreRT/tools/rt_set_streams.jl):

  • GaussLegQuad — half-space Gauss-Legendre on ; the solar zenith angle and viewing zenith angles are appended as zero-weight nodes. Cheap; interpolation between roots loses accuracy at off-quadrature angles.

  • RadauQuad — block-Radau composite on   with as a true quadrature node carrying non-zero weight (Sanghavi 2014 App. B, Eqs. B.1–B.2). Required when the solar direction does not coincide with a Gauss root, because vSmartMOM excludes solar SFI from (see callout below) and therefore relies on direct ray-tracing through .

Design choice — no solar SFI in Sanghavi 2014 §2.2: vSmartMOM does not include single scattering of the direct solar beam in the source term . Two consequences:

  1. The matrix-operator equations reduce from Sanghavi 2014 Eqs. (23)–(33) to Eqs. (23)–(28) for solar bands (no term). Doubling and adding are correspondingly cheaper.

  2. For thermal-only ( ) , the source is isotropic, so its computation can be reused.

The price is needing Block-Radau quadrature for off-quadrature solar/viewing angles. This is a deliberate departure from VLIDORT, SCIATRAN, and the Plass/Fell/Hollstein matrix-operator formulations. :::

The supermatrix form (per Fourier moment)

For each Fourier moment , the RTE becomes (Sanghavi 2014, Eq. 12):

After quadrature in , the upwelling and downwelling streams are stacked into supervectors of length (Sanghavi 2014, Eqs. 14–15):

  and   are diagonal supermatrices of quadrature angles and weights. The four blocks are the directional pairings of the phase matrix.

The / superscripts on that you see throughout the source code (r⁻⁺, t⁺⁺, R⁻⁺, etc.) are exactly this convention: is downward, is upward.

The four arrays each layer reduces to

For each Fourier moment , every atmospheric layer is captured by four arrays of shape (NquadN, NquadN, nSpec) where   :

SymbolMeaningShape
optical depth (per spectral point)(nSpec,)
single-scattering albedo(nSpec,)
forward (downward → downward) Fourier-moment phase matrix(NquadN, NquadN, nSpec)
backscatter (downward → upward) Fourier-moment phase matrix(NquadN, NquadN, nSpec)

The third dimension is spectral. That is the whole batched-matmul story — see Concepts/07. Once these four arrays exist, the elemental/doubling/adding kernels of Concepts/04 are pure linear algebra over the spectral axis.

Polarization is a type, not a runtime branch

Stokes_I ( , scalar intensity), Stokes_IQ ( ), Stokes_IQU ( ), Stokes_IQUV ( ) each subtype AbstractPolarizationType in src/Scattering/types.jl:92–143. Each carries the polarization symmetry vector (used in Concepts/04) and the unit solar Stokes vector .

The kernels specialize on the type at compile time. There is no if pol_type.n == 4 branch inside the inner loops — Julia's multiple dispatch emits separate code paths for each polarization at JIT time. Adding Stokes_IQ took adding one struct definition; the kernels picked it up automatically.

Code anchors

ConceptSource
Top-level RT loopsrc/CoreRT/rt_run.jl:53–329
Fourier-moment dispatchsrc/CoreRT/rt_run.jl:208 (for m = 0:max_m-1)
Quadrature constructionsrc/CoreRT/tools/rt_set_streams.jl:24–110
QuadPoints structsrc/CoreRT/types.jl::QuadPoints
Polarization typessrc/Scattering/types.jl:92–143
Phase matrix Z from Greeksrc/Scattering/compute_Z_matrices.jl::compute_Z_moments
Supermatrix layer typessrc/CoreRT/types.jl::CompositeLayer/AddedLayer

Hands-on tutorials

Runnable examples with Plotly figures:

References

  • Sanghavi et al. (2014), vSmartMOM, JQSRT 133:412–433, doi:10.1016/j.jqsrt.2013.09.004. Eqs. (2), (4), (6), (8)–(10), (12), (14)–(15); App. B (Block-Radau). Primary reference for everything on this page.

  • Chandrasekhar (1950), Radiative Transfer. (Original vector RTE in this form.)

  • Hovenier (1971), Symmetry relationships for scattering of polarized light…, J. Atmos. Sci. 28:488. (Block diagonal/anti-diagonal structure of .)

  • Full crib sheet: docs/dev_notes/theory_references.md §A, §E.