← Octo-TIGER / All test problemsNonequilibrium radiation diffusion

Marshak wave

Boundary illumination heats initially cold stationary material. Radiation and material energy profiles are compared with the semi-analytic Marshak diffusion solution.

The experiment

Physical setup

The normalized cube spans −1 ≤ x,y,z ≤ 1. The material density and momentum are fixed; radiation exchanges energy with the material through a linear equilibrium law Eeq = e. The propagation speed is c = 1, absorption and extinction are both χ = 100, and the default final time is t = 1. Equal material and radiation energy capacities correspond to ε = 1. Gravity and hydrodynamic transport are disabled.

This is a heating test in normalized units, rather than an ideal-gas temperature benchmark. Plots show E/4 and e/4. The saved metadata and source record the actual boundary treatment and reference quadrature.

This case uses normalized units with c = 1, opacity 100, and equal material/radiation energy capacities. Hydrodynamics and gravity are disabled. The default cube spans −1 to +1, and the default final time is 1. Adaptive runs retain their actual leaf spacing in exported data.

The comparison

Comparison and interpretation

The Su–Olson reference is the solution of a linear, gray, nonequilibrium diffusion problem in a semi-infinite stationary medium:

∂tE − (1/3χ)∂x²E = χ(e − E)
∂te = χ(E − e)

With ξ = √3χ(x + 1), τ = χt, u = E/4, and v = e/4, the reference has u − (2/√3)∂ξu = 1 at the illuminated boundary and zero initial u and v. The existing integral quadrature supplies radiation energy, material energy, and diffusion flux. Reference values are evaluated at cell centers; the exported comparison does not apply an additional cell-volume average.

What this test measures

Inspect the penetration depth and delay between radiation and material heating. The full-volume L1, L2, and maximum norms include material energy as well as all four radiation moments. Nonlinear M1 has finite propagation speed and a closure different from diffusion near the illuminated boundary and early front. Consequently these norms include model and boundary differences; they are not an exact-M1 order test. A zero transverse-flux error has no logarithmic convergence order. The numerical boundary combines a prescribed incoming half-isotropic bath with outgoing moments from the interior M1 intensity. Its finite-grid ghost-cell/HLL realization and early causal transient can differ from the diffusion boundary solution.

Measured results

All resolutions

LevelGridBatch status
232³Ready

Convergence requires at least two completed, comparable resolutions.

Full-volume norms and orders

GridFieldNormErrorOrderError units
32³ cellserL11.516612e-01—code units
32³ cellserL22.362060e-01—code units
32³ cellserLinf5.011435e-01—code units
32³ cellsfxL11.678540e-03—code units
32³ cellsfxL23.526129e-03—code units
32³ cellsfxLinf1.695846e-02—code units
32³ cellsfyL10.000000e+00—code units
32³ cellsfyL20.000000e+00—code units
32³ cellsfyLinf0.000000e+00—code units
32³ cellsfzL10.000000e+00—code units
32³ cellsfzL20.000000e+00—code units
32³ cellsfzLinf0.000000e+00—code units
32³ cellsmaterial_egasL11.516600e-01—code units
32³ cellsmaterial_egasL22.362131e-01—code units
32³ cellsmaterial_egasLinf5.011783e-01—code units

32³ cells / level 2

Recorded setup and run provenance
OriginOcto-TIGER application
DomainCube of side 2.000000e+00 (code length); centered at the origin
Resolution32³ cells; finest Δx = 6.250000e-02 (code length)
Final comparison time6.000000e+01 (code time)
Radiation modelM1
Radiation integrationM1 RK1 (first order in time)
Flux semanticsEvolved M1 flux moment
Radiation backendcpu
Physical light speed1.000000e+00 (code speed)
Reduced-speed factor α1.000000e+00
Radiation transport speed1.000000e+00 (code speed)
Periodic boundariesoff
Hydrodynamicsoff
Gravityoff
Implicit radiation sourceon
Uniform gridon
CFL0.4
Timestep cap0.050000000000000003 (code time)
Output interval0.5 (code time)
{
  "background": 0.0,
  "c": 1.0,
  "case": "marshak",
  "cells": 32,
  "command": [
    "/ddnB/work/dmarce1/octotiger-m1-fld-20261014/campaign/control/launch-solver",
    "/ddnB/work/dmarce1/octotiger-m1-fld-20261014/build/octotiger",
    "--config_file=/ddnB/work/dmarce1/octotiger-m1-fld-20261014/campaign/runs/l2-marshak-m1/marshak/l2/run.ini",
    "--hpx:threads=64"
  ],
  "comparison_signature": {
    "atomic_mass": "1",
    "atomic_number": "1",
    "cfl": "0.4",
    "code_to_cm": "1",
    "code_to_g": "1",
    "code_to_s": "1",
    "correct_am_hydro": "off",
    "disable_analytic": "off",
    "disable_diagnostics": "on",
    "eos": "IDEAL",
    "gravity": "off",
    "hard_dt": "0.050000000000000003",
    "hydro": "off",
    "hydro_device_kernel_type": "OFF",
    "hydro_host_kernel_type": "LEGACY",
    "monopole_device_kernel_type": "OFF",
    "multipole_device_kernel_type": "OFF",
    "n_species": "1",
    "omega": "0",
    "periodic": "off",
    "problem": "MARSHAK",
    "rad_implicit": "on",
    "radiation": "on",
    "radiation_algorithm_version": 1,
    "radiation_integrator": "euler",
    "radiation_model": "m1",
    "radiation_reduced_c_factor": "1",
    "stop_step": "1000000",
    "stop_time": "60",
    "unigrid": "on",
    "xscale": "1"
  },
  "config": {
    "atomic_mass": "1",
    "atomic_number": "1",
    "cfl": "0.4",
    "code_to_cm": "1",
    "code_to_g": "1",
    "code_to_s": "1",
    "correct_am_hydro": "off",
    "datadir": "/ddnB/work/dmarce1/octotiger-m1-fld-20261014/campaign/runs/l2-marshak-m1/marshak/l2/",
    "disable_analytic": "off",
    "disable_diagnostics": "on",
    "disable_output": "off",
    "eos": "IDEAL",
    "gravity": "off",
    "hard_dt": "0.050000000000000003",
    "hydro": "off",
    "hydro_device_kernel_type": "OFF",
    "hydro_host_kernel_type": "LEGACY",
    "max_level": "2",
    "min_level": "2",
    "monopole_device_kernel_type": "OFF",
    "multipole_device_kernel_type": "OFF",
    "n_species": "1",
    "odt": "0.5",
    "omega": "0",
    "periodic": "off",
    "problem": "MARSHAK",
    "rad_implicit": "on",
    "radiation": "on",
    "radiation_backend": "cpu",
    "radiation_integrator": "euler",
    "radiation_model": "m1",
    "radiation_reduced_c_factor": "1",
    "stop_step": "1000000",
    "stop_time": "60",
    "unigrid": "on",
    "xscale": "1"
  },
  "conservation": [
    {
      "boundary": -13.624558525091505,
      "boundary_fraction": -0.42576745390910953,
      "field": "er",
      "final": 6.812521669367501,
      "final_fraction": 0.21289130216773441,
      "history_end_time": 60.0,
      "history_scope": "full_run",
      "history_start_time": 0.0,
      "initial": 0.0,
      "initial_fraction": 0.0,
      "integral_units": "code units",
      "max_abs_residual": 5.3290705182007514e-15,
      "max_normalized_error": 1.6653345369377348e-16,
      "normalization_basis": "bath_energy * volume",
      "normalization_scale": 32.0,
      "normalization_time": 0.0,
      "normalized_error": 0.0,
      "normalized_residual": 0.0,
      "raw_change": 6.812521669367501,
      "residual": 0.0,
      "source": -6.812036855724004,
      "source_fraction": -0.2128761517413751,
      "zero_initial_flux_relative_tolerance": 1.4210854715202004e-14
    },
    {
      "boundary": -311.43366853255327,
      "boundary_fraction": -9.73230214164229,
      "field": "fx",
      "final": 0.05268788866313581,
      "final_fraction": 0.0016464965207229941,
      "history_end_time": 60.0,
      "history_scope": "full_run",
      "history_start_time": 0.0,
      "initial": 0.0,
      "initial_fraction": 0.0,
      "integral_units": "code units",
      "max_abs_residual": 5.684341886080802e-14,
      "max_normalized_error": 1.7763568394002505e-15,
      "normalization_basis": "c * bath_energy * volume",
      "normalization_scale": 32.0,
      "normalization_time": 0.0,
      "normalized_error": 0.0,
      "normalized_residual": 0.0,
      "raw_change": 0.05268788866313581,
      "residual": 0.0,
      "source": -311.38098064389015,
      "source_fraction": -9.730655645121567,
      "zero_initial_flux_relative_tolerance": 1.4210854715202004e-14
    },
    {
      "boundary": -5.694003734421379e-21,
      "boundary_fraction": -1.779376167006681e-22,
      "field": "fy",
      "final": 0.0,
      "final_fraction": 0.0,
      "history_end_time": 60.0,
      "history_scope": "full_run",
      "history_start_time": 0.0,
      "initial": 0.0,
      "initial_fraction": 0.0,
      "integral_units": "code units",
      "max_abs_residual": 7.506356456525063e-21,
      "max_normalized_error": 2.345736392664082e-22,
      "normalization_basis": "c * bath_energy * volume",
      "normalization_scale": 32.0,
      "normalization_time": 0.0,
      "normalized_error": 1.779376167006681e-22,
      "normalized_residual": -1.779376167006681e-22,
      "raw_change": 0.0,
      "residual": -5.694003734421379e-21,
      "source": 0.0,
      "source_fraction": 0.0,
      "zero_initial_flux_relative_tolerance": 1.4210854715202004e-14
    },
    {
      "boundary": 6.310868495560712e-21,
      "boundary_fraction": 1.9721464048627224e-22,
      "field": "fz",
      "final": 0.0,
      "final_fraction": 0.0,
      "history_end_time": 60.0,
      "history_scope": "full_run",
      "history_start_time": 0.0,
      "initial": 0.0,
      "initial_fraction": 0.0,
      "integral_units": "code units",
      "max_abs_residual": 1.4002821011691316e-20,
      "max_normalized_error": 4.375881566153536e-22,
      "normalization_basis": "c * bath_energy * volume",
      "normalization_scale": 32.0,
      "normalization_time": 0.0,
      "normalized_error": 1.9721464048627224e-22,
      "normalized_residual": 1.9721464048627224e-22,
      "raw_change": 0.0,
      "residual": 6.310868495560712e-21,
      "source": 0.0,
      "source_fraction": 0.0,
      "zero_initial_flux_relative_tolerance": 1.4210854715202004e-14
    }
  ],
  "dx": 0.0625,
  "executable_sha256": "ca38f12e338eed0b7074dbb15c70c089a8d31d1c12aa87f5b044c42824721f0b",
  "inx": 8,
  "length": 2.0,
  "level": 2,
  "marshak_capacity_ratio": 1.0,
  "marshak_opacity": 100.0,
  "movie_capture": {
    "hard_dt": 0.05,
    "odt": 0.5,
    "requested_snapshots": 121,
    "steps_per_output_check": 1
  },
  "norms": {
    "er": {
      "L1": 0.1516612,
      "L2": 0.236206,
      "Linf": 0.5011435
    },
    "fx": {
      "L1": 0.00167854,
      "L2": 0.003526129,
      "Linf": 0.01695846
    },
    "fy": {
      "L1": 0.0,
      "L2": 0.0,
      "Linf": 0.0
    },
    "fz": {
      "L1": 0.0,
      "L2": 0.0,
      "Linf": 0.0
    },
    "material_egas": {
      "L1": 0.15166,
      "L2": 0.2362131,
      "Linf": 0.5011783
    }
  },
  "origin": "Octo-TIGER application",
  "physical_c": 1.0,
  "radiation_algorithm_version": 1,
  "radiation_backend": "cpu",
  "radiation_equations": "physical-c",
  "radiation_flux_semantics": "evolved_m1_flux",
  "radiation_integrator": "euler",
  "radiation_model": "m1",
  "radiation_reduced_c_factor": 1.0,
  "reference_kind": "Su-Olson nonequilibrium diffusion reference",
  "status": "complete",
  "time": 60.0,
  "transport_c": 1.0,
  "units": {
    "c": 1.0,
    "energy_density": "normalized energy density",
    "flux": "normalized flux",
    "length": "normalized length",
    "system": "dimensionless",
    "time": "normalized time"
  },
  "validation_scope": "Normalized stationary material, equal radiation/material energy capacities. Analytic diffusion comparison is not an exact solution of nonlinear M1 transport."
}

Final fields and profiles

Central cell layer at z=dx/2; ray at y=z=dx/2.

Numerical, reference, and signed-error energy slices. Central cell layer at z=dx/2; ray at y=z=dx/2.
Numerical, reference, and signed-error energy slices. Central cell layer at z=dx/2; ray at y=z=dx/2. · PDF
Numerical/reference line profiles and signed errors for all four radiation fields
Numerical/reference line profiles and signed errors for all four radiation fields · PDF
Normalized radiation and material energies with their semi-analytic diffusion references
Normalized radiation and material energies with their semi-analytic diffusion references · PDF
All three flux maps
fx: numerical, reference, and signed-error slice
fx: numerical, reference, and signed-error slice · PDF
fy: numerical, reference, and signed-error slice
fy: numerical, reference, and signed-error slice · PDF
fz: numerical, reference, and signed-error slice
fz: numerical, reference, and signed-error slice · PDF

Saved-state movie

Marshak · L2 · M1 RK1

121 stored states; 121 shown. Physical time: 0–60 normalized time. Playback: 60 s.

Zero-order hold; no solution interpolation. Playback holds each saved numerical state over its time interval and includes endpoint holds.

All 121 snapshot times

Snapshot times: 0, 0.505208333, 1, 1.50520833, 2.00520833, 2.5, 3, 3.5, 4, 4.5, 5.00520833, 5.50520833, 6.00520833, 6.50520833, 7.00520833, 7.50520833, 8.00520833, 8.50520833, 9.00520833, 9.50520833, 10, 10.5, 11, 11.5, 12, 12.5, 13, 13.5, 14, 14.5, 15, 15.5, 16, 16.5, 17, 17.5, 18, 18.5, 19, 19.5052083, 20.0052083, 20.5052083, 21.0052083, 21.5052083, 22.0052083, 22.5052083, 23.0052083, 23.5052083, 24.0052083, 24.5052083, 25.0052083, 25.5052083, 26.0052083, 26.5052083, 27.0052083, 27.5052083, 28.0052083, 28.5052083, 29.0052083, 29.5052083, 30.0052083, 30.5052083, 31.0052083, 31.5052083, 32.0052083, 32.5052083, 33.0052083, 33.5052083, 34.0052083, 34.5052083, 35.0052083, 35.5052083, 36.0052083, 36.5052083, 37.0052083, 37.5052083, 38.0052083, 38.5, 39, 39.5, 40, 40.5, 41, 41.5, 42, 42.5, 43, 43.5, 44, 44.5, 45, 45.5, 46, 46.5, 47, 47.5, 48, 48.5, 49, 49.5, 50, 50.5, 51, 51.5, 52, 52.5, 53, 53.5, 54, 54.5, 55, 55.5, 56, 56.5, 57, 57.5, 58, 58.5, 59, 59.5, 60 normalized time.

Arrows show the independently evolved M1 radiation flux. Arrows have equal length and show direction only.

Movies show saved numerical states; model-specific analytic comparisons are in the final-time profiles.

Central z=0+ radiation energy slice with faint cell mesh and equal-length flux direction arrows. Arrows show the independently evolved M1 radiation flux. Fixed color limits; actual stored Silo times held without solution interpolation. Movies show saved numerical states; model-specific analytic comparisons are in the final-time profiles. The maps show the full simulation domain without cropping.

Download MP4 · Rendering and timing metadata · Full run report

All campaign movies

Conservation

All values below are dimensionless fractions of the fixed normalization scale at t_ref = 0.

FieldScaleQ(t_ref)/scaleQ(t)/scaleBoundary B/scaleSource S/scaleResidual R/scaleMaximum |R|/scale
erEbath V0.000000e+002.128913e-01-4.257675e-01-2.128762e-010.000000e+001.665335e-16
fxc Ebath V0.000000e+001.646497e-03-9.732302e+00-9.730656e+000.000000e+001.776357e-15
fyc Ebath V0.000000e+000.000000e+00-1.779376e-220.000000e+00-1.779376e-222.345736e-22
fzc Ebath V0.000000e+000.000000e+001.972146e-220.000000e+001.972146e-224.375882e-22
Integrated radiation budgets and their transport/source-corrected residuals
Integrated radiation budgets and their transport/source-corrected residuals · PDF

Data and reproducibility

Momentum plots are omitted for this fixed-material problem. The prescribed medium acts as an external momentum reservoir, so radiation momentum exchange does not imply that a compensating gas momentum was evolved.

Reading the results

Errors, budgets, and movies

Every norm uses the full three-dimensional domain, with cell-volume weighting. Here e = numerical − reference and V is the domain volume:

L1 = Σ |e| ΔV / V,   L2 = √(Σ e² ΔV / V),   L∞ = max |e|
p = ln(errorcoarse / errorfine) / ln(Δxcoarse / Δxfine)

Orders require two completed runs with matching physical time, configuration, units, executable, and origin. An exact zero error has no logarithmic order. The central slice and line profiles are visual diagnostics: the report slice is the cell layer at z = Δx/2, and its line profile is at y = z = Δx/2.

Conservation is assessed using the measured domain integrals and the signed cumulative boundary and source budgets:

A labeled checkpoint continuation uses its checkpoint as t_ref, with counters and fixed normalization starting there. Its reported residuals cover that resumed segment only; original histories remain separate.

R(t) = Q(t) − Q(0) + B(t) − S(t)

B is outward transport and S is the source change. Conservation plots and this page's budget tables show fractions of the initial totals: Q(t)/D, [Q(0) − B(t) + S(t)]/D, and the signed residual R(t)/D. The denominator D = |Q(0)| is fixed for the entire history. It never grows with the current total or accumulated transport/source budgets. A negative initial flux therefore starts at −1, and the plotted fraction retains its sign.

A component whose initial net flux is zero or negligible cannot be divided by its own initial total. It instead uses the fixed scale D = c |E₀|, where E₀ is the initial integrated radiation energy. This fallback is labeled on its panels and in the table. “Negligible” means |Q(0)|/(c |E₀|) ≤ 64ε ≈ 1.421085 × 10−14, with ε the double-precision machine epsilon. This prevents roundoff left by cancellation in a symmetric field from becoming the denominator. For the initially cold Marshak case, the prescribed bath sets D = Ebath V for energy and D = c Ebath V for flux, with Ebath = 4/c. These scales remain fixed even though radiation enters the initially empty domain.

Raw histories and dimensional summary columns remain in physical units: energy integrals in erg, physical-flux integrals in erg cm/s for CGS runs. Dividing a flux integral by c² gives radiation momentum. Fraction-history CSV files are dimensionless except for their time column; summary files record each denominator and its basis. Maximum fractional residuals use the same fixed denominator as the final residual.

For runs labeled RSLA, these budgets check the modified radiation equations. They do not establish conservation of physical gas-plus-radiation energy. With reduced-speed factor α, the closed-system coupled invariants are Egas + Er/α and gas momentum + F/(αc²), with physical c in the closure. Compare multiple α values before interpreting an RSLA result as an approximation to physical-c evolution.

Movies show the numerical Silo states with a fixed color range for each movie. Recorded states are held for their actual simulation-time gaps; intermediate solutions are not interpolated. The 3D option shows the domain's exterior surface. Movies and plots retain their recorded physical units.

Sources & provenance

Where this problem comes from

Su & Olson (1996), “Benchmark results for the non-equilibrium Marshak diffusion problem,” JQSRT 56, 337–351, gives the nonequilibrium diffusion benchmark. The local reference is implemented in src/test_problems/marshak/marshak.cpp; the run records source hashes for reproducibility.

Radiation method context

M1 RK1 evolves radiation energy and flux moments with the retained first-order radiation update. FLD uses backward-Euler radiation diffusion/material exchange with the Levermore–Pomraning limiter. Moving cases add first-order conservative gas transport and radiation force/work; their reconstructed lab-frame flux includes (E I+P) v as well as −D∇E. The scalar implicit transport solve uses composite-grid multigrid and an iterative linear solver. Both time integrators are first order.

The suite includes Gaussian scattering, Marshak heating, a Gaussian emitter, explicit thin positive-opacity wave/front comparisons, and moving-gas shock/blast problems on CPU at physical light speed. A configured case does not establish that a solve converged; only validated completed runs appear as results. Flux diagnostics in FLD have no separate M1 momentum conservation test. A difference between paired M1 and FLD solutions is a model comparison, not an analytic error. ADER2 has been removed from the current solver; historical result files retain their original provenance.

Skinner & Ostriker (2013) and Hanawa & Audit (2014) provide context for the M1 equations. The implementation and saved executable/source hashes identify the actual numerical methods used in this batch.

Reproducing this batch

Descriptions document the reviewed source implementation. Older or edited solver versions can use different profile constants or averaging rules; these are not all encoded in run.ini. Use the run metadata, executable hash, and saved batch source hashes to identify the implementation behind a result. A source path or descriptive formula alone is not proof of the executable's contents.