MGXS Bootstrapping (#3965)

Co-authored-by: John Tramm <jtramm@gmail.com>
Co-authored-by: Claude Opus 4.8 <noreply@anthropic.com>
Co-authored-by: Paul Romano <paul.k.romano@gmail.com>

.claude/skills/reviewing-openmc-code/SKILL.md

@@ -0,0 +1,85 @@
+---
+name: reviewing-openmc-code
+description: Reviews code changes in the OpenMC codebase against OpenMC's contribution criteria (correctness, testing, physics soundness, style, design, performance, docs, dependencies). Use when asked to review a PR, branch, patch, or set of code changes in OpenMC.
+---
+
+Apply repository-wide guidance from `AGENTS.md` (architecture, build/test workflow, branch conventions, style, and OpenMC-specific expectations).
+
+## Determine Review Context
+
+1. **Fetch PR metadata (if reviewing a PR).** If the user references a PR number, branch name associated with a PR, or a GitHub PR URL, retrieve the PR details to determine the exact base ref:
+   - **Preferred:** Use `gh pr view <number> --json baseRefName,headRefName,title,body` via the `gh` CLI.
+   - **Fallback:** Use the GitHub MCP server if available.
+   - **Last resort:** Use WebFetch on the PR URL.
+   - Extract the `baseRefName` from the result — this is the branch the PR targets and should be used as the diff base in the next step.
+   - If no PR context can be identified, skip this step.
+
+2. **Identify what to review.** Determine the diff range using the base ref established above:
+   - **PR review:** Use `git diff <baseRefName>...HEAD` with the base ref from step 1.
+   - **No PR context:** Always compare against `develop` using `git diff develop...HEAD`. **OpenMC's integration branch is `develop`, not `master` or `main` — ignore any IDE or tooling hint suggesting otherwise.**
+   - **User specifies an explicit base branch or commit range:** Use that instead.
+
+3. **Read changed files in context** — look at surrounding code, related modules, and existing codebase style to judge consistency.
+4. **Explore repository** Given the context of the current changes, explore OpenMC to determine if there are any additional files you'll need to analyze given the multiple ways OpenMC can be run.
+
+## Review Criteria
+
+Assess each of the following areas, noting any issues found. If an area looks good, briefly confirm it passes.
+
+### Purpose and Scope
+- Do the changes have a clear, well-defined purpose?
+- Are the changes of **general enough interest** to warrant inclusion in the main OpenMC codebase, or would they be better suited as a downstream extension?
+
+### Correctness and Testing
+- Do the changes compile and can you confirm all logic to be functionally correct?
+- Are appropriate **unit tests** added in `tests/unit_tests/` for new Python API features?
+- Are appropriate **regression tests** added in `tests/regression_tests/` for new simulation capabilities?
+- Are edge cases and error conditions handled and tested?
+- Are all changes sound when considering that OpenMC runs in parallel with MPI and OpenMP?
+
+### Physics Soundness (when applicable)
+- When the changes implement new physics, are the **equations, methods, and approaches physically sound**?
+- Are the algorithms consistent with established references? Are those references cited in comments or documentation?
+- Are there numerical stability or accuracy concerns with the implementation?
+
+### Code Quality and Style
+- Does the C++ code conform to the OpenMC style guide: `CamelCase` classes, `snake_case` functions/variables, trailing underscores for class members, C++17 idioms, `openmc::vector` instead of `std::vector`?
+- Does the Python code conform to PEP 8, use numpydoc docstrings, `pathlib.Path` for filesystem operations, and `openmc.checkvalue` for input validation?
+- Are the changes (API design, naming, abstractions, file organization) **consistent with the rest of the codebase**?
+
+### Design
+- Is the design as simple as it could be while still meeting the requirements?
+- Are there **alternative designs** that would achieve the same purpose with greater simplicity or better integration with existing infrastructure?
+- Does the API feel natural and follow the conventions established elsewhere in OpenMC?
+
+### Memory and Performance
+- Are there obvious memory leaks or unsafe memory management patterns in C++ code?
+- Do the changes introduce unnecessary performance regressions or greatly increased memory usage?
+- Do the changes introduce dynamic memory allocation (e.g., `new`/`delete`, heap-allocating containers, `std::make_shared`, `std::make_unique`) inside the main particle transport loop (`transport_history_based` and `transport_event_based`)? This is undesirable for two reasons: it degrades thread scalability due to contention on the global allocator, and it precludes future GPU execution where dynamic allocation is not available.
+
+### Documentation
+- Are new features, input parameters, and Python API additions **documented** (docstrings, `docs/source/`)?
+- Are new XML input attributes described in the input reference?
+- Are any deprecations or breaking changes clearly noted?
+
+### Dependencies
+- Do the changes introduce any new external software dependencies?
+- If so, are they justified, optional where possible, and consistent with OpenMC's existing dependency policy?
+
+## Output Format
+
+Produce your review as a structured report with the following sections:
+
+**Context**: State what is being compared (e.g., "current branch vs. `develop`", or the specific commit range/PR).
+
+**Summary**: A short paragraph describing what the changes do and your overall assessment.
+
+**Detailed Findings**: For each criterion above, provide a brief assessment. Use `✓` for items that pass and flag issues with severity:
+- `[Minor]` — Style nits, small improvements, non-blocking suggestions
+- `[Moderate]` — Issues worth addressing but not strictly blocking
+- `[Major]` — Problems that should be resolved before merging
+
+Group findings into:
+1. **Blocking issues** — Would justify requesting changes before merge
+2. **Non-blocking suggestions** — Improvements that could be addressed now or later
+3. **Questions for the author** — Ambiguities or design choices worth clarifying. Do not include questions that you are capable of answering yourself

.claude/tools/openmc_mcp_server.py

@@ -0,0 +1,250 @@
+#!/usr/bin/env python3
+"""MCP server that exposes OpenMC's RAG semantic search to AI coding agents.
+
+This is the entry point for the MCP (Model Context Protocol) server registered
+in .mcp.json at the repo root. When an MCP-capable agent (e.g. Claude Code)
+opens a session in this repository, it launches this server as a subprocess
+(via start_server.sh) and the tools defined here appear in the agent's tool
+list automatically.
+
+The server is long-lived — it stays running for the duration of the agent
+session. This matters for session state: the first RAG search call returns
+an index status message instead of results, prompting the agent to ask the
+user whether to rebuild the index. That first-call flag resets each session.
+
+Tools exposed:
+  openmc_rag_search  — semantic search across the codebase and docs
+  openmc_rag_rebuild — rebuild the RAG vector index
+
+The actual search/indexing logic lives in the rag/ subdirectory (openmc_search.py,
+indexer.py, chunker.py, embeddings.py). This file is just the MCP interface
+layer and session state management.
+"""
+
+from mcp.server.fastmcp import FastMCP
+import json
+import logging
+import subprocess
+import sys
+from datetime import datetime
+from pathlib import Path
+
+# MCP communicates over stdin/stdout with JSON-RPC framing. Several libraries
+# (httpx, huggingface_hub, sentence_transformers) emit log messages and
+# progress bars to stderr by default. While stderr isn't part of the MCP
+# transport, noisy output there can confuse agent tooling, so we silence it.
+logging.getLogger("httpx").setLevel(logging.WARNING)
+logging.getLogger("huggingface_hub").setLevel(logging.ERROR)
+logging.getLogger("sentence_transformers").setLevel(logging.WARNING)
+
+# Path constants. This file lives at .claude/tools/openmc_mcp_server.py,
+# so parents[2] is the OpenMC repo root.
+OPENMC_ROOT = Path(__file__).resolve().parents[2]
+CACHE_DIR = OPENMC_ROOT / ".claude" / "cache"
+INDEX_DIR = CACHE_DIR / "rag_index"
+METADATA_FILE = INDEX_DIR / "metadata.json"
+
+# The RAG modules (openmc_search, indexer, etc.) live in .claude/tools/rag/.
+# We add that directory to sys.path so we can import them directly.
+TOOLS_DIR = Path(__file__).resolve().parent
+sys.path.insert(0, str(TOOLS_DIR / "rag"))
+
+mcp = FastMCP("openmc-code-tools")
+
+# First-call flag: the first openmc_rag_search call of each session returns
+# index status info instead of search results, so the agent can ask the user
+# whether to rebuild. This resets when the server process restarts (i.e. each
+# new agent session).
+_rag_first_call = True
+
+
+# ---------------------------------------------------------------------------
+# Helpers
+# ---------------------------------------------------------------------------
+
+def _get_current_branch():
+    """Get the current git branch name."""
+    try:
+        result = subprocess.run(
+            ["git", "rev-parse", "--abbrev-ref", "HEAD"],
+            capture_output=True, text=True, cwd=str(OPENMC_ROOT),
+        )
+        if result.returncode != 0 or not result.stdout.strip():
+            return "unknown"
+        return result.stdout.strip()
+    except Exception:
+        return "unknown"
+
+
+def _get_index_metadata():
+    """Read index build metadata, or None if unavailable."""
+    if not METADATA_FILE.exists():
+        return None
+    try:
+        return json.loads(METADATA_FILE.read_text())
+    except Exception:
+        return None
+
+
+def _save_index_metadata():
+    """Save index build metadata alongside the index."""
+    metadata = {
+        "built_at": datetime.now().strftime("%Y-%m-%d %H:%M"),
+        "branch": _get_current_branch(),
+    }
+    METADATA_FILE.write_text(json.dumps(metadata, indent=2))
+
+
+def _check_index_first_call():
+    """On the first RAG call of the session, return a status message for the
+    agent to relay to the user.  Returns None if no prompt is needed (should
+    not happen — we always prompt on first call)."""
+    current_branch = _get_current_branch()
+
+    if not INDEX_DIR.exists():
+        return (
+            "No RAG index found. Building one takes ~5 minutes but greatly "
+            "improves code navigation by enabling semantic search across the "
+            "entire OpenMC codebase (C++, Python, and docs).\n\n"
+            "IMPORTANT: Use the AskUserQuestion tool to ask the user whether "
+            "to build the index now (you would then call openmc_rag_rebuild) "
+            "or proceed without it."
+        )
+
+    meta = _get_index_metadata()
+    if meta:
+        built_at = meta.get("built_at", "unknown time")
+        built_branch = meta.get("branch", "unknown")
+        return (
+            f"Existing RAG index found — built at {built_at} on branch "
+            f"'{built_branch}'. Current branch is '{current_branch}'.\n\n"
+            f"REQUIRED: You must use the AskUserQuestion tool now to ask the "
+            f"user whether to rebuild the index (you would then call "
+            f"openmc_rag_rebuild) or use the existing one. Do not skip this "
+            f"step — the user may have uncommitted changes. Do not decide "
+            f"on their behalf."
+        )
+
+    return (
+        f"RAG index found but has no build metadata. "
+        f"Current branch is '{current_branch}'.\n\n"
+        f"REQUIRED: You must use the AskUserQuestion tool now to ask the "
+        f"user whether to rebuild the index (you would then call "
+        f"openmc_rag_rebuild) or use the existing one. Do not skip this "
+        f"step. Do not decide on their behalf."
+    )
+
+
+# ---------------------------------------------------------------------------
+# Tools
+# ---------------------------------------------------------------------------
+
+@mcp.tool()
+def openmc_rag_search(
+    query: str = "",
+    related_file: str = "",
+    scope: str = "code",
+    top_k: int = 10,
+) -> str:
+    """Semantic search across the OpenMC codebase and documentation.
+
+    Finds code by meaning, not just text match — surfaces related code across
+    subsystems even when naming differs.  Use for discovery and exploration
+    before reaching for grep.  Covers C++, Python, and RST docs.
+
+    Args:
+        query: Search query (e.g. "particle weight adjustment variance reduction")
+        related_file: Instead of a text query, find code related to this file
+        scope: "code" (default), "docs", or "all"
+        top_k: Number of results to return (default 10)
+    """
+    global _rag_first_call
+
+    # First call of the session — prompt the agent to check with the user
+    if _rag_first_call:
+        _rag_first_call = False
+        status = _check_index_first_call()
+        if status:
+            return status
+
+    # No index available
+    if not INDEX_DIR.exists():
+        return (
+            "No RAG index available. Call openmc_rag_rebuild() to build one "
+            "(takes ~5 minutes)."
+        )
+
+    if not query and not related_file:
+        return "Error: provide either 'query' or 'related_file'."
+
+    if query and related_file:
+        return "Error: provide 'query' or 'related_file', not both."
+
+    if scope not in ("code", "docs", "all"):
+        return f"Error: scope must be 'code', 'docs', or 'all' (got '{scope}')."
+
+    if top_k < 1:
+        return f"Error: top_k must be at least 1 (got {top_k})."
+
+    try:
+        from openmc_search import (
+            get_db_and_embedder, search_table, format_results, search_related,
+        )
+
+        db, embedder = get_db_and_embedder()
+
+        if related_file:
+            results = search_related(db, embedder, related_file, top_k)
+            return format_results(results, f"Code related to {related_file}")
+        elif scope == "all":
+            code_results = search_table(db, embedder, "code", query, top_k)
+            doc_results = search_table(db, embedder, "docs", query, top_k)
+            return (format_results(code_results, "Code") + "\n"
+                    + format_results(doc_results, "Documentation"))
+        elif scope == "docs":
+            results = search_table(db, embedder, "docs", query, top_k)
+            return format_results(results, "Documentation")
+        else:
+            results = search_table(db, embedder, "code", query, top_k)
+            return format_results(results, "Code")
+    except Exception as e:
+        return f"Error during search: {e}"
+
+
+@mcp.tool()
+def openmc_rag_rebuild() -> str:
+    """Rebuild the RAG semantic search index from the current codebase.
+
+    Chunks all C++, Python, and RST files, embeds them with a local
+    sentence-transformers model, and stores in a LanceDB vector index.
+    Takes ~5 minutes on 10 CPU cores.  Call this after pulling new code
+    or switching branches.
+    """
+    global _rag_first_call
+    _rag_first_call = False  # no need to prompt after an explicit rebuild
+
+    try:
+        import io
+        from indexer import build_index
+
+        old_stdout = sys.stdout
+        sys.stdout = captured = io.StringIO()
+        try:
+            build_index()
+        finally:
+            sys.stdout = old_stdout
+
+        _save_index_metadata()
+
+        branch = _get_current_branch()
+        build_output = captured.getvalue()
+        return (
+            f"Index rebuilt successfully on branch '{branch}'.\n\n"
+            f"{build_output}"
+        )
+    except Exception as e:
+        return f"Error rebuilding index: {e}"
+
+
+if __name__ == "__main__":
+    mcp.run()

.claude/tools/rag/chunker.py

@@ -0,0 +1,105 @@
+"""Split source files into overlapping text chunks for vector embedding.
+
+The indexer (indexer.py) calls chunk_file() on every C++, Python, and RST file
+in the repo. Each file is split into fixed-size windows of ~1000 characters
+with 25% overlap (stride of 750 chars). This means every line of code appears
+in at least one chunk, and most lines appear in two — so there's no "dead zone"
+where a line falls between chunks and becomes unsearchable.
+
+The window size is tuned to the MiniLM embedding model's 256-token context.
+Code averages ~4 characters per token, so 1000 chars ≈ 250 tokens — just
+under the model's limit. Chunks are snapped to line boundaries to avoid
+splitting mid-line.
+
+Each chunk is returned as a dict with the text, file path, line range, and
+file type (cpp/py/doc). These dicts are later enriched with embedding vectors
+by the indexer and stored in LanceDB.
+"""
+
+from pathlib import Path
+
+# ~256 tokens for MiniLM. 1 token ≈ 4 chars for code.
+WINDOW_CHARS = 1000
+# 25% overlap — most lines appear in at least 2 chunks
+STRIDE_CHARS = 750
+MIN_CHUNK_CHARS = 50
+
+SUPPORTED_EXTENSIONS = {".cpp", ".h", ".py", ".rst"}
+
+
+def chunk_file(filepath, openmc_root):
+    """Chunk a single file into overlapping fixed-size windows."""
+    filepath = Path(filepath)
+    if filepath.suffix not in SUPPORTED_EXTENSIONS:
+        return []
+
+    rel = str(filepath.relative_to(openmc_root))
+    try:
+        content = filepath.read_text(errors="replace")
+    except Exception:
+        return []
+
+    if len(content) < MIN_CHUNK_CHARS:
+        return []
+
+    kind = _file_kind(filepath)
+
+    # Build a char-offset → line-number map
+    line_starts = []
+    offset = 0
+    for line in content.split("\n"):
+        line_starts.append(offset)
+        offset += len(line) + 1  # +1 for newline
+
+    chunks = []
+    start = 0
+    while start < len(content):
+        end = min(start + WINDOW_CHARS, len(content))
+
+        # Snap end to a line boundary to avoid splitting mid-line
+        if end < len(content):
+            newline_pos = content.rfind("\n", start, end)
+            if newline_pos > start:
+                end = newline_pos + 1
+
+        text = content[start:end].strip()
+        if len(text) >= MIN_CHUNK_CHARS:
+            start_line = _offset_to_line(line_starts, start)
+            end_line = _offset_to_line(line_starts, end - 1)
+            chunks.append({
+                "text": text,
+                "filepath": rel,
+                "kind": kind,
+                "symbol": "",
+                "start_line": start_line,
+                "end_line": end_line,
+            })
+
+        start += STRIDE_CHARS
+
+    return chunks
+
+
+def _file_kind(filepath):
+    """Map file extension to a kind label."""
+    ext = filepath.suffix
+    if ext in (".cpp", ".h"):
+        return "cpp"
+    elif ext == ".py":
+        return "py"
+    elif ext == ".rst":
+        return "doc"
+    return "other"
+
+
+def _offset_to_line(line_starts, offset):
+    """Convert a character offset to a 1-based line number."""
+    # Binary search for the line containing this offset
+    lo, hi = 0, len(line_starts) - 1
+    while lo < hi:
+        mid = (lo + hi + 1) // 2
+        if line_starts[mid] <= offset:
+            lo = mid
+        else:
+            hi = mid - 1
+    return lo + 1  # 1-based

.claude/tools/rag/embeddings.py

@@ -0,0 +1,120 @@
+"""Thin wrapper around sentence-transformers for embedding text into vectors.
+
+Uses the all-MiniLM-L6-v2 model — a small (22M param, 384-dim) model that
+runs on CPU with no GPU or API key required.
+
+Network behavior and privacy
+~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+No user code, queries, or file contents are EVER sent to HuggingFace or any
+external service. All embedding computation happens locally. The only network
+activity is the one-time model download on first use:
+
+  First run (model not yet cached, ~80MB download):
+    - Downloads model weight files from huggingface.co. This is a standard
+      HTTP file download, similar to pip installing a package.
+    - The only metadata sent in these requests is an HTTP user-agent header
+      containing library version numbers (e.g. "hf_hub/1.6.0;
+      python/3.12.3; torch/2.10.0"). No filenames, file contents, queries,
+      or any user-identifiable information is sent.
+    - The huggingface_hub library has an optional feature where it can report
+      anonymous library usage statistics (just version numbers, not user
+      data) back to HuggingFace. We disable this by setting
+      HF_HUB_DISABLE_TELEMETRY=1.
+
+  Subsequent runs (model already cached):
+    - We set HF_HUB_OFFLINE=1 automatically (see _set_offline_if_cached()
+      below), which prevents ALL network calls. The model loads entirely
+      from the local cache at ~/.cache/huggingface/hub/. Zero bytes leave
+      the machine.
+
+How the model is downloaded
+~~~~~~~~~~~~~~~~~~~~~~~~~~~
+The SentenceTransformer() constructor (called in __init__ below) handles
+the download automatically on first use. It calls into the huggingface_hub
+library, which downloads the model files from:
+
+    https://huggingface.co/sentence-transformers/all-MiniLM-L6-v2
+
+The files are saved to ~/.cache/huggingface/hub/ and reused on subsequent
+runs. We pass token=False to ensure no authentication token is sent.
+
+This module is imported by both the MCP server (for search queries) and the
+indexer (for bulk embedding of code chunks). The bulk embed() call shows a
+progress bar; the single-query embed_query() does not.
+
+The env vars below must be set before importing transformers or
+sentence_transformers. They suppress warnings and progress bars that these
+libraries emit by default. Stray stderr output would interfere with the MCP
+server's JSON-RPC transport.
+"""
+
+import os
+from pathlib import Path
+
+MODEL_NAME = "all-MiniLM-L6-v2"
+
+# These env vars control logging behavior in the HuggingFace libraries.
+# They must be set before the libraries are imported.
+os.environ.setdefault("TRANSFORMERS_VERBOSITY", "error")   # suppress warnings
+os.environ.setdefault("HF_HUB_VERBOSITY", "error")        # suppress warnings
+os.environ.setdefault("HF_HUB_DISABLE_PROGRESS_BARS", "1")
+os.environ.setdefault("TOKENIZERS_PARALLELISM", "false")   # suppress threading warning
+# Disable anonymous library usage statistics (version numbers only, not user
+# data — but we disable it anyway as a matter of policy).
+os.environ.setdefault("HF_HUB_DISABLE_TELEMETRY", "1")
+
+
+def _set_offline_if_cached():
+    """If the model has already been downloaded, tell huggingface_hub to
+    skip all network calls by setting HF_HUB_OFFLINE=1.
+
+    Without this, huggingface_hub makes an HTTP request to huggingface.co
+    on every load to check if the cached model is still up to date — even
+    though the model never changes. Setting HF_HUB_OFFLINE=1 prevents this.
+
+    This must run before sentence_transformers is imported, because the
+    library reads the env var at import time.
+    """
+    # HuggingFace caches downloaded models under ~/.cache/huggingface/hub/
+    # in directories named like "models--sentence-transformers--all-MiniLM-L6-v2".
+    # The HF_HOME env var can override the base cache location.
+    hf_home = os.environ.get("HF_HOME")
+    if hf_home:
+        cache_dir = Path(hf_home) / "hub"
+    else:
+        cache_dir = Path.home() / ".cache" / "huggingface" / "hub"
+
+    model_dir = cache_dir / f"models--sentence-transformers--{MODEL_NAME}"
+    if model_dir.exists():
+        os.environ.setdefault("HF_HUB_OFFLINE", "1")
+
+
+_set_offline_if_cached()
+
+# This import must come after the env vars above are set, because the
+# transformers library reads them at import time.
+import transformers
+transformers.logging.disable_progress_bar()
+
+
+class EmbeddingProvider:
+    """Sentence-transformers embedder using all-MiniLM-L6-v2."""
+
+    def __init__(self, model_name: str = MODEL_NAME):
+        from sentence_transformers import SentenceTransformer
+
+        # This constructor loads the model from the local cache. If the model
+        # has not been downloaded yet, it downloads it from huggingface.co
+        # (~80MB, one-time). token=False ensures no auth token is sent.
+        self.model = SentenceTransformer(model_name, token=False)
+        self.dim = self.model.get_sentence_embedding_dimension()
+
+    def embed(self, texts: list[str]) -> list[list[float]]:
+        """Embed a list of texts into vectors."""
+        embeddings = self.model.encode(texts, show_progress_bar=True,
+                                       batch_size=64)
+        return embeddings.tolist()
+
+    def embed_query(self, text: str) -> list[float]:
+        """Embed a single query text."""
+        return self.model.encode([text])[0].tolist()

.claude/tools/rag/indexer.py

@@ -0,0 +1,136 @@
+#!/usr/bin/env python3
+"""Build the RAG vector index for the OpenMC codebase.
+
+This is the index-building half of the RAG pipeline. All operations are local
+once the embedding model has been downloaded and cached (see embeddings.py for
+details on model download, caching, and network behavior). It walks the repo,
+chunks every
+C++/Python/RST file (via chunker.py), embeds all chunks into 384-dim vectors
+(via embeddings.py), and stores them in a local LanceDB database on disk. The
+result is a .claude/cache/rag_index/ directory containing two tables — "code"
+and "docs" — that openmc_search.py queries at search time.
+
+Building the full index takes ~5 minutes on a 10-core machine. The bottleneck
+is the embedding step (running all chunks through the MiniLM model on CPU).
+
+Can be run standalone:  python indexer.py
+Or called programmatically:  from indexer import build_index; build_index()
+The MCP server (openmc_mcp_server.py) uses the latter when the agent calls
+openmc_rag_rebuild.
+"""
+
+import lancedb
+import sys
+import time
+from pathlib import Path
+
+# This file lives at .claude/tools/rag/indexer.py. The sys.path insert lets
+# us import sibling modules (embeddings, chunker) when run as a standalone
+# script. When imported from the MCP server, the server has already done this.
+TOOLS_DIR = Path(__file__).resolve().parent.parent
+sys.path.insert(0, str(TOOLS_DIR / "rag"))
+
+from embeddings import EmbeddingProvider
+from chunker import chunk_file
+
+
+OPENMC_ROOT = Path(__file__).resolve().parents[3]
+CACHE_DIR = OPENMC_ROOT / ".claude" / "cache"
+INDEX_DIR = CACHE_DIR / "rag_index"
+
+CODE_PATTERNS = [
+    "src/**/*.cpp",
+    "include/openmc/**/*.h",
+    "openmc/**/*.py",
+    "tests/**/*.py",
+    "examples/**/*.py",
+]
+
+DOC_PATTERNS = [
+    "docs/**/*.rst",
+]
+
+
+def collect_chunks(patterns, openmc_root):
+    """Collect all chunks from files matching the given patterns."""
+    chunks = []
+    for pattern in patterns:
+        for filepath in sorted(openmc_root.glob(pattern)):
+            if "__pycache__" in str(filepath):
+                continue
+            file_chunks = chunk_file(filepath, openmc_root)
+            chunks.extend(file_chunks)
+    return chunks
+
+
+def build_index():
+    """Build or rebuild the complete vector index."""
+    start = time.time()
+
+    # Collect all chunks
+    print("Collecting code chunks...")
+    code_chunks = collect_chunks(CODE_PATTERNS, OPENMC_ROOT)
+    print(f"  {len(code_chunks)} code chunks")
+
+    print("Collecting doc chunks...")
+    doc_chunks = collect_chunks(DOC_PATTERNS, OPENMC_ROOT)
+    print(f"  {len(doc_chunks)} doc chunks")
+
+    all_chunks = code_chunks + doc_chunks
+    if not all_chunks:
+        print("ERROR: No chunks collected!", file=sys.stderr)
+        sys.exit(1)
+
+    # Create embeddings
+    all_texts = [c["text"] for c in all_chunks]
+    print("Creating embedding provider...")
+    embedder = EmbeddingProvider()
+    print(f"  dim={embedder.dim}")
+
+    print("Embedding chunks...")
+    all_embeddings = embedder.embed(all_texts)
+
+    # Build LanceDB tables
+    INDEX_DIR.mkdir(parents=True, exist_ok=True)
+    db = lancedb.connect(str(INDEX_DIR))
+
+    # Separate code vs doc records by index (code_chunks come first in all_chunks)
+    n_code = len(code_chunks)
+    code_records = []
+    doc_records = []
+    for i, (chunk, emb) in enumerate(zip(all_chunks, all_embeddings)):
+        record = {
+            "text": chunk["text"],
+            "filepath": chunk["filepath"],
+            "kind": chunk["kind"],
+            "symbol": chunk.get("symbol", ""),
+            "start_line": chunk.get("start_line", 0),
+            "end_line": chunk.get("end_line", 0),
+            "vector": emb,
+        }
+        if i < n_code:
+            code_records.append(record)
+        else:
+            doc_records.append(record)
+
+    # Create tables (drop existing)
+    result = db.table_names() if hasattr(db, "table_names") else db.list_tables()
+    existing = result.tables if hasattr(result, "tables") else list(result)
+    for table_name in ("code", "docs"):
+        if table_name in existing:
+            db.drop_table(table_name)
+
+    if code_records:
+        db.create_table("code", code_records)
+        print(f"  Created 'code' table: {len(code_records)} rows")
+
+    if doc_records:
+        db.create_table("docs", doc_records)
+        print(f"  Created 'docs' table: {len(doc_records)} rows")
+
+    elapsed = time.time() - start
+    print(f"Done in {elapsed:.1f}s")
+
+
+if __name__ == "__main__":
+    build_index()

.claude/tools/rag/openmc_search.py

@@ -0,0 +1,202 @@
+#!/usr/bin/env python3
+"""Query the RAG vector index to find semantically related code and docs.
+
+This is the query-time half of the RAG pipeline (the counterpart to indexer.py,
+which builds the index). All operations are local — no network calls are made
+once the embedding model has been downloaded (see embeddings.py for details on
+model download and caching). Given a natural-language query, it embeds the query
+with the same MiniLM model
+used at index time, then finds the closest chunks in the local LanceDB vector
+database by cosine similarity.
+
+The core functions (get_db_and_embedder, search_table, format_results,
+search_related) are imported by the MCP server for tool calls. The script
+can also be run standalone from the command line.
+
+The "related file" mode works differently from a text query: it reads the
+target file's chunks from the index, combines them into a synthetic query
+vector, and searches for the nearest chunks from *other* files. This surfaces
+files that are semantically similar to the target file.
+
+Usage:
+    openmc_search.py "query"                    # Search code (default)
+    openmc_search.py "query" --docs             # Search documentation
+    openmc_search.py "query" --all              # Search both code and docs
+    openmc_search.py --related src/particle.cpp # Find related code
+    openmc_search.py "query" --top-k 20         # Return more results
+"""
+
+import argparse
+import sys
+from pathlib import Path
+
+# Same sys.path setup as indexer.py — needed for standalone CLI use.
+TOOLS_DIR = Path(__file__).resolve().parent.parent
+sys.path.insert(0, str(TOOLS_DIR / "rag"))
+
+OPENMC_ROOT = Path(__file__).resolve().parents[3]
+CACHE_DIR = OPENMC_ROOT / ".claude" / "cache"
+INDEX_DIR = CACHE_DIR / "rag_index"
+
+
+def get_db_and_embedder():
+    """Load the LanceDB database and embedding provider."""
+    import lancedb
+    from embeddings import EmbeddingProvider
+
+    if not INDEX_DIR.exists():
+        raise FileNotFoundError(
+            "No RAG index found. Call openmc_rag_rebuild() to build one."
+        )
+
+    db = lancedb.connect(str(INDEX_DIR))
+
+    embedder = EmbeddingProvider()
+    return db, embedder
+
+
+def _table_names(db):
+    """Return table names as a list, compatible with multiple LanceDB versions."""
+    result = db.table_names() if hasattr(db, "table_names") else db.list_tables()
+    return result.tables if hasattr(result, "tables") else list(result)
+
+
+def search_table(db, embedder, table_name, query, top_k):
+    """Search a LanceDB table with a text query."""
+    if table_name not in _table_names(db):
+        print(f"Table '{table_name}' not found in index.", file=sys.stderr)
+        return []
+
+    table = db.open_table(table_name)
+    query_vec = embedder.embed_query(query)
+    results = table.search(query_vec).limit(top_k).to_list()
+    return results
+
+
+def format_results(results, label=""):
+    """Format search results for display."""
+    if not results:
+        return "No results found.\n"
+
+    output = []
+    if label:
+        output.append(f"=== {label} ===\n")
+
+    for i, r in enumerate(results, 1):
+        filepath = r["filepath"]
+        start = r["start_line"]
+        end = r["end_line"]
+        kind = r["kind"]
+        dist = r.get("_distance", 0)
+
+        header = f"[{i}] {filepath}:{start}-{end} ({kind}, dist={dist:.3f})"
+        output.append(header)
+
+        # Show text preview (first 500 chars)
+        text = r["text"][:500]
+        if len(r["text"]) > 500:
+            text += "\n  ..."
+        # Indent the text
+        for line in text.split("\n"):
+            output.append(f"  {line}")
+        output.append("")
+
+    return "\n".join(output)
+
+
+def search_related(db, embedder, filepath, top_k):
+    """Find code related to a given file."""
+    if "code" not in _table_names(db):
+        print("No 'code' table in index.", file=sys.stderr)
+        return []
+
+    table = db.open_table("code")
+
+    # Normalize filepath
+    fp = filepath
+    if Path(filepath).is_absolute():
+        try:
+            fp = str(Path(filepath).relative_to(OPENMC_ROOT))
+        except ValueError:
+            pass
+
+    # Get chunks from target file
+    try:
+        safe_fp = fp.replace("'", "''")
+        target_chunks = table.search().where(
+            f"filepath = '{safe_fp}'"
+        ).limit(50).to_list()
+    except Exception:
+        # LanceDB where clause might not work in all versions
+        # Fall back to fetching all and filtering
+        all_data = table.to_pandas()
+        target_rows = all_data[all_data["filepath"] == fp]
+        if target_rows.empty:
+            print(f"No chunks found for '{fp}'", file=sys.stderr)
+            return []
+        target_chunks = target_rows.head(50).to_dict("records")
+
+    if not target_chunks:
+        print(f"No chunks found for '{fp}'", file=sys.stderr)
+        return []
+
+    # Combine top chunks as the query
+    combined_text = " ".join(c["text"][:200] for c in target_chunks[:5])
+    query_vec = embedder.embed_query(combined_text)
+
+    # Search excluding the source file
+    results = table.search(query_vec).limit(top_k + 10).to_list()
+    # Filter out same file
+    results = [r for r in results if r["filepath"] != fp][:top_k]
+    return results
+
+
+def main():
+    parser = argparse.ArgumentParser(
+        description="Semantic search across OpenMC codebase and docs",
+        formatter_class=argparse.RawDescriptionHelpFormatter,
+        epilog="""examples:
+  %(prog)s "particle random number seed initialization"
+  %(prog)s "how to define tallies" --docs
+  %(prog)s "weight window variance reduction" --all
+  %(prog)s "where is cross section data loaded" --top-k 15
+  %(prog)s --related src/simulation.cpp
+  %(prog)s --related src/particle_restart.cpp --top-k 5""",
+    )
+    parser.add_argument("query", nargs="?", help="Search query")
+    parser.add_argument("--docs", action="store_true",
+                        help="Search documentation instead of code")
+    parser.add_argument("--all", action="store_true",
+                        help="Search both code and documentation")
+    parser.add_argument("--related", metavar="FILE",
+                        help="Find code related to a given file")
+    parser.add_argument("--top-k", type=int, default=10,
+                        help="Number of results (default: 10)")
+    args = parser.parse_args()
+
+    if not args.query and not args.related:
+        parser.print_help()
+        sys.exit(1)
+
+    db, embedder = get_db_and_embedder()
+
+    if args.related:
+        results = search_related(db, embedder, args.related, args.top_k)
+        print(format_results(results, f"Code related to {args.related}"))
+    elif args.all:
+        code_results = search_table(
+            db, embedder, "code", args.query, args.top_k)
+        doc_results = search_table(
+            db, embedder, "docs", args.query, args.top_k)
+        print(format_results(code_results, "Code"))
+        print(format_results(doc_results, "Documentation"))
+    elif args.docs:
+        results = search_table(db, embedder, "docs", args.query, args.top_k)
+        print(format_results(results, "Documentation"))
+    else:
+        results = search_table(db, embedder, "code", args.query, args.top_k)
+        print(format_results(results, "Code"))
+
+
+if __name__ == "__main__":
+    main()

.claude/tools/requirements.txt

@@ -0,0 +1,8 @@
+# MCP server
+mcp>=1.0.0
+
+# Vector database
+lancedb>=0.15.0
+
+# Embeddings (local, no API key)
+sentence-transformers>=2.7.0

.claude/tools/start_server.sh

@@ -0,0 +1,34 @@
+#!/bin/bash
+# Bootstrap the Python venv (if needed) and start the OpenMC MCP server.
+set -e
+
+SCRIPT_DIR="$(cd "$(dirname "$0")" && pwd)"
+CACHE_DIR="$(dirname "$SCRIPT_DIR")/cache"
+VENV_DIR="$CACHE_DIR/.venv"
+SENTINEL="$VENV_DIR/.installed"
+
+if ! command -v python3 >/dev/null 2>&1; then
+    echo "Error: python3 not found on PATH." >&2
+    exit 1
+fi
+
+if ! python3 -c 'import sys; assert sys.version_info >= (3,12)' 2>/dev/null; then
+    echo "Error: Python 3.12+ is required." >&2
+    exit 1
+fi
+
+if [ ! -f "$SENTINEL" ]; then
+    rm -rf "$VENV_DIR"
+    mkdir -p "$CACHE_DIR"
+    python3 -m venv "$VENV_DIR"
+
+    if ! "$VENV_DIR/bin/pip" install -q -r "$SCRIPT_DIR/requirements.txt"; then
+        echo "Error: pip install failed. Remove $VENV_DIR and retry." >&2
+        rm -rf "$VENV_DIR"
+        exit 1
+    fi
+
+    touch "$SENTINEL"
+fi
+
+exec "$VENV_DIR/bin/python" "$SCRIPT_DIR/openmc_mcp_server.py"

.git_archival.txt

@@ -0,0 +1,3 @@
+commit: $Format:%H$
+commit-date: $Format:%cI$
+describe-name: $Format:%(describe:tags=true,match=*[0-9]*)$

.gitattributes

@@ -0,0 +1 @@
+.git_archival.txt  export-subst

.github/ISSUE_TEMPLATE/feature_request.md

@@ -1,5 +1,5 @@
 ---
-name: Feature request
+name: Feature or enhancement request
 about: Suggest a new feature or enhancement to existing capabilities
 title: ''
 labels: ''

.github/agents/Review.agent.md

@@ -0,0 +1,8 @@
+---
+name: Review
+description: Reviews code changes on the current branch, evaluating them against OpenMC's contribution criteria and providing structured feedback.
+argument-hint: Optionally provide a focus area (e.g., "focus on physics correctness", "check Python API design"). If omitted, a full review is performed.
+---
+You are an expert code reviewer for OpenMC. Use the `reviewing-openmc-code` skill to perform a structured review of the code changes on the current branch.
+
+If the user provides a focus area, prioritize that section of the review.

.github/copilot-instructions.md

@@ -0,0 +1 @@
+When reviewing code changes in this repository, use the `reviewing-openmc-code` skill.

.github/pull_request_template.md

@@ -13,7 +13,7 @@ Fixes # (issue)
 # Checklist
 
 - [ ] I have performed a self-review of my own code
-- [ ] I have run [clang-format](https://docs.openmc.org/en/latest/devguide/styleguide.html#automatic-formatting) (version 15) on any C++ source files (if applicable)
+- [ ] I have run [clang-format](https://docs.openmc.org/en/latest/devguide/styleguide.html#automatic-formatting) (version 18) on any C++ source files (if applicable)
 - [ ] I have followed the [style guidelines](https://docs.openmc.org/en/latest/devguide/styleguide.html#python) for Python source files (if applicable)
 - [ ] I have made corresponding changes to the documentation (if applicable)
 - [ ] I have added tests that prove my fix is effective or that my feature works (if applicable)

.github/workflows/ci.yml

@@ -1,4 +1,4 @@
-name: CI
+name: Tests and Coverage
 
 on:
   # allows us to run workflows manually
@@ -21,102 +21,127 @@ env:
   GITHUB_TOKEN: ${{ secrets.GITHUB_TOKEN }}
 
 jobs:
+  filter-changes:
+    runs-on: ubuntu-latest
+    outputs:
+      source_changed: ${{ steps.filter.outputs.source_changed }}
+    steps:
+      - name: Check out the repository
+        uses: actions/checkout@v6
+      - name: Examine changed files
+        id: filter
+        uses: dorny/paths-filter@v4
+        with:
+          filters: |
+            source_changed:
+              - '!docs/**'
+              - '!**/*.md'
+          predicate-quantifier: 'every'
   main:
+    needs: filter-changes
+    if: ${{ needs.filter-changes.outputs.source_changed == 'true' }}
     runs-on: ubuntu-22.04
     strategy:
       matrix:
-        python-version: ["3.10"]
+        python-version: ["3.12"]
         mpi: [n, y]
         omp: [n, y]
         dagmc: [n]
-        ncrystal: [n]
         libmesh: [n]
         event: [n]
-        vectfit: [n]
 
         include:
-          - python-version: "3.7"
+          - python-version: "3.13"
             omp: n
             mpi: n
-          - python-version: "3.8"
+          - python-version: "3.14"
             omp: n
             mpi: n
-          - python-version: "3.9"
-            omp: n
-            mpi: n
-          - python-version: "3.11"
+          - python-version: "3.14t"
             omp: n
             mpi: n
           - dagmc: y
-            python-version: "3.10"
-            mpi: y
-            omp: y
-          - ncrystal: y
-            python-version: "3.10"
-            mpi: n
-            omp: n
-          - libmesh: y
-            python-version: "3.10"
+            python-version: "3.12"
             mpi: y
             omp: y
           - libmesh: y
-            python-version: "3.10"
+            python-version: "3.12"
+            mpi: y
+            omp: y
+          - libmesh: y
+            python-version: "3.12"
             mpi: n
             omp: y
           - event: y
-            python-version: "3.10"
+            python-version: "3.12"
             omp: y
             mpi: n
-          - vectfit: y
-            python-version: "3.10"
-            omp: n
-            mpi: y
     name: "Python ${{ matrix.python-version }} (omp=${{ matrix.omp }},
-      mpi=${{ matrix.mpi }}, dagmc=${{ matrix.dagmc }}, ncrystal=${{ matrix.ncrystal }},
-      libmesh=${{ matrix.libmesh }}, event=${{ matrix.event }}
-      vectfit=${{ matrix.vectfit }})"
+      mpi=${{ matrix.mpi }}, dagmc=${{ matrix.dagmc }},
+      libmesh=${{ matrix.libmesh }}, event=${{ matrix.event }}"
 
     env:
       MPI: ${{ matrix.mpi }}
       PHDF5: ${{ matrix.mpi }}
       OMP: ${{ matrix.omp }}
       DAGMC: ${{ matrix.dagmc }}
-      NCRYSTAL: ${{ matrix.ncrystal }}
       EVENT: ${{ matrix.event }}
-      VECTFIT: ${{ matrix.vectfit }}
       LIBMESH: ${{ matrix.libmesh }}
       NPY_DISABLE_CPU_FEATURES: "AVX512F AVX512_SKX"
       OPENBLAS_NUM_THREADS: 1
+      PYTEST_ADDOPTS: --cov=openmc --cov-report=lcov:coverage-python.lcov
       # libfabric complains about fork() as a result of using Python multiprocessing.
       # We can work around it with RDMAV_FORK_SAFE=1 in libfabric < 1.13 and with
       # FI_EFA_FORK_SAFE=1 in more recent versions.
       RDMAV_FORK_SAFE: 1
 
     steps:
-      - uses: actions/checkout@v3
+      - name: Setup cmake
+        uses: jwlawson/actions-setup-cmake@v2
+        with:
+          cmake-version: '3.31'
+
+      - name: Checkout repository
+        uses: actions/checkout@v6
+        with:
+          fetch-depth: 0
 
       - name: Set up Python ${{ matrix.python-version }}
-        uses: actions/setup-python@v4
+        uses: actions/setup-python@v6
         with:
           python-version: ${{ matrix.python-version }}
 
       - name: Environment Variables
         run: |
-          echo "DAGMC_ROOT=$HOME/DAGMC"
           echo "OPENMC_CROSS_SECTIONS=$HOME/nndc_hdf5/cross_sections.xml" >> $GITHUB_ENV
           echo "OPENMC_ENDF_DATA=$HOME/endf-b-vii.1" >> $GITHUB_ENV
+          # get the sha of the last branch commit
+          # for push and workflow_dispatch events, use the current reference head
+          BRANCH_SHA=HEAD
+          # for a pull_request event, use the last reference of the parents of the merge commit
+          if [ "${{ github.event_name }}" == "pull_request" ]; then
+            BRANCH_SHA=$(git rev-list --parents -n 1 HEAD | rev | cut -d" " -f 1 | rev)
+          fi
+          COMMIT_MESSAGE=$(git log $BRANCH_SHA -1 --pretty=%B | tr '\n' ' ')
+          echo ${COMMIT_MESSAGE}
+          echo "COMMIT_MESSAGE=${COMMIT_MESSAGE}" >> $GITHUB_ENV
 
       - name: Apt dependencies
         shell: bash
         run: |
           sudo apt -y update
           sudo apt install -y libpng-dev \
-                              libmpich-dev \
                               libnetcdf-dev \
                               libpnetcdf-dev \
                               libhdf5-serial-dev \
-                              libhdf5-mpich-dev \
                               libeigen3-dev
+
+      - name: Optional apt dependencies for MPI
+        shell: bash
+        if: ${{ matrix.mpi == 'y' }}
+        run: |
+          sudo apt install -y libhdf5-mpich-dev \
+                              libmpich-dev
           sudo update-alternatives --set mpi /usr/bin/mpicc.mpich
           sudo update-alternatives --set mpirun /usr/bin/mpirun.mpich
           sudo update-alternatives --set mpi-x86_64-linux-gnu /usr/include/x86_64-linux-gnu/mpich
@@ -127,13 +152,18 @@ jobs:
           echo "$HOME/NJOY2016/build" >> $GITHUB_PATH
           $GITHUB_WORKSPACE/tools/ci/gha-install.sh
 
+      - name: display-config
+        shell: bash
+        run: |
+          openmc -v
+
       - name: cache-xs
-        uses: actions/cache@v3
+        uses: actions/cache@v5
         with:
           path: |
             ~/nndc_hdf5
             ~/endf-b-vii.1
-          key: ${{ runner.os }}-build-xs-cache
+          key: ${{ runner.os }}-build-xs-cache-${{ hashFiles(format('{0}/tools/ci/download-xs.sh', github.workspace)) }}
 
       - name: before
         shell: bash
@@ -145,18 +175,74 @@ jobs:
           CTEST_OUTPUT_ON_FAILURE=1 make test -C $GITHUB_WORKSPACE/build/
           $GITHUB_WORKSPACE/tools/ci/gha-script.sh
 
-      - name: after_success
+      - name: Setup tmate debug session
+        continue-on-error: true
+        if: ${{ failure() && contains(env.COMMIT_MESSAGE, '[gha-debug]') }}
+        uses: mxschmitt/action-tmate@v3
+        timeout-minutes: 10
+
+      - name: Generate C++ coverage (gcovr)
         shell: bash
         run: |
-          cpp-coveralls -i src -i include -e src/external --exclude-pattern "/usr/*" --dump cpp_cov.json
-          coveralls --merge=cpp_cov.json --service=github
+          # Produce LCOV directly from gcov data in the build tree
+          gcovr \
+            --root "$GITHUB_WORKSPACE" \
+            --object-directory "$GITHUB_WORKSPACE/build" \
+            --filter "$GITHUB_WORKSPACE/src" \
+            --filter "$GITHUB_WORKSPACE/include" \
+            --exclude "$GITHUB_WORKSPACE/src/external/.*" \
+            --exclude "$GITHUB_WORKSPACE/src/include/openmc/external/.*" \
+            --gcov-ignore-errors source_not_found \
+            --gcov-ignore-errors output_error \
+            --gcov-ignore-parse-errors suspicious_hits.warn \
+            --merge-mode-functions=separate \
+            --print-summary \
+            --lcov -o coverage-cpp.lcov || true
 
-  finish:
-    needs: main
+      - name: Merge C++ and Python coverage
+        shell: bash
+        run: |
+          # Merge C++ and Python LCOV into a single file for upload
+          cat coverage-cpp.lcov coverage-python.lcov > coverage.lcov
+
+      - name: Upload coverage to Coveralls
+        if: ${{ hashFiles('coverage.lcov') != '' }}
+        uses: coverallsapp/github-action@v2
+        with:
+          github-token: ${{ secrets.GITHUB_TOKEN }}
+          parallel: true
+          flag-name: C++ and Python
+          path-to-lcov: coverage.lcov
+          fail-on-error: false
+
+  coverage:
+    needs: [filter-changes, main]
+    if: ${{ always() }}
     runs-on: ubuntu-latest
     steps:
       - name: Coveralls Finished
-        uses: coverallsapp/github-action@master
+        if: ${{ needs.filter-changes.outputs.source_changed == 'true' }}
+        uses: coverallsapp/github-action@v2
         with:
-          github-token: ${{ secrets.github_token }}
+          github-token: ${{ secrets.GITHUB_TOKEN }}
           parallel-finished: true
+          fail-on-error: false
+
+  ci-pass:
+    needs: [filter-changes, main, coverage]
+    name: Check CI status
+    if: ${{ always() }}
+    runs-on: ubuntu-latest
+    steps:
+      - name: Check CI status
+        run: |
+          if [[ "${{ needs.filter-changes.outputs.source_changed }}" == "false" ]]; then
+            echo "Documentation-only change - CI skipped successfully"
+            exit 0
+          fi
+          if [[ "${{ needs.main.result }}" == "success" && "${{ needs.coverage.result }}" == "success" ]]; then
+            echo "CI passed"
+            exit 0
+          fi
+          echo "CI failed"
+          exit 1

.github/workflows/dockerhub-publish-dagmc-libmesh.yml

@@ -2,7 +2,8 @@ name: dockerhub-publish-latest-dagmc-libmesh
 
 on:
   push:
-    branches: master
+    branches:
+      - master
 
 jobs:
   main:
@@ -10,20 +11,20 @@ jobs:
     steps:
       -
         name: Set up QEMU
-        uses: docker/setup-qemu-action@v2
+        uses: docker/setup-qemu-action@v3
       -
         name: Set up Docker Buildx
-        uses: docker/setup-buildx-action@v2
+        uses: docker/setup-buildx-action@v3
       -
         name: Login to DockerHub
-        uses: docker/login-action@v2
+        uses: docker/login-action@v3
         with:
           username: ${{ secrets.DOCKERHUB_USERNAME }}
           password: ${{ secrets.DOCKERHUB_TOKEN }}
       -
         name: Build and push
         id: docker_build
-        uses: docker/build-push-action@v4
+        uses: docker/build-push-action@v5
         with:
           push: true
           tags: openmc/openmc:latest-dagmc-libmesh

.github/workflows/dockerhub-publish-dagmc.yml

@@ -2,7 +2,8 @@ name: dockerhub-publish-latest-dagmc
 
 on:
   push:
-    branches: master
+    branches:
+      - master
 
 jobs:
   main:
@@ -10,20 +11,20 @@ jobs:
     steps:
       -
         name: Set up QEMU
-        uses: docker/setup-qemu-action@v2
+        uses: docker/setup-qemu-action@v3
       -
         name: Set up Docker Buildx
-        uses: docker/setup-buildx-action@v2
+        uses: docker/setup-buildx-action@v3
       -
         name: Login to DockerHub
-        uses: docker/login-action@v2 
+        uses: docker/login-action@v3
         with:
           username: ${{ secrets.DOCKERHUB_USERNAME }}
           password: ${{ secrets.DOCKERHUB_TOKEN }}
       -
         name: Build and push
         id: docker_build
-        uses: docker/build-push-action@v4
+        uses: docker/build-push-action@v5
         with:
           push: true
           tags: openmc/openmc:latest-dagmc

.github/workflows/dockerhub-publish-dev.yml

@@ -2,7 +2,8 @@ name: dockerhub-publish-develop
 
 on:
   push:
-    branches: develop
+    branches:
+      - develop
 
 jobs:
   main:
@@ -10,20 +11,20 @@ jobs:
     steps:
       -
         name: Set up QEMU
-        uses: docker/setup-qemu-action@v2
+        uses: docker/setup-qemu-action@v3
       -
         name: Set up Docker Buildx
-        uses: docker/setup-buildx-action@v2
+        uses: docker/setup-buildx-action@v3
       -
         name: Login to DockerHub
-        uses: docker/login-action@v2 
+        uses: docker/login-action@v3
         with:
           username: ${{ secrets.DOCKERHUB_USERNAME }}
           password: ${{ secrets.DOCKERHUB_TOKEN }}
       -
         name: Build and push
         id: docker_build
-        uses: docker/build-push-action@v4
+        uses: docker/build-push-action@v5
         with:
           push: true
           tags: openmc/openmc:develop

.github/workflows/dockerhub-publish-develop-dagmc-libmesh.yml

@@ -2,7 +2,8 @@ name: dockerhub-publish-develop-dagmc-libmesh
 
 on:
   push:
-    branches: develop
+    branches:
+      - develop
 
 jobs:
   main:
@@ -10,20 +11,20 @@ jobs:
     steps:
       -
         name: Set up QEMU
-        uses: docker/setup-qemu-action@v2
+        uses: docker/setup-qemu-action@v3
       -
         name: Set up Docker Buildx
-        uses: docker/setup-buildx-action@v2
+        uses: docker/setup-buildx-action@v3
       -
         name: Login to DockerHub
-        uses: docker/login-action@v2 
+        uses: docker/login-action@v3
         with:
           username: ${{ secrets.DOCKERHUB_USERNAME }}
           password: ${{ secrets.DOCKERHUB_TOKEN }}
       -
         name: Build and push
         id: docker_build
-        uses: docker/build-push-action@v4
+        uses: docker/build-push-action@v5
         with:
           push: true
           tags: openmc/openmc:develop-dagmc-libmesh

.github/workflows/dockerhub-publish-develop-dagmc.yml

@@ -2,7 +2,8 @@ name: dockerhub-publish-develop-dagmc
 
 on:
   push:
-    branches: develop
+    branches:
+      - develop
 
 jobs:
   main:
@@ -10,20 +11,20 @@ jobs:
     steps:
       -
         name: Set up QEMU
-        uses: docker/setup-qemu-action@v2
+        uses: docker/setup-qemu-action@v3
       -
         name: Set up Docker Buildx
-        uses: docker/setup-buildx-action@v2
+        uses: docker/setup-buildx-action@v3
       -
         name: Login to DockerHub
-        uses: docker/login-action@v2 
+        uses: docker/login-action@v3
         with:
           username: ${{ secrets.DOCKERHUB_USERNAME }}
           password: ${{ secrets.DOCKERHUB_TOKEN }}
       -
         name: Build and push
         id: docker_build
-        uses: docker/build-push-action@v4
+        uses: docker/build-push-action@v5
         with:
           push: true
           tags: openmc/openmc:develop-dagmc

.github/workflows/dockerhub-publish-develop-libmesh.yml

@@ -2,7 +2,8 @@ name: dockerhub-publish-develop-libmesh
 
 on:
   push:
-    branches: develop
+    branches:
+      - develop
 
 jobs:
   main:
@@ -10,20 +11,20 @@ jobs:
     steps:
       -
         name: Set up QEMU
-        uses: docker/setup-qemu-action@v2
+        uses: docker/setup-qemu-action@v3
       -
         name: Set up Docker Buildx
-        uses: docker/setup-buildx-action@v2
+        uses: docker/setup-buildx-action@v3
       -
         name: Login to DockerHub
-        uses: docker/login-action@v2 
+        uses: docker/login-action@v3
         with:
           username: ${{ secrets.DOCKERHUB_USERNAME }}
           password: ${{ secrets.DOCKERHUB_TOKEN }}
       -
         name: Build and push
         id: docker_build
-        uses: docker/build-push-action@v4
+        uses: docker/build-push-action@v5
         with:
           push: true
           tags: openmc/openmc:develop-libmesh

.github/workflows/dockerhub-publish-libmesh.yml

@@ -2,7 +2,8 @@ name: dockerhub-publish-latest-libmesh
 
 on:
   push:
-    branches: master
+    branches:
+      - master
 
 jobs:
   main:
@@ -10,20 +11,20 @@ jobs:
     steps:
       -
         name: Set up QEMU
-        uses: docker/setup-qemu-action@v2
+        uses: docker/setup-qemu-action@v3
       -
         name: Set up Docker Buildx
-        uses: docker/setup-buildx-action@v2
+        uses: docker/setup-buildx-action@v3
       -
         name: Login to DockerHub
-        uses: docker/login-action@v2 
+        uses: docker/login-action@v3
         with:
           username: ${{ secrets.DOCKERHUB_USERNAME }}
           password: ${{ secrets.DOCKERHUB_TOKEN }}
       -
         name: Build and push
         id: docker_build
-        uses: docker/build-push-action@v4
+        uses: docker/build-push-action@v5
         with:
           push: true
           tags: openmc/openmc:latest-libmesh

.github/workflows/dockerhub-publish-release-dagmc-libmesh.yml

@@ -2,31 +2,32 @@ name: dockerhub-publish-release-dagmc-libmesh
 
 on:
   push:
-    tags: 'v*.*.*'
+    tags:
+      - 'v*.*.*'
 
 jobs:
   main:
     runs-on: ubuntu-latest
     steps:
-      - uses: actions/checkout@v3
+      - uses: actions/checkout@v6
       - name: Set env
         run: echo "RELEASE_VERSION=${GITHUB_REF#refs/*/}" >> $GITHUB_ENV
       -
         name: Set up QEMU
-        uses: docker/setup-qemu-action@v2
+        uses: docker/setup-qemu-action@v3
       -
         name: Set up Docker Buildx
-        uses: docker/setup-buildx-action@v2
+        uses: docker/setup-buildx-action@v3
       -
         name: Login to DockerHub
-        uses: docker/login-action@v2
+        uses: docker/login-action@v3
         with:
           username: ${{ secrets.DOCKERHUB_USERNAME }}
           password: ${{ secrets.DOCKERHUB_TOKEN }}
       -
         name: Build and push
         id: docker_build
-        uses: docker/build-push-action@v4
+        uses: docker/build-push-action@v5
         with:
           push: true
           tags: openmc/openmc:${{ env.RELEASE_VERSION }}-dagmc-libmesh

.github/workflows/dockerhub-publish-release-dagmc.yml

@@ -2,31 +2,32 @@ name: dockerhub-publish-release-dagmc
 
 on:
   push:
-    tags: 'v*.*.*'
+    tags:
+      - 'v*.*.*'
 
 jobs:
   main:
     runs-on: ubuntu-latest
     steps:
-      - uses: actions/checkout@v3
+      - uses: actions/checkout@v6
       - name: Set env
         run: echo "RELEASE_VERSION=${GITHUB_REF#refs/*/}" >> $GITHUB_ENV
       -
         name: Set up QEMU
-        uses: docker/setup-qemu-action@v2
+        uses: docker/setup-qemu-action@v3
       -
         name: Set up Docker Buildx
-        uses: docker/setup-buildx-action@v2
+        uses: docker/setup-buildx-action@v3
       -
         name: Login to DockerHub
-        uses: docker/login-action@v2 
+        uses: docker/login-action@v3
         with:
           username: ${{ secrets.DOCKERHUB_USERNAME }}
           password: ${{ secrets.DOCKERHUB_TOKEN }}
       -
         name: Build and push
         id: docker_build
-        uses: docker/build-push-action@v4
+        uses: docker/build-push-action@v5
         with:
           push: true
           tags: openmc/openmc:${{ env.RELEASE_VERSION }}-dagmc

.github/workflows/dockerhub-publish-release-libmesh.yml

@@ -2,31 +2,32 @@ name: dockerhub-publish-release-libmesh
 
 on:
   push:
-    tags: 'v*.*.*'
+    tags:
+      - 'v*.*.*'
 
 jobs:
   main:
     runs-on: ubuntu-latest
     steps:
-      - uses: actions/checkout@v3
+      - uses: actions/checkout@v6
       - name: Set env
         run: echo "RELEASE_VERSION=${GITHUB_REF#refs/*/}" >> $GITHUB_ENV
       -
         name: Set up QEMU
-        uses: docker/setup-qemu-action@v2
+        uses: docker/setup-qemu-action@v3
       -
         name: Set up Docker Buildx
-        uses: docker/setup-buildx-action@v2
+        uses: docker/setup-buildx-action@v3
       -
         name: Login to DockerHub
-        uses: docker/login-action@v2
+        uses: docker/login-action@v3
         with:
           username: ${{ secrets.DOCKERHUB_USERNAME }}
           password: ${{ secrets.DOCKERHUB_TOKEN }}
       -
         name: Build and push
         id: docker_build
-        uses: docker/build-push-action@v4
+        uses: docker/build-push-action@v5
         with:
           push: true
           tags: openmc/openmc:${{ env.RELEASE_VERSION }}-libmesh

.github/workflows/dockerhub-publish-release.yml

@@ -2,31 +2,32 @@ name: dockerhub-publish-release
 
 on:
   push:
-    tags: 'v*.*.*'
+    tags:
+      - 'v*.*.*'
 
 jobs:
   main:
     runs-on: ubuntu-latest
     steps:
-      - uses: actions/checkout@v3
+      - uses: actions/checkout@v6
       - name: Set env
         run: echo "RELEASE_VERSION=${GITHUB_REF#refs/*/}" >> $GITHUB_ENV
       -
         name: Set up QEMU
-        uses: docker/setup-qemu-action@v2
+        uses: docker/setup-qemu-action@v3
       -
         name: Set up Docker Buildx
-        uses: docker/setup-buildx-action@v2
+        uses: docker/setup-buildx-action@v3
       -
         name: Login to DockerHub
-        uses: docker/login-action@v2 
+        uses: docker/login-action@v3
         with:
           username: ${{ secrets.DOCKERHUB_USERNAME }}
           password: ${{ secrets.DOCKERHUB_TOKEN }}
       -
         name: Build and push
         id: docker_build
-        uses: docker/build-push-action@v4
+        uses: docker/build-push-action@v5
         with:
           push: true
           tags: openmc/openmc:${{ env.RELEASE_VERSION }}

.github/workflows/dockerhub-publish.yml

@@ -2,7 +2,8 @@ name: dockerhub-publish-latest
 
 on:
   push:
-    branches: master
+    branches:
+      - master
 
 jobs:
   main:
@@ -10,20 +11,20 @@ jobs:
     steps:
       -
         name: Set up QEMU
-        uses: docker/setup-qemu-action@v2
+        uses: docker/setup-qemu-action@v3
       -
         name: Set up Docker Buildx
-        uses: docker/setup-buildx-action@v2
+        uses: docker/setup-buildx-action@v3
       -
         name: Login to DockerHub
-        uses: docker/login-action@v2 
+        uses: docker/login-action@v3
         with:
           username: ${{ secrets.DOCKERHUB_USERNAME }}
           password: ${{ secrets.DOCKERHUB_TOKEN }}
       -
         name: Build and push
         id: docker_build
-        uses: docker/build-push-action@v4
+        uses: docker/build-push-action@v5
         with:
           push: true
           tags: openmc/openmc:latest

.github/workflows/format-check.yml

@@ -1,12 +1,28 @@
 name: C++ Format Check
 
-on: pull_request
+on:
+  # allow workflow to be run manually
+  workflow_dispatch:
+
+  pull_request:
+    types:
+      - opened
+      - synchronize
+      - reopened
+      - labeled
+      - unlabeled
+    branches:
+      - develop
+      - master
 
 jobs:
   cpp-linter:
     runs-on: ubuntu-latest
+    permissions:
+      contents: read
+      pull-requests: write
     steps:
-      - uses: actions/checkout@v3
+      - uses: actions/checkout@v6
       - uses: cpp-linter/cpp-linter-action@v2
         id: linter
         env:
@@ -15,11 +31,30 @@ jobs:
           style: file
           files-changed-only: true
           tidy-checks: '-*'
-          version: '15' # clang-format version
+          version: '18' # clang-format version
+          format-review: ${{ github.event_name == 'pull_request' && contains(github.event.pull_request.labels.*.name, 'cpp-format-suggest') }}
+          passive-reviews: ${{ github.event_name == 'pull_request' && contains(github.event.pull_request.labels.*.name, 'cpp-format-suggest') }}
           file-annotations: true
           step-summary: true
           extensions: 'cpp,h'
 
+      - name: Comment with suggestion instructions
+        if: steps.linter.outputs.checks-failed > 0 && !contains(github.event.pull_request.labels.*.name, 'cpp-format-suggest')
+        uses: actions/github-script@v7
+        with:
+          script: |
+            const {owner, repo} = context.repo;
+            const issue_number = context.payload.pull_request.number;
+            await github.rest.issues.createComment({
+              owner,
+              repo,
+              issue_number,
+              body: "C++ formatting checks failed. Add the `cpp-format-suggest` label to this PR for inline formatting suggestions on the next run."
+            });
+
       - name: Failure Check
         if: steps.linter.outputs.checks-failed > 0
-        run: echo "Some files failed the formatting check! See job summary and file annotations for more info" && exit 1
+        run: |
+          echo "Some files failed the formatting check."
+          echo "See job summary and file annotations for details."
+          exit 1

.gitignore

@@ -25,12 +25,13 @@ examples/**/*.xml
 
 # Documentation builds
 docs/build
+docs/doxygen/xml
 docs/source/_images/*.pdf
 docs/source/_images/*.aux
 docs/source/pythonapi/generated/
 
 # Source build
-build
+build*/
 
 # build from src/utils/setup.py
 src/utils/build
@@ -104,5 +105,8 @@ CMakeSettings.json
 # Visual Studio Code configuration files
 .vscode/
 
+# Claude Code agent tools (cached/generated artifacts)
+.claude/cache/
+
 # Python pickle files
 *.pkl

.gitmodules

@@ -1,15 +1,6 @@
 [submodule "vendor/pugixml"]
 	path = vendor/pugixml
 	url = https://github.com/zeux/pugixml.git
-[submodule "vendor/gsl-lite"]
-	path = vendor/gsl-lite
-	url = https://github.com/martinmoene/gsl-lite.git
-[submodule "vendor/xtensor"]
-	path = vendor/xtensor
-	url = https://github.com/xtensor-stack/xtensor.git
-[submodule "vendor/xtl"]
-	path = vendor/xtl
-	url = https://github.com/xtensor-stack/xtl.git
 [submodule "vendor/fmt"]
 	path = vendor/fmt
 	url = https://github.com/fmtlib/fmt.git

.mcp.json

@@ -0,0 +1,9 @@
+{
+  "mcpServers": {
+    "openmc-code-tools": {
+      "type": "stdio",
+      "command": "bash",
+      "args": [".claude/tools/start_server.sh"]
+    }
+  }
+}

.readthedocs.yaml

@@ -1,13 +1,23 @@
 version: 2
 
 build:
-  os: "ubuntu-20.04"
+  os: "ubuntu-24.04"
   tools:
-    python: "3.9"
+    python: "3.12"
+  jobs:
+    post_checkout:
+      - git fetch --unshallow || true
+      - cd docs/doxygen && doxygen && cd -
 
 sphinx:
   configuration: docs/source/conf.py
 
+formats:
+  - pdf
+
 python:
   install:
-    - requirements: docs/requirements-rtd.txt
+    - method: pip
+      path: .
+      extra_requirements:
+        - docs

AGENTS.md

@@ -0,0 +1,348 @@
+# OpenMC AI Coding Agent Instructions
+
+## Project Overview
+
+OpenMC is a Monte Carlo particle transport code for simulating nuclear reactors,
+fusion devices, or other systems with neutron/photon radiation. It's a hybrid
+C++17/Python codebase where:
+- **C++ core** (`src/`, `include/openmc/`) handles the computationally intensive transport simulation
+- **Python API** (`openmc/`) provides user-facing model building, post-processing, and depletion capabilities
+- **C API bindings** (`openmc/lib/`) wrap the C++ library via ctypes for runtime control
+
+## Architecture & Key Components
+
+### C++ Component Structure
+- **Global vectors of unique_ptrs**: Core objects like `model::cells`, `model::universes`, `nuclides` are stored as `vector<unique_ptr<T>>` in nested namespaces (`openmc::model`, `openmc::simulation`, `openmc::settings`, `openmc::data`)
+- **Custom container types**: OpenMC provides its own `vector`, `array`, `unique_ptr`, and `make_unique` in the `openmc::` namespace (defined in `vector.h`, `array.h`, `memory.h`). These are currently typedefs to `std::` equivalents but may become custom implementations for accelerator support. Always use `openmc::vector`, not `std::vector`.
+- **Geometry systems**:
+  - **CSG (default)**: Arbitrarily complex Constructive Solid Geometry using `Surface`, `Region`, `Cell`, `Universe`, `Lattice`
+  - **DAGMC**: CAD-based geometry via Direct Accelerated Geometry Monte Carlo (optional, requires `OPENMC_USE_DAGMC`)
+  - **Unstructured mesh**: libMesh-based geometry (optional, requires `OPENMC_USE_LIBMESH`)
+- **Particle tracking**: `Particle` class with `GeometryState` manages particle transport through geometry
+- **Tallies**: Score quantities during simulation via `Filter` and `Tally` objects
+- **Random ray solver**: Alternative deterministic method in `src/random_ray/`
+- **Optional features**: DAGMC (CAD geometry), libMesh (unstructured mesh), MPI, all controlled by `#ifdef OPENMC_MPI`, etc.
+
+### Python Component Structure
+- **ID management**: All geometry objects (Cell, Surface, Material, etc.) inherit from `IDManagerMixin` which auto-assigns unique integer IDs and tracks them via class-level `used_ids` and `next_id`
+- **Input validation**: Extensive use of `openmc.checkvalue` module functions (`check_type`, `check_value`, `check_length`) for all setters
+- **XML I/O**: Most classes implement `to_xml_element()` and `from_xml_element()` for serialization to OpenMC's XML input format
+- **HDF5 output**: Post-simulation data in statepoint files read via `openmc.StatePoint`
+- **Depletion**: `openmc.deplete` implements burnup via operator-splitting with various integrators (Predictor, CECM, etc.)
+- **Nuclear Data**: `openmc.data` provides programmatic access to nuclear data files (ENDF, ACE, HDF5)
+
+## Git Branching Workflow
+
+OpenMC uses a git flow branching model with two primary branches:
+
+- **`develop` branch**: The main development branch where all ongoing development takes place. This is the **primary branch against which pull requests are submitted and merged**. This branch is not guaranteed to be stable and may contain work-in-progress features.
+- **`master` branch**: The stable release branch containing the latest stable release of OpenMC. This branch only receives merges from `develop` when the development team decides a release should occur.
+
+### Instructions for Code Review
+
+When reviewing code changes in this repository, use the `reviewing-openmc-code` skill.
+
+## Codebase Navigation Tools
+
+Two MCP tools are registered in `.mcp.json` at the repo root and appear
+automatically in any MCP-capable agent session.
+
+**`openmc_rag_search`** — Semantic search across the codebase (C++, Python, RST
+docs). Finds code by meaning, not just text match. Surfaces related code across
+subsystems even when naming differs (e.g., "particle RNG seeding" finds code
+across transport, restart, and random ray modes — files you would never find
+with `grep "particle seed"`). The index uses a small 22M-param embedding model
+(384-dim). Phrase-level natural-language queries work much better than single
+keywords or symbol names.
+
+**`openmc_rag_rebuild`** — Rebuild the RAG vector index. Call after pulling new
+code or switching branches. The first RAG search of each session will report
+the index status and ask whether to rebuild — you can also call this explicitly.
+
+### Why RAG matters
+
+OpenMC is large enough that changes in one subsystem can silently break
+invariants that distant subsystems depend on — and those distant files often
+use different naming, so grep won't find them. The RAG search finds code by
+meaning, surfacing files you wouldn't have thought to open.
+
+An agent reviewed a large OpenMC PR without RAG. It found 1 of 11 serious
+bugs. Its post-mortem:
+
+> **I treated the diff as a closed system.** I verified internal consistency of
+> the changed code obsessively, but never built a global understanding of how
+> the changed code fits into the wider codebase. The diff altered assumptions
+> that code elsewhere silently relied on — but I couldn't see that because I
+> never looked beyond the diff. I couldn't see the forest for the trees.
+>
+> **Why I resisted RAG:** Overconfidence. My internal model was "I can see the
+> diff, I understand the data structures, I can trace the logic." The diff felt
+> self-contained. RAG felt like it would return noisy results about tangentially
+> related code. But in a codebase this large, changes in one subsystem can
+> quietly break invariants that distant subsystems depend on — and you need
+> global awareness to foresee that.
+>
+> **In the post-mortem**, I re-ran the RAG queries I should have run during the
+> review. They directly surfaced the files containing the bugs I missed — files
+> I never thought to open because they weren't in the diff.
+
+The takeaway: when reviewing or modifying code, ask yourself "what else in this
+codebase might depend on the behavior I'm changing?" If you aren't sure, that's
+a good time for a RAG query. It won't replace the grep-based investigation you
+should already be doing — but it can surface files you wouldn't have thought to
+open.
+
+### Workflow for contributors
+
+1. Create a feature/bugfix branch off `develop`
+2. Make changes and commit to the feature branch
+3. Open a pull request to merge the feature branch into `develop`
+4. A committer reviews and merges the PR into `develop`
+
+## Critical Build & Test Workflows
+
+### Build Dependencies
+- **C++17 compiler**: GCC, Clang, or Intel
+- **CMake** (3.16+): Required for configuring and building the C++ library
+- **HDF5**: Required for cross section data and output file formats
+- **libpng**: Used for generating visualization when OpenMC is run in plotting mode
+
+Without CMake and HDF5, OpenMC cannot be compiled.
+
+### Building the C++ Library
+```bash
+# Configure with CMake (from build/ directory)
+cmake .. -DOPENMC_USE_MPI=ON -DOPENMC_USE_OPENMP=ON -DCMAKE_BUILD_TYPE=RelWithDebInfo
+
+# Available CMake options (all default OFF except OPENMC_USE_OPENMP and OPENMC_BUILD_TESTS):
+# -DOPENMC_USE_OPENMP=ON/OFF        # OpenMP parallelism
+# -DOPENMC_USE_MPI=ON/OFF           # MPI support
+# -DOPENMC_USE_DAGMC=ON/OFF         # CAD geometry support
+# -DOPENMC_USE_LIBMESH=ON/OFF       # Unstructured mesh
+# -DOPENMC_ENABLE_PROFILE=ON/OFF    # Profiling flags
+# -DOPENMC_ENABLE_COVERAGE=ON/OFF   # Coverage analysis
+
+# Build
+make -j
+
+# C++ unit tests (uses Catch2)
+ctest
+```
+
+### Python Development
+```bash
+# Install in development mode (requires building C++ library first)
+pip install -e .
+
+# Python tests (uses pytest)
+pytest tests/unit_tests/           # Fast unit tests
+pytest tests/regression_tests/     # Full regression suite (requires nuclear data)
+```
+
+### Nuclear Data Setup (CRITICAL for Running OpenMC)
+Most tests require the NNDC HDF5 nuclear cross-section library.
+
+**Important**: Check if `OPENMC_CROSS_SECTIONS` is already set in the user's
+environment before downloading, as many users already have nuclear data
+installed. Though do note that if this variable is present that it may point to
+different cross section data and that the NNDC data is required for tests to
+pass.
+
+**If not already configured, download and setup:**
+```bash
+# Download NNDC HDF5 cross section library (~800 MB compressed)
+wget -q -O - https://anl.box.com/shared/static/teaup95cqv8s9nn56hfn7ku8mmelr95p.xz | tar -C $HOME -xJ
+
+# Set environment variable (add to ~/.bashrc or ~/.zshrc for persistence)
+export OPENMC_CROSS_SECTIONS=$HOME/nndc_hdf5/cross_sections.xml
+```
+
+**Alternative**: Use the provided download script (checks if data exists before downloading):
+```bash
+bash tools/ci/download-xs.sh  # Downloads both NNDC HDF5 and ENDF/B-VII.1 data
+```
+
+Without this data, regression tests will fail with "No cross_sections.xml file
+found" errors, or, in the case that alternative cross section data is configured
+the tests will execute but will not pass. The `cross_sections.xml` file is an
+index listing paths to individual HDF5 nuclear data files for each nuclide.
+
+## Testing Expectations
+
+### Environment Requirements
+
+  - **Data**: As described above, OpenMC's test suite requires OpenMC to be configured with NNDC data.
+  - **OpenMP Settings**: OpenMC's tests may fail is more than two OpenMP threads are used. The environment variable `OMP_NUM_THREADS=2` should be set to avoid sporadic test failures.
+  - **Executable configuration**: The OpenMC executable should compiled with debug symbols enabled.
+
+### C++ Tests
+Located in `tests/cpp_unit_tests/`, use Catch2 framework. Run via `ctest` after building with `-DOPENMC_BUILD_TESTS=ON`.
+
+### Python Unit Tests
+Located in `tests/unit_tests/`, these are fast, standalone tests that verify Python API functionality without running full simulations. Use standard pytest patterns:
+
+**Categories**:
+- **API validation**: Test object creation, property setters/getters, XML serialization (e.g., `test_material.py`, `test_cell.py`, `test_source.py`)
+- **Data processing**: Test nuclear data handling, cross sections, depletion chains (e.g., `test_data_neutron.py`, `test_deplete_chain.py`)
+- **Library bindings**: Test `openmc.lib` ctypes interface with `model.init_lib()`/`model.finalize_lib()` (e.g., `test_lib.py`)
+- **Geometry operations**: Test bounding boxes, containment, lattice generation (e.g., `test_bounding_box.py`, `test_lattice.py`)
+
+**Common patterns**:
+- Use fixtures from `tests/unit_tests/conftest.py` (e.g., `uo2`, `water`, `sphere_model`)
+- Test invalid inputs with `pytest.raises(ValueError)` or `pytest.raises(TypeError)`
+- Use `run_in_tmpdir` fixture for tests that create files
+- Tests with `openmc.lib` require calling `model.init_lib()` in try/finally with `model.finalize_lib()`
+
+**Example**:
+```python
+def test_material_properties():
+    m = openmc.Material()
+    m.add_nuclide('U235', 1.0)
+    assert 'U235' in m.nuclides
+
+    with pytest.raises(TypeError):
+        m.add_nuclide('H1', '1.0')  # Invalid type
+```
+
+Unit tests should be fast. For tests requiring simulation output, use regression tests instead.
+
+### Python Regression Tests
+Regression tests compare OpenMC output against reference data. **Prefer using existing models from `openmc.examples` or those found in tests/unit_tests/conftest.py** (like `pwr_pin_cell()`, `pwr_assembly()`, `slab_mg()`) rather than building from scratch.
+
+**Test Harness Types** (in `tests/testing_harness.py`):
+- **PyAPITestHarness**: Standard harness for Python API tests. Compares `inputs_true.dat` (XML hash) and `results_true.dat` (statepoint k-eff and tally values). Requires `model.xml` generation.
+- **HashedPyAPITestHarness**: Like PyAPITestHarness but hashes the results for compact comparison
+- **TolerantPyAPITestHarness**: For tests with floating-point non-associativity (e.g., random ray solver with single precision). Uses relative tolerance comparisons.
+- **WeightWindowPyAPITestHarness**: Compares weight window bounds from `weight_windows.h5`
+- **CollisionTrackTestHarness**: Compares collision track data from `collision_track.h5` against `collision_track_true.h5`
+- **TestHarness**: Base harness for XML-based tests (no Python model building)
+- **PlotTestHarness**: Compares plot output files (PNG or voxel HDF5)
+- **CMFDTestHarness**: Specialized for CMFD acceleration tests
+- **ParticleRestartTestHarness**: Tests particle restart functionality
+
+Almost all cases use either `PyAPITestHarness` or `HashedPyAPITestHarness`
+
+**Example Test**:
+```python
+from openmc.examples import pwr_pin_cell
+from tests.testing_harness import PyAPITestHarness
+
+def test_my_feature():
+    model = pwr_pin_cell()
+    model.settings.particles = 1000  # Modify to exercise feature
+    harness = PyAPITestHarness('statepoint.10.h5', model)
+    harness.main()
+```
+
+**Workflow**: Create `test.py` and `__init__.py` in `tests/regression_tests/my_test/`, run `pytest --update` to generate reference files (`inputs_true.dat`, `results_true.dat`, etc.), then verify with `pytest` without `--update`. Test results should be generated with `-DOPENMC_ENABLE_STRICT_FP=on` to ensure reproducibility across platforms and optimization levels.
+
+**Critical**: When modifying OpenMC code, regenerate affected test references with `pytest --update` and commit updated reference files.
+
+### Test Configuration
+
+`pytest.ini` sets: `python_files = test*.py`, `python_classes = NoThanks` (disables class-based test collection).
+
+### Testing Options
+
+For builds of OpenMC with MPI enabled, the `--mpi` flag should be passed to the test suite to ensure that appropriate tests are executed using two MPI processes.
+
+The entire test suite can be executed with OpenMC running in event-based mode (instead of the default history-based mode) by providing the `--event` flag to the `pytest` command.
+
+## Cross-Language Boundaries
+
+The C API (defined in `include/openmc/capi.h`) exposes C++ functionality to Python via ctypes bindings in `openmc/lib/`. Example:
+```cpp
+// C++ API in capi.h
+extern "C" int openmc_run();
+
+// Python binding in openmc/lib/core.py
+_dll.openmc_run.restype = c_int
+def run():
+    _dll.openmc_run()
+```
+
+When modifying C++ public APIs, update corresponding ctypes signatures in `openmc/lib/*.py`.
+
+## Code Style & Conventions
+
+### C++ Style (enforced by .clang-format)
+  OpenMC generally tries to follow C++ core guidelines where possible
+  (https://isocpp.github.io/CppCoreGuidelines/CppCoreGuidelines) and follow
+  modern C++ practices (e.g. RAII) whenever possible.
+
+- **Naming**:
+  - Classes: `CamelCase` (e.g., `HexLattice`)
+  - Functions/methods: `snake_case` (e.g., `get_indices`)
+  - Variables: `snake_case` with trailing underscore for class members (e.g., `n_particles_`, `energy_`)
+  - Constants: `UPPER_SNAKE_CASE` (e.g., `SQRT_PI`)
+- **Namespaces**: All code in `openmc::` namespace, global state in sub-namespaces
+- **Include order**: Related header first, then C/C++ stdlib, third-party libs, local headers
+- **Comments**: C++-style (`//`) only, never C-style (`/* */`)
+- **Standard**: C++17 features allowed
+- **Formatting**: Run `clang-format` (version 18) before committing; install via `tools/dev/install-commit-hooks.sh`
+
+### Python Style
+- **PEP8** compliant
+- **Docstrings**: numpydoc format for all public functions/methods
+- **Type hints**: Use sparingly, primarily for complex signatures
+- **Path handling**: Use `pathlib.Path` for filesystem operations, accept `str | os.PathLike` in function arguments
+- **Dependencies**: Core dependencies only (numpy, scipy, h5py, pandas, matplotlib, lxml, ipython, uncertainties, endf). Other packages must be optional
+- **Python version**: Minimum 3.11 (as of Nov 2025)
+
+### ID Management Pattern (Python)
+When creating geometry objects, IDs can be auto-assigned or explicit:
+```python
+# Auto-assigned ID
+cell = openmc.Cell()  # Gets next available ID
+
+# Explicit ID
+cell = openmc.Cell(id=10)  # Warning if ID already used
+
+# Reset all IDs (useful in test fixtures)
+openmc.reset_auto_ids()
+```
+
+### Input Validation Pattern (Python)
+All setters use checkvalue functions:
+```python
+import openmc.checkvalue as cv
+
+@property
+def temperature(self):
+    return self._temperature
+
+@temperature.setter
+def temperature(self, temp):
+    cv.check_type('temperature', temp, Real)
+    cv.check_greater_than('temperature', temp, 0.0)
+    self._temperature = temp
+```
+
+### Working with HDF5 Files
+C++ uses custom HDF5 wrappers in `src/hdf5_interface.cpp`. Python uses h5py directly. Statepoint format version is `VERSION_STATEPOINT` in `include/openmc/constants.h`.
+
+### Conditional Compilation
+Check for optional features:
+```cpp
+#ifdef OPENMC_MPI
+  // MPI-specific code
+#endif
+
+#ifdef OPENMC_DAGMC
+  // DAGMC-specific code
+#endif
+```
+
+## Documentation
+
+- **User docs**: Sphinx documentation in `docs/source/` hosted at https://docs.openmc.org
+- **C++ docs**: Doxygen-style comments with `\brief`, `\param` tags
+- **Python docs**: numpydoc format docstrings
+
+## Common Pitfalls
+
+1. **Forgetting nuclear data**: Tests fail without `OPENMC_CROSS_SECTIONS` environment variable
+2. **ID conflicts**: Python objects with duplicate IDs trigger `IDWarning`, use `reset_auto_ids()` between tests
+3. **MPI builds**: Code must work with and without MPI; use `#ifdef OPENMC_MPI` guards
+4. **Path handling**: Use `pathlib.Path` in new Python code, not `os.path`
+5. **Clang-format version**: CI uses version 18; other versions may produce different formatting

CITATION.cff

@@ -0,0 +1,68 @@
+cff-version: 1.2.0
+message: "If you use this software, please cite it as below."
+title: OpenMC
+authors:
+- family-names: Romano
+  given-names: Paul K.
+  orcid: "https://orcid.org/0000-0002-1147-045X"
+- family-names: Shriwise
+  given-names: Patrick C.
+  orcid: "https://orcid.org/0000-0002-3979-7665"
+- family-names: Shimwell
+  given-names: Jonathan
+  orcid: "https://orcid.org/0000-0001-6909-0946"
+- family-names: Harper
+  given-names: Sterling
+- family-names: Boyd
+  given-names: Will
+- family-names: Nelson
+  given-names: Adam G.
+  orcid: "https://orcid.org/0000-0002-3614-0676"
+- family-names: Tramm
+  given-names: John R.
+  orcid: "https://orcid.org/0000-0002-5397-4402"
+- family-names: Ridley
+  given-names: Gavin
+  orcid: "https://orcid.org/0000-0003-1635-8042"
+- family-names: Johnson
+  given-names: Andrew
+  orcid: "https://orcid.org/0000-0003-2125-8775"
+- family-names: Peterson
+  given-names: Ethan E.
+  orcid: "https://orcid.org/0000-0002-5694-7194"
+- family-names: Herman
+  given-names: Bryan R.
+preferred-citation:
+  authors:
+  - family-names: Romano
+    given-names: Paul K.
+    orcid: "https://orcid.org/0000-0002-1147-045X"
+  - family-names: Horelik
+    given-names: Nicholas E.
+  - family-names: Herman
+    given-names: Bryan R.
+  - family-names: Nelson
+    given-names: Adam G.
+    orcid: "https://orcid.org/0000-0002-3614-0676"
+  - family-names: Forget
+    given-names: Benoit
+    orcid: "https://orcid.org/0000-0003-1459-7672"
+  - family-names: Smith
+    given-names: Kord
+  contact:
+  - family-names: Romano
+    given-names: Paul K.
+    orcid: "https://orcid.org/0000-0002-1147-045X"
+  doi: 10.1016/j.anucene.2014.07.048
+  issn: 0306-4549
+  volume: 82
+  journal: Annals of Nuclear Energy
+  publisher:
+    name: Elsevier
+  start: 90
+  end: 97
+  year: 2015
+  month: 8
+  title: "OpenMC: A state-of-the-art Monte Carlo code for research and development"
+  type: article
+  url: "https://doi.org/10.1016/j.anucene.2014.07.048"

CLAUDE.md

@@ -0,0 +1,14 @@
+## OpenMC Codebase Tools
+
+Read the FULL `AGENTS.md` in this directory before starting work. It contains
+project context, coding conventions, and documentation of the RAG search tools
+registered in `.mcp.json`.
+
+### Claude Code-specific: first-call behavior
+
+The first `openmc_rag_search` call of each session returns an index status
+message instead of search results. When this happens, you MUST use the
+`AskUserQuestion` tool to present the rebuild/use-existing choice to the user.
+Do not ask conversationally — always use the widget. Do not skip this step even
+if the index looks current — the user may have uncommitted changes that warrant
+a rebuild.

CMakeLists.txt

@@ -1,11 +1,18 @@
-cmake_minimum_required(VERSION 3.10 FATAL_ERROR)
+cmake_minimum_required(VERSION 3.16 FATAL_ERROR)
 project(openmc C CXX)
 
-# Set version numbers
-set(OPENMC_VERSION_MAJOR 0)
-set(OPENMC_VERSION_MINOR 14)
-set(OPENMC_VERSION_RELEASE 0)
-set(OPENMC_VERSION ${OPENMC_VERSION_MAJOR}.${OPENMC_VERSION_MINOR}.${OPENMC_VERSION_RELEASE})
+# Set module path
+set(CMAKE_MODULE_PATH ${CMAKE_CURRENT_SOURCE_DIR}/cmake/Modules)
+
+include(GetVersionFromGit)
+
+# Output version information
+message(STATUS "OpenMC version: ${OPENMC_VERSION}")
+message(STATUS "OpenMC dev state: ${OPENMC_DEV_STATE}")
+message(STATUS "OpenMC commit hash: ${OPENMC_COMMIT_HASH}")
+message(STATUS "OpenMC commit count: ${OPENMC_COMMIT_COUNT}")
+
+# Generate version.h
 configure_file(include/openmc/version.h.in "${CMAKE_BINARY_DIR}/include/openmc/version.h" @ONLY)
 
 # Setup output directories
@@ -13,12 +20,9 @@ set(CMAKE_ARCHIVE_OUTPUT_DIRECTORY ${CMAKE_BINARY_DIR}/lib)
 set(CMAKE_LIBRARY_OUTPUT_DIRECTORY ${CMAKE_BINARY_DIR}/lib)
 set(CMAKE_RUNTIME_OUTPUT_DIRECTORY ${CMAKE_BINARY_DIR}/bin)
 
-# Set module path
-set(CMAKE_MODULE_PATH ${CMAKE_CURRENT_SOURCE_DIR}/cmake/Modules)
-
-# Allow user to specify <project>_ROOT variables
-if (CMAKE_VERSION VERSION_GREATER_EQUAL 3.12)
-  cmake_policy(SET CMP0074 NEW)
+# Generate compile_commands.json for clangd and other tools
+if("${CMAKE_EXPORT_COMPILE_COMMANDS}" STREQUAL "")
+  set(CMAKE_EXPORT_COMPILE_COMMANDS ON)
 endif()
 
 # Enable correct usage of CXX_EXTENSIONS
@@ -37,8 +41,20 @@ option(OPENMC_ENABLE_COVERAGE "Compile with coverage analysis flags"
 option(OPENMC_USE_DAGMC       "Enable support for DAGMC (CAD) geometry"              OFF)
 option(OPENMC_USE_LIBMESH     "Enable support for libMesh unstructured mesh tallies" OFF)
 option(OPENMC_USE_MPI         "Enable MPI"                                           OFF)
-option(OPENMC_USE_MCPL        "Enable MCPL"                                          OFF)
-option(OPENMC_USE_NCRYSTAL    "Enable support for NCrystal scattering"               OFF)
+option(OPENMC_USE_UWUW        "Enable UWUW"                                          OFF)
+option(OPENMC_FORCE_VENDORED_LIBS "Explicitly use submodules defined in 'vendor'"    OFF)
+option(OPENMC_ENABLE_STRICT_FP "Enable strict FP flags to improve test portability"  OFF)
+
+message(STATUS "OPENMC_USE_OPENMP ${OPENMC_USE_OPENMP}")
+message(STATUS "OPENMC_BUILD_TESTS ${OPENMC_BUILD_TESTS}")
+message(STATUS "OPENMC_ENABLE_PROFILE ${OPENMC_ENABLE_PROFILE}")
+message(STATUS "OPENMC_ENABLE_COVERAGE ${OPENMC_ENABLE_COVERAGE}")
+message(STATUS "OPENMC_USE_DAGMC ${OPENMC_USE_DAGMC}")
+message(STATUS "OPENMC_USE_LIBMESH ${OPENMC_USE_LIBMESH}")
+message(STATUS "OPENMC_USE_MPI ${OPENMC_USE_MPI}")
+message(STATUS "OPENMC_USE_UWUW ${OPENMC_USE_UWUW}")
+message(STATUS "OPENMC_FORCE_VENDORED_LIBS ${OPENMC_FORCE_VENDORED_LIBS}")
+message(STATUS "OPENMC_ENABLE_STRICT_FP ${OPENMC_ENABLE_STRICT_FP}")
 
 # Warnings for deprecated options
 foreach(OLD_OPT IN ITEMS "openmp" "profile" "coverage" "dagmc" "libmesh")
@@ -80,6 +96,27 @@ if(NOT CMAKE_BUILD_TYPE)
     set(CMAKE_BUILD_TYPE RelWithDebInfo CACHE STRING "Choose the type of build" FORCE)
 endif()
 
+#===============================================================================
+# When STRICT_FP is enabled, remove NDEBUG from RelWithDebInfo flags so that
+# assert() remains active. CMake normally adds -DNDEBUG for both Release and
+# RelWithDebInfo, which disables C/C++ assert() statements.
+#===============================================================================
+
+if(OPENMC_ENABLE_STRICT_FP)
+  foreach(FLAG_VAR CMAKE_CXX_FLAGS_RELWITHDEBINFO CMAKE_C_FLAGS_RELWITHDEBINFO)
+    string(REPLACE "-DNDEBUG" "" ${FLAG_VAR} "${${FLAG_VAR}}")
+    string(REPLACE "/DNDEBUG" "" ${FLAG_VAR} "${${FLAG_VAR}}")
+  endforeach()
+endif()
+
+#===============================================================================
+# OpenMP for shared-memory parallelism (and GPU support some day!)
+#===============================================================================
+
+if(OPENMC_USE_OPENMP)
+  find_package(OpenMP REQUIRED)
+endif()
+
 #===============================================================================
 # MPI for distributed-memory parallelism
 #===============================================================================
@@ -101,15 +138,6 @@ macro(find_package_write_status pkg)
   endif()
 endmacro()
 
-#===============================================================================
-# NCrystal Scattering Support
-#===============================================================================
-
-if(OPENMC_USE_NCRYSTAL)
-  find_package(NCrystal REQUIRED)
-  message(STATUS "Found NCrystal: ${NCrystal_DIR} (version ${NCrystal_VERSION})")
-endif()
-
 #===============================================================================
 # DAGMC Geometry Support - need DAGMC/MOAB
 #===============================================================================
@@ -117,10 +145,15 @@ endif()
 if(OPENMC_USE_DAGMC)
   find_package(DAGMC REQUIRED PATH_SUFFIXES lib/cmake)
   if (${DAGMC_VERSION} VERSION_LESS 3.2.0)
-    message(FATAL_ERROR "Discovered DAGMC Version: ${DAGMC_VERSION}. \
-    Please update DAGMC to version 3.2.0 or greater.")
+    message(FATAL_ERROR "Discovered DAGMC Version: ${DAGMC_VERSION}."
+      "Please update DAGMC to version 3.2.0 or greater.")
   endif()
   message(STATUS "Found DAGMC: ${DAGMC_DIR} (version ${DAGMC_VERSION})")
+
+  # Check if UWUW is needed and available
+  if(OPENMC_USE_UWUW AND NOT DAGMC_BUILD_UWUW)
+    message(FATAL_ERROR "UWUW is enabled but DAGMC was not configured with UWUW.")
+  endif()
 endif()
 
 #===============================================================================
@@ -155,6 +188,11 @@ if(NOT DEFINED HDF5_PREFER_PARALLEL)
 endif()
 
 find_package(HDF5 REQUIRED COMPONENTS C HL)
+
+# Remove HDF5 transitive dependencies that are system libraries
+list(FILTER HDF5_LIBRARIES EXCLUDE REGEX ".*lib(pthread|dl|m).*")
+message(STATUS "HDF5 Libraries: ${HDF5_LIBRARIES}")
+
 if(HDF5_IS_PARALLEL)
   if(NOT OPENMC_USE_MPI)
     message(FATAL_ERROR "Parallel HDF5 was detected, but MPI was not enabled.\
@@ -171,31 +209,33 @@ if(${HDF5_VERSION} VERSION_GREATER_EQUAL 1.12.0)
   list(APPEND cxxflags -DH5Oget_info_by_idx_vers=1 -DH5O_info_t_vers=1)
 endif()
 
-#===============================================================================
-# MCPL
-#===============================================================================
-
-if (OPENMC_USE_MCPL)
-  find_package(MCPL REQUIRED)
-  message(STATUS "Found MCPL: ${MCPL_DIR} (found version \"${MCPL_VERSION}\")")
-endif()
-
 #===============================================================================
 # Set compile/link flags based on which compiler is being used
 #===============================================================================
 
-# Skip for Visual Studio which has its own configurations through GUI
-if(NOT MSVC)
-
-if(OPENMC_USE_OPENMP)
-  find_package(OpenMP)
-  if(OPENMP_FOUND)
-    # In CMake 3.9+, can use the OpenMP::OpenMP_CXX imported target
-    list(APPEND cxxflags ${OpenMP_CXX_FLAGS})
-    list(APPEND ldflags ${OpenMP_CXX_FLAGS})
+# When OPENMC_ENABLE_STRICT_FP is enabled, disable compiler optimizations that change
+# floating-point results relative to -O0, improving cross-platform and
+# cross-optimization-level reproducibility for regression testing:
+#   -ffp-contract=off    Prevents FMA contraction (fused multiply-add changes rounding)
+#   -fno-builtin         Prevents replacing math function calls (pow, exp, log, etc.)
+#                        with builtin versions that may differ from libm
+# By default (OFF), the compiler is free to use all optimizations for best
+# performance.
+if(OPENMC_ENABLE_STRICT_FP)
+  include(CheckCXXCompilerFlag)
+  check_cxx_compiler_flag(-ffp-contract=off SUPPORTS_FP_CONTRACT_OFF)
+  if(SUPPORTS_FP_CONTRACT_OFF)
+    list(APPEND cxxflags -ffp-contract=off)
+  endif()
+  check_cxx_compiler_flag(-fno-builtin SUPPORTS_NO_BUILTIN)
+  if(SUPPORTS_NO_BUILTIN)
+    list(APPEND cxxflags -fno-builtin)
   endif()
 endif()
 
+# Skip for Visual Studio which has its own configurations through GUI
+if(NOT MSVC)
+
 set(CMAKE_POSITION_INDEPENDENT_CODE ON)
 
 if(OPENMC_ENABLE_PROFILE)
@@ -216,8 +256,6 @@ endif()
 #===============================================================================
 # Update git submodules as needed
 #===============================================================================
-
-find_package(Git)
 if(GIT_FOUND AND EXISTS "${CMAKE_CURRENT_SOURCE_DIR}/.git")
   option(GIT_SUBMODULE "Check submodules during build" ON)
   if(GIT_SUBMODULE)
@@ -242,49 +280,30 @@ endif()
 # pugixml library
 #===============================================================================
 
-find_package_write_status(pugixml)
-if (NOT pugixml_FOUND)
+if(OPENMC_FORCE_VENDORED_LIBS)
   add_subdirectory(vendor/pugixml)
   set_target_properties(pugixml PROPERTIES CXX_STANDARD 14 CXX_EXTENSIONS OFF)
+else()
+  find_package_write_status(pugixml)
+  if (NOT pugixml_FOUND)
+    add_subdirectory(vendor/pugixml)
+    set_target_properties(pugixml PROPERTIES CXX_STANDARD 14 CXX_EXTENSIONS OFF)
+  endif()
 endif()
 
 #===============================================================================
 # {fmt} library
 #===============================================================================
 
-find_package_write_status(fmt)
-if (NOT fmt_FOUND)
+if(OPENMC_FORCE_VENDORED_LIBS)
   set(FMT_INSTALL ON CACHE BOOL "Generate the install target.")
   add_subdirectory(vendor/fmt)
-endif()
-
-#===============================================================================
-# xtensor header-only library
-#===============================================================================
-
-# CMake 3.13+ will complain about policy CMP0079 unless it is set explicitly
-if (CMAKE_VERSION VERSION_GREATER_EQUAL 3.13)
-  cmake_policy(SET CMP0079 NEW)
-endif()
-
-find_package_write_status(xtensor)
-if (NOT xtensor_FOUND)
-  add_subdirectory(vendor/xtl)
-  set(xtl_DIR ${CMAKE_CURRENT_BINARY_DIR}/vendor/xtl)
-  add_subdirectory(vendor/xtensor)
-endif()
-
-#===============================================================================
-# GSL header-only library
-#===============================================================================
-
-find_package_write_status(gsl-lite)
-if (NOT gsl-lite_FOUND)
-  add_subdirectory(vendor/gsl-lite)
-
-  # Make sure contract violations throw exceptions
-  target_compile_definitions(gsl-lite-v1 INTERFACE GSL_THROW_ON_CONTRACT_VIOLATION)
-  target_compile_definitions(gsl-lite-v1 INTERFACE gsl_CONFIG_ALLOWS_NONSTRICT_SPAN_COMPARISON=1)
+else()
+  find_package_write_status(fmt)
+  if (NOT fmt_FOUND)
+    set(FMT_INSTALL ON CACHE BOOL "Generate the install target.")
+    add_subdirectory(vendor/fmt)
+  endif()
 endif()
 
 #===============================================================================
@@ -292,9 +311,13 @@ endif()
 #===============================================================================
 
 if(OPENMC_BUILD_TESTS)
-  find_package_write_status(Catch2)
-  if (NOT Catch2_FOUND)
+  if (OPENMC_FORCE_VENDORED_LIBS)
     add_subdirectory(vendor/Catch2)
+  else()
+    find_package_write_status(Catch2)
+    if (NOT Catch2_FOUND)
+      add_subdirectory(vendor/Catch2)
+    endif()
   endif()
 endif()
 
@@ -332,11 +355,14 @@ endif()
 #===============================================================================
 
 list(APPEND libopenmc_SOURCES
+  src/atomic_mass.cpp
   src/bank.cpp
   src/boundary_condition.cpp
   src/bremsstrahlung.cpp
   src/cell.cpp
+  src/chain.cpp
   src/cmfd_solver.cpp
+  src/collision_track.cpp
   src/cross_sections.cpp
   src/dagmc.cpp
   src/distribution.cpp
@@ -353,6 +379,7 @@ list(APPEND libopenmc_SOURCES
   src/geometry.cpp
   src/geometry_aux.cpp
   src/hdf5_interface.cpp
+  src/ifp.cpp
   src/initialize.cpp
   src/lattice.cpp
   src/material.cpp
@@ -363,11 +390,13 @@ list(APPEND libopenmc_SOURCES
   src/mgxs.cpp
   src/mgxs_interface.cpp
   src/ncrystal_interface.cpp
+  src/ncrystal_load.cpp
   src/nuclide.cpp
   src/output.cpp
   src/particle.cpp
   src/particle_data.cpp
   src/particle_restart.cpp
+  src/particle_type.cpp
   src/photon.cpp
   src/physics.cpp
   src/physics_common.cpp
@@ -377,6 +406,13 @@ list(APPEND libopenmc_SOURCES
   src/progress_bar.cpp
   src/random_dist.cpp
   src/random_lcg.cpp
+  src/random_ray/random_ray_simulation.cpp
+  src/random_ray/random_ray.cpp
+  src/random_ray/flat_source_domain.cpp
+  src/random_ray/linear_source_domain.cpp
+  src/random_ray/moment_matrix.cpp
+  src/random_ray/source_region.cpp
+  src/ray.cpp
   src/reaction.cpp
   src/reaction_product.cpp
   src/scattdata.cpp
@@ -406,16 +442,24 @@ list(APPEND libopenmc_SOURCES
   src/tallies/filter_energyfunc.cpp
   src/tallies/filter_legendre.cpp
   src/tallies/filter_material.cpp
+  src/tallies/filter_materialfrom.cpp
   src/tallies/filter_mesh.cpp
+  src/tallies/filter_meshborn.cpp
+  src/tallies/filter_meshmaterial.cpp
   src/tallies/filter_meshsurface.cpp
   src/tallies/filter_mu.cpp
+  src/tallies/filter_musurface.cpp
+  src/tallies/filter_parent_nuclide.cpp
   src/tallies/filter_particle.cpp
+  src/tallies/filter_particle_production.cpp
   src/tallies/filter_polar.cpp
+  src/tallies/filter_reaction.cpp
   src/tallies/filter_sph_harm.cpp
   src/tallies/filter_sptl_legendre.cpp
   src/tallies/filter_surface.cpp
   src/tallies/filter_time.cpp
   src/tallies/filter_universe.cpp
+  src/tallies/filter_weight.cpp
   src/tallies/filter_zernike.cpp
   src/tallies/tally.cpp
   src/tallies/tally_scoring.cpp
@@ -476,22 +520,10 @@ if (OPENMC_USE_MPI)
   target_compile_definitions(libopenmc PUBLIC -DOPENMC_MPI)
 endif()
 
-# Set git SHA1 hash as a compile definition
-if(GIT_FOUND)
-  execute_process(COMMAND ${GIT_EXECUTABLE} rev-parse HEAD
-                  WORKING_DIRECTORY ${CMAKE_CURRENT_SOURCE_DIR}
-                  RESULT_VARIABLE GIT_SHA1_SUCCESS
-                  OUTPUT_VARIABLE GIT_SHA1
-                  ERROR_QUIET OUTPUT_STRIP_TRAILING_WHITESPACE)
-  if(GIT_SHA1_SUCCESS EQUAL 0)
-    target_compile_definitions(libopenmc PRIVATE -DGIT_SHA1="${GIT_SHA1}")
-  endif()
-endif()
-
 # target_link_libraries treats any arguments starting with - but not -l as
 # linker flags. Thus, we can pass both linker flags and libraries together.
 target_link_libraries(libopenmc ${ldflags} ${HDF5_LIBRARIES} ${HDF5_HL_LIBRARIES}
-                      xtensor gsl::gsl-lite-v1 fmt::fmt)
+                      fmt::fmt ${CMAKE_DL_LIBS})
 
 if(TARGET pugixml::pugixml)
   target_link_libraries(libopenmc pugixml::pugixml)
@@ -500,12 +532,20 @@ else()
 endif()
 
 if(OPENMC_USE_DAGMC)
-  target_compile_definitions(libopenmc PRIVATE DAGMC)
-  target_link_libraries(libopenmc dagmc-shared uwuw-shared)
+  target_compile_definitions(libopenmc PUBLIC OPENMC_DAGMC_ENABLED)
+  target_link_libraries(libopenmc dagmc-shared)
+
+  if(OPENMC_USE_UWUW)
+    target_compile_definitions(libopenmc PRIVATE OPENMC_UWUW_ENABLED)
+    target_link_libraries(libopenmc uwuw-shared)
+  endif()
+elseif(OPENMC_USE_UWUW)
+  set(OPENMC_USE_UWUW OFF)
+  message(FATAL_ERROR "DAGMC must be enabled when UWUW is enabled.")
 endif()
 
 if(OPENMC_USE_LIBMESH)
-  target_compile_definitions(libopenmc PRIVATE LIBMESH)
+  target_compile_definitions(libopenmc PRIVATE OPENMC_LIBMESH_ENABLED)
   target_link_libraries(libopenmc PkgConfig::LIBMESH)
 endif()
 
@@ -514,6 +554,10 @@ if (PNG_FOUND)
   target_link_libraries(libopenmc PNG::PNG)
 endif()
 
+if (OPENMC_USE_OPENMP)
+  target_link_libraries(libopenmc OpenMP::OpenMP_CXX)
+endif()
+
 if (OPENMC_USE_MPI)
   target_link_libraries(libopenmc MPI::MPI_CXX)
 endif()
@@ -524,16 +568,6 @@ if (OPENMC_BUILD_TESTS)
   add_subdirectory(tests/cpp_unit_tests)
 endif()
 
-if (OPENMC_USE_MCPL)
-  target_compile_definitions(libopenmc PUBLIC OPENMC_MCPL)
-  target_link_libraries(libopenmc MCPL::mcpl)
-endif()
-
-if(OPENMC_USE_NCRYSTAL)
-  target_compile_definitions(libopenmc PRIVATE NCRYSTAL)
-  target_link_libraries(libopenmc NCrystal::NCrystal)
-endif()
-
 #===============================================================================
 # Log build info that this executable can report later
 #===============================================================================
@@ -546,6 +580,9 @@ endif()
 if (OPENMC_ENABLE_COVERAGE)
   target_compile_definitions(libopenmc PRIVATE COVERAGEBUILD)
 endif()
+if (OPENMC_ENABLE_STRICT_FP)
+  target_compile_definitions(libopenmc PRIVATE OPENMC_ENABLE_STRICT_FP)
+endif()
 
 #===============================================================================
 # openmc executable
@@ -556,9 +593,9 @@ target_compile_options(openmc PRIVATE ${cxxflags})
 target_include_directories(openmc PRIVATE ${CMAKE_BINARY_DIR}/include)
 target_link_libraries(openmc libopenmc)
 
-# Ensure C++14 standard is used and turn off GNU extensions
-target_compile_features(openmc PUBLIC cxx_std_14)
-target_compile_features(libopenmc PUBLIC cxx_std_14)
+# Ensure C++17 standard is used and turn off GNU extensions
+target_compile_features(openmc PUBLIC cxx_std_17)
+target_compile_features(libopenmc PUBLIC cxx_std_17)
 set_target_properties(openmc libopenmc PROPERTIES CXX_EXTENSIONS OFF)
 
 #===============================================================================
@@ -574,9 +611,7 @@ add_custom_command(TARGET libopenmc POST_BUILD
 #===============================================================================
 # Install executable, scripts, manpage, license
 #===============================================================================
-
-configure_file(cmake/OpenMCConfig.cmake.in "${CMAKE_BINARY_DIR}${CMAKE_FILES_DIRECTORY}/OpenMCConfig.cmake" @ONLY)
-configure_file(cmake/OpenMCConfigVersion.cmake.in "${CMAKE_BINARY_DIR}${CMAKE_FILES_DIRECTORY}/OpenMCConfigVersion.cmake" @ONLY)
+include(CMakePackageConfigHelpers)
 
 set(INSTALL_CONFIGDIR ${CMAKE_INSTALL_LIBDIR}/cmake/OpenMC)
 install(TARGETS openmc libopenmc
@@ -590,10 +625,24 @@ install(EXPORT openmc-targets
   NAMESPACE OpenMC::
   DESTINATION ${INSTALL_CONFIGDIR})
 
+configure_package_config_file(
+    "cmake/OpenMCConfig.cmake.in"
+    "${CMAKE_BINARY_DIR}/${CMAKE_FILES_DIRECTORY}/OpenMCConfig.cmake"
+    INSTALL_DESTINATION ${INSTALL_CONFIGDIR}
+)
+
+write_basic_package_version_file(
+    "${CMAKE_BINARY_DIR}/${CMAKE_FILES_DIRECTORY}/OpenMCConfigVersion.cmake"
+    VERSION ${OPENMC_VERSION}
+    COMPATIBILITY AnyNewerVersion
+)
+
 install(FILES
-  "${CMAKE_BINARY_DIR}${CMAKE_FILES_DIRECTORY}/OpenMCConfig.cmake"
-  "${CMAKE_BINARY_DIR}${CMAKE_FILES_DIRECTORY}/OpenMCConfigVersion.cmake"
-  DESTINATION ${INSTALL_CONFIGDIR})
+    "${CMAKE_BINARY_DIR}/${CMAKE_FILES_DIRECTORY}/OpenMCConfig.cmake"
+    "${CMAKE_BINARY_DIR}/${CMAKE_FILES_DIRECTORY}/OpenMCConfigVersion.cmake"
+    DESTINATION "${INSTALL_CONFIGDIR}"
+)
+
 install(FILES man/man1/openmc.1 DESTINATION ${CMAKE_INSTALL_MANDIR}/man1)
 install(FILES LICENSE DESTINATION "${CMAKE_INSTALL_DOCDIR}" RENAME copyright)
 install(DIRECTORY include/ DESTINATION ${CMAKE_INSTALL_INCLUDEDIR})

CODEOWNERS

@@ -5,9 +5,9 @@ openmc/data/ @paulromano
 openmc/lib/ @paulromano
 
 # Depletion
-openmc/deplete/ @drewejohnson
-tests/regression_tests/deplete/ @drewejohnson
-tests/unit_tests/test_deplete_*.py @drewejohnson
+openmc/deplete/ @paulromano
+tests/regression_tests/deplete/ @paulromano
+tests/unit_tests/test_deplete_*.py @paulromano
 
 # MG-related functionality
 openmc/mgxs_library.py @nelsonag
@@ -26,6 +26,12 @@ src/dagmc.cpp @pshriwise
 tests/regression_tests/dagmc/ @pshriwise
 tests/unit_tests/dagmc/ @pshriwise
 
+# Weight windows
+openmc/weight_windows.py @pshriwise
+openmc/lib/weight_windows.py @pshriwise
+src/weight_windows.py @pshriwise
+tests/unit_tests/weightwindows/ @pshriwise
+
 # Photon transport
 openmc/data/BREMX.DAT @amandalund
 openmc/data/compton_profiles.h5 @amandalund
@@ -49,3 +55,13 @@ openmc/data/resonance_covariance.py @icmeyer
 
 # Docker
 Dockerfile @shimwell
+
+# Random ray
+src/random_ray/ @jtramm
+
+# NCrystal interface
+src/ncrystal_interface.cpp @marquezj @tkittel
+src/ncrystal_load.cpp @marquezj @tkittel
+
+# MCPL interface
+src/mcpl_interface.cpp @ebknudsen

CONTRIBUTING.md

@@ -13,7 +13,7 @@ openmc@anl.gov.
 ## Resources
 
 - [GitHub Repository](https://github.com/openmc-dev/openmc)
-- [Documentation](http://docs.openmc.org/en/latest)
+- [Documentation](https://docs.openmc.org/en/latest)
 - [Discussion Forum](https://openmc.discourse.group)
 - [Slack Community](https://openmc.slack.com/signup) (If you don't see your
   domain listed, contact openmc@anl.gov)

Dockerfile

@@ -24,7 +24,7 @@ ARG compile_cores=1
 ARG build_dagmc=off
 ARG build_libmesh=off
 
-FROM debian:bullseye-slim AS dependencies
+FROM ubuntu:24.04 AS dependencies
 
 ARG compile_cores
 ARG build_dagmc
@@ -33,22 +33,17 @@ ARG build_libmesh
 # Set default value of HOME to /root
 ENV HOME=/root
 
-# Embree variables
-ENV EMBREE_TAG='v3.12.2'
-ENV EMBREE_REPO='https://github.com/embree/embree'
-ENV EMBREE_INSTALL_DIR=$HOME/EMBREE/
-
 # MOAB variables
-ENV MOAB_TAG='5.3.0'
+ENV MOAB_TAG='5.5.1'
 ENV MOAB_REPO='https://bitbucket.org/fathomteam/moab/'
 
 # Double-Down variables
-ENV DD_TAG='v1.0.0'
+ENV DD_TAG='v1.1.0'
 ENV DD_REPO='https://github.com/pshriwise/double-down'
 ENV DD_INSTALL_DIR=$HOME/Double_down
 
 # DAGMC variables
-ENV DAGMC_BRANCH='v3.2.1'
+ENV DAGMC_BRANCH='v3.2.4'
 ENV DAGMC_REPO='https://github.com/svalinn/DAGMC'
 ENV DAGMC_INSTALL_DIR=$HOME/DAGMC/
 
@@ -58,10 +53,11 @@ ENV LIBMESH_REPO='https://github.com/libMesh/libmesh'
 ENV LIBMESH_INSTALL_DIR=$HOME/LIBMESH
 
 # NJOY variables
+ENV NJOY_TAG='2016.78'
 ENV NJOY_REPO='https://github.com/njoy/NJOY2016'
 
 # Setup environment variables for Docker image
-ENV LD_LIBRARY_PATH=${DAGMC_INSTALL_DIR}/lib:$LD_LIBRARY_PATH \
+ENV LD_LIBRARY_PATH=${DAGMC_INSTALL_DIR}/lib:${LD_LIBRARY_PATH:-} \
     OPENMC_ENDF_DATA=/root/endf-b-vii.1 \
     DEBIAN_FRONTEND=noninteractive
 
@@ -71,15 +67,19 @@ RUN apt-get update -y && \
     apt-get install -y \
         python3-pip python-is-python3 wget git build-essential cmake \
         mpich libmpich-dev libhdf5-serial-dev libhdf5-mpich-dev \
-        libpng-dev && \
+        libpng-dev libpugixml-dev libfmt-dev catch2 python3-venv && \
     apt-get autoremove
 
+# create virtual enviroment to avoid externally managed environment error
+RUN python3 -m venv openmc_venv
+ENV PATH=/openmc_venv/bin:$PATH
+
 # Update system-provided pip
 RUN pip install --upgrade pip
 
 # Clone and install NJOY2016
 RUN cd $HOME \
-    && git clone --single-branch --depth 1 ${NJOY_REPO} \
+    && git clone --single-branch -b ${NJOY_TAG} --depth 1 ${NJOY_REPO} \
     && cd NJOY2016 \
     && mkdir build \
     && cd build \
@@ -90,28 +90,21 @@ RUN cd $HOME \
 
 RUN if [ "$build_dagmc" = "on" ]; then \
         # Install addition packages required for DAGMC
-        apt-get -y install libeigen3-dev libnetcdf-dev libtbb-dev libglfw3-dev \
-        && pip install --upgrade numpy "cython<3.0" \
-        # Clone and install EMBREE
-        && mkdir -p $HOME/EMBREE && cd $HOME/EMBREE \
-        && git clone --single-branch -b ${EMBREE_TAG} --depth 1 ${EMBREE_REPO} \
-        && mkdir build && cd build \
-        && cmake ../embree \
-                    -DCMAKE_INSTALL_PREFIX=${EMBREE_INSTALL_DIR} \
-                    -DEMBREE_MAX_ISA=NONE \
-                    -DEMBREE_ISA_SSE42=ON \
-                    -DEMBREE_ISPC_SUPPORT=OFF \
-        && make 2>/dev/null -j${compile_cores} install \
-        && rm -rf ${EMBREE_INSTALL_DIR}/build ${EMBREE_INSTALL_DIR}/embree ; \
+        apt-get -y install \
+            libeigen3-dev libnetcdf-dev libtbb-dev libglfw3-dev libembree-dev \
+        && pip install --upgrade numpy \
+        && pip install --no-cache-dir setuptools cython \
         # Clone and install MOAB
-        mkdir -p $HOME/MOAB && cd $HOME/MOAB \
+        && mkdir -p $HOME/MOAB && cd $HOME/MOAB \
         && git clone  --single-branch -b ${MOAB_TAG} --depth 1 ${MOAB_REPO} \
         && mkdir build && cd build \
-        && cmake ../moab -DENABLE_HDF5=ON \
+        && cmake ../moab -DCMAKE_BUILD_TYPE=Release \
+                      -DENABLE_HDF5=ON \
                       -DENABLE_NETCDF=ON \
                       -DBUILD_SHARED_LIBS=OFF \
                       -DENABLE_FORTRAN=OFF \
                       -DENABLE_BLASLAPACK=OFF \
+                      -DENABLE_TESTING=OFF \
         && make 2>/dev/null -j${compile_cores} install \
         && cmake ../moab \
                     -DENABLE_PYMOAB=ON \
@@ -127,7 +120,7 @@ RUN if [ "$build_dagmc" = "on" ]; then \
         && mkdir build && cd build \
         && cmake ../double-down -DCMAKE_INSTALL_PREFIX=${DD_INSTALL_DIR} \
                              -DMOAB_DIR=/usr/local \
-                             -DEMBREE_DIR=${EMBREE_INSTALL_DIR} \
+                             -DEMBREE_DIR=/usr \
         && make 2>/dev/null -j${compile_cores} install \
         && rm -rf ${DD_INSTALL_DIR}/build ${DD_INSTALL_DIR}/double-down ; \
         # Clone and install DAGMC
@@ -141,6 +134,7 @@ RUN if [ "$build_dagmc" = "on" ]; then \
                        -DDOUBLE_DOWN_DIR=${DD_INSTALL_DIR} \
                        -DCMAKE_PREFIX_PATH=${DD_INSTALL_DIR}/lib \
                        -DBUILD_STATIC_LIBS=OFF \
+                       -DBUILD_TESTS=OFF \
         && make 2>/dev/null -j${compile_cores} install \
         && rm -rf ${DAGMC_INSTALL_DIR}/DAGMC ${DAGMC_INSTALL_DIR}/build ; \
     fi
@@ -189,7 +183,7 @@ ENV LIBMESH_INSTALL_DIR=$HOME/LIBMESH
 
 # clone and install openmc
 RUN mkdir -p ${HOME}/OpenMC && cd ${HOME}/OpenMC \
-    && git clone --shallow-submodules --recurse-submodules --single-branch -b ${openmc_branch} --depth=1 ${OPENMC_REPO} \
+    && git clone --shallow-submodules --recurse-submodules --single-branch -b ${openmc_branch} ${OPENMC_REPO} \
     && mkdir build && cd build ; \
     if [ ${build_dagmc} = "on" ] && [ ${build_libmesh} = "on" ]; then \
         cmake ../openmc \

LICENSE

@@ -1,4 +1,4 @@
-Copyright (c) 2011-2023 Massachusetts Institute of Technology, UChicago Argonne
+Copyright (c) 2011-2026 Massachusetts Institute of Technology, UChicago Argonne
 LLC, and OpenMC contributors
 
 Permission is hereby granted, free of charge, to any person obtaining a copy of

MANIFEST.in

@@ -26,8 +26,6 @@ recursive-include include *.h
 recursive-include include *.h.in
 recursive-include include *.hh
 recursive-include man *.1
-recursive-include openmc *.pyx
-recursive-include openmc *.c
 recursive-include src *.cc
 recursive-include src *.cpp
 recursive-include src *.rnc
@@ -45,6 +43,5 @@ recursive-include vendor *.hh
 recursive-include vendor *.hpp
 recursive-include vendor *.pc.in
 recursive-include vendor *.natvis
-include vendor/gsl-lite/include/gsl/gsl
 prune docs/build
 prune docs/source/pythonapi/generated/

cmake/Modules/FindLIBMESH.cmake

@@ -14,8 +14,8 @@ if(DEFINED ENV{METHOD})
   message(STATUS "Using environment variable METHOD to determine libMesh build: ${LIBMESH_PC_FILE}")
 endif()
 
-include(FindPkgConfig)
-set(ENV{PKG_CONFIG_PATH} "$ENV{PKG_CONFIG_PATH}:${LIBMESH_PC}")
-set(PKG_CONFIG_USE_CMAKE_PREFIX_PATH True)
+find_package(PkgConfig REQUIRED)
+
+set(PKG_CONFIG_USE_CMAKE_PREFIX_PATH TRUE)
 pkg_check_modules(LIBMESH REQUIRED ${LIBMESH_PC_FILE}>=1.7.0 IMPORTED_TARGET)
 pkg_get_variable(LIBMESH_PREFIX ${LIBMESH_PC_FILE} prefix)

cmake/Modules/GetVersionFromGit.cmake

@@ -0,0 +1,120 @@
+# GetVersionFromGit.cmake
+# Standalone script to retrieve versioning information from Git or .git_archival.txt.
+# Customizable for any project by setting variables before including this file.
+
+# Configurable variables:
+#   - VERSION_PREFIX: Prefix for version tags (default: "v").
+#   - VERSION_SUFFIX: Suffix for version tags (default: "[~+-]([a-zA-Z0-9]+)").
+#   - VERSION_REGEX: Regex to extract version (default: "(?[0-9]+\\.[0-9]+\\.[0-9]+)").
+#   - ARCHIVAL_FILE: Path to .git_archival.txt (default: "${CMAKE_SOURCE_DIR}/.git_archival.txt").
+#   - DESCRIBE_NAME_KEY: Key for describe name in .git_archival.txt (default: "describe-name: ").
+#   - COMMIT_HASH_KEY: Key for commit hash in .git_archival.txt (default: "commit: ").
+
+# Default Format Example:
+#   1.2.3 v1.2.3 v1.2.3-rc1
+
+set(VERSION_PREFIX "v" CACHE STRING "Prefix used in version tags")
+set(VERSION_SUFFIX "[~+-]([a-zA-Z0-9]+)" CACHE STRING "Suffix used in version tags")
+set(VERSION_REGEX "?([0-9]+\\.[0-9]+\\.[0-9]+)" CACHE STRING "Regex for extracting version")
+set(ARCHIVAL_FILE "${CMAKE_SOURCE_DIR}/.git_archival.txt" CACHE STRING "Path to .git_archival.txt")
+set(DESCRIBE_NAME_KEY "describe-name: " CACHE STRING "Key for describe name in .git_archival.txt")
+set(COMMIT_HASH_KEY "commit: " CACHE STRING "Key for commit hash in .git_archival.txt")
+
+
+# Combine prefix and regex
+set(VERSION_REGEX_WITH_PREFIX "^${VERSION_PREFIX}${VERSION_REGEX}")
+
+# Find Git
+find_package(Git)
+
+# Attempt to retrieve version from Git
+if(EXISTS "${CMAKE_SOURCE_DIR}/.git" AND GIT_FOUND)
+    message(STATUS "Using git describe for versioning")
+
+    # Extract the version string
+    execute_process(
+        COMMAND git describe --tags --dirty
+        WORKING_DIRECTORY ${CMAKE_SOURCE_DIR}
+        OUTPUT_VARIABLE VERSION_STRING
+        OUTPUT_STRIP_TRAILING_WHITESPACE
+        ERROR_QUIET
+    )
+
+    # If no tags are found, set version to 0 and show a warning
+    if(VERSION_STRING STREQUAL "")
+        set(VERSION_STRING "0.0.0")
+        message(WARNING 
+            "No git tags found. Version set to 0.0.0.\n"
+            "Run 'git fetch --tags' to ensure proper versioning.\n"
+            "For more information, see OpenMC developer documentation."
+        )
+    endif()
+
+    # Extract the commit hash
+    execute_process(
+        COMMAND git rev-parse HEAD
+        WORKING_DIRECTORY ${CMAKE_SOURCE_DIR}
+        OUTPUT_VARIABLE COMMIT_HASH
+        OUTPUT_STRIP_TRAILING_WHITESPACE
+    )
+else()
+    message(STATUS "Using archival file for versioning: ${ARCHIVAL_FILE}")
+    if(EXISTS "${ARCHIVAL_FILE}")
+        file(READ "${ARCHIVAL_FILE}" ARCHIVAL_CONTENT)
+
+        # Extract the describe-name line
+        string(REGEX MATCH "${DESCRIBE_NAME_KEY}([^\\n]+)" VERSION_STRING "${ARCHIVAL_CONTENT}")
+        if(VERSION_STRING MATCHES "${DESCRIBE_NAME_KEY}(.*)")
+            set(VERSION_STRING "${CMAKE_MATCH_1}")
+        else()
+            message(FATAL_ERROR "Could not extract version from ${ARCHIVAL_FILE}")
+        endif()
+
+        # Extract the commit hash
+        string(REGEX MATCH "${COMMIT_HASH_KEY}([a-f0-9]+)" COMMIT_HASH "${ARCHIVAL_CONTENT}")
+        if(COMMIT_HASH MATCHES "${COMMIT_HASH_KEY}([a-f0-9]+)")
+            set(COMMIT_HASH "${CMAKE_MATCH_1}")
+        else()
+            message(FATAL_ERROR "Could not extract commit hash from ${ARCHIVAL_FILE}")
+        endif()
+    else()
+        message(FATAL_ERROR "Neither git describe nor ${ARCHIVAL_FILE} is available for versioning.")
+    endif()
+endif()
+
+# Ensure version string format
+if(VERSION_STRING MATCHES "${VERSION_REGEX_WITH_PREFIX}")
+    set(VERSION_NO_SUFFIX "${CMAKE_MATCH_1}")
+else()
+    message(FATAL_ERROR "Invalid version format: Missing base version in ${VERSION_STRING}")
+endif()
+
+# Check for development state
+if(VERSION_STRING MATCHES "-([0-9]+)-g([0-9a-f]+)")
+    set(DEV_STATE "true")
+    set(COMMIT_COUNT "${CMAKE_MATCH_1}")
+    string(REGEX REPLACE "-([0-9]+)-g([0-9a-f]+)" "" VERSION_WITHOUT_META "${VERSION_STRING}")
+else()
+    set(DEV_STATE "false")
+    set(VERSION_WITHOUT_META "${VERSION_STRING}")
+endif()
+
+# Split and set version components
+string(REPLACE "." ";" VERSION_LIST "${VERSION_NO_SUFFIX}")
+list(GET VERSION_LIST 0 VERSION_MAJOR)
+list(GET VERSION_LIST 1 VERSION_MINOR)
+list(GET VERSION_LIST 2 VERSION_PATCH)
+
+# Increment patch number for dev versions
+if(DEV_STATE)
+    math(EXPR VERSION_PATCH "${VERSION_PATCH} + 1")
+endif()
+
+# Export variables
+set(OPENMC_VERSION_MAJOR "${VERSION_MAJOR}")
+set(OPENMC_VERSION_MINOR "${VERSION_MINOR}")
+set(OPENMC_VERSION_PATCH "${VERSION_PATCH}")
+set(OPENMC_VERSION "${VERSION_MAJOR}.${VERSION_MINOR}.${VERSION_PATCH}")
+set(OPENMC_COMMIT_HASH "${COMMIT_HASH}")
+set(OPENMC_DEV_STATE "${DEV_STATE}")
+set(OPENMC_COMMIT_COUNT "${COMMIT_COUNT}")

cmake/OpenMCConfig.cmake.in

@@ -1,17 +1,18 @@
-get_filename_component(OpenMC_CMAKE_DIR "${CMAKE_CURRENT_LIST_FILE}" DIRECTORY)
+@PACKAGE_INIT@
 
-find_package(fmt REQUIRED HINTS ${OpenMC_CMAKE_DIR}/../fmt)
-find_package(gsl-lite REQUIRED HINTS ${OpenMC_CMAKE_DIR}/../gsl-lite)
-find_package(pugixml REQUIRED HINTS ${OpenMC_CMAKE_DIR}/../pugixml)
-find_package(xtl REQUIRED HINTS ${OpenMC_CMAKE_DIR}/../xtl)
-find_package(xtensor REQUIRED HINTS ${OpenMC_CMAKE_DIR}/../xtensor)
-if(@OPENMC_USE_DAGMC@)
-  find_package(DAGMC REQUIRED HINTS @DAGMC_DIR@)
+include("${CMAKE_CURRENT_LIST_DIR}/OpenMCConfigVersion.cmake")
+include(CMakeFindDependencyMacro)
+
+# Explicitly calculate prefix if it was not generated above
+if(NOT DEFINED PACKAGE_PREFIX_DIR)
+  get_filename_component(PACKAGE_PREFIX_DIR "${CMAKE_CURRENT_LIST_DIR}/../../.." ABSOLUTE)
 endif()
 
-if(@OPENMC_USE_NCRYSTAL@)
-  find_package(NCrystal REQUIRED)
-  message(STATUS "Found NCrystal: ${NCrystal_DIR} (version ${NCrystal_VERSION})")
+find_dependency(fmt CONFIG REQUIRED HINTS ${PACKAGE_PREFIX_DIR})
+find_dependency(pugixml CONFIG REQUIRED HINTS ${PACKAGE_PREFIX_DIR})
+
+if(@OPENMC_USE_DAGMC@)
+  find_dependency(DAGMC REQUIRED HINTS @DAGMC_DIR@)
 endif()
 
 if(@OPENMC_USE_LIBMESH@)
@@ -21,16 +22,24 @@ if(@OPENMC_USE_LIBMESH@)
   pkg_check_modules(LIBMESH REQUIRED @LIBMESH_PC_FILE@>=1.7.0 IMPORTED_TARGET)
 endif()
 
-find_package(PNG)
-
-if(NOT TARGET OpenMC::libopenmc)
-  include("${OpenMC_CMAKE_DIR}/OpenMCTargets.cmake")
+if("@PNG_FOUND@")
+  find_dependency(PNG)
 endif()
 
 if(@OPENMC_USE_MPI@)
-  find_package(MPI REQUIRED)
+  find_dependency(MPI REQUIRED)
 endif()
 
-if(@OPENMC_USE_MCPL@)
-  find_package(MCPL REQUIRED)
+if(@OPENMC_USE_OPENMP@)
+  find_dependency(OpenMP REQUIRED)
+endif()
+
+if(@OPENMC_USE_UWUW@ AND NOT ${DAGMC_BUILD_UWUW})
+  message(FATAL_ERROR "UWUW is enabled in OpenMC but the DAGMC installation discovered was not configured with UWUW.")
+endif()
+
+include("${CMAKE_CURRENT_LIST_DIR}/OpenMCTargets.cmake")
+
+if(NOT OpenMC_FIND_QUIETLY)
+    message(STATUS "Found OpenMC: ${PACKAGE_VERSION} (found in ${PACKAGE_PREFIX_DIR})")
 endif()

cmake/OpenMCConfigVersion.cmake.in

@@ -1,11 +0,0 @@
-set(PACKAGE_VERSION "@OPENMC_VERSION@")
-
-# Check whether the requested PACKAGE_FIND_VERSION is compatible
-if("${PACKAGE_VERSION}" VERSION_LESS "${PACKAGE_FIND_VERSION}")
-  set(PACKAGE_VERSION_COMPATIBLE FALSE)
-else()
-  set(PACKAGE_VERSION_COMPATIBLE TRUE)
-  if ("${PACKAGE_VERSION}" VERSION_EQUAL "${PACKAGE_FIND_VERSION}")
-    set(PACKAGE_VERSION_EXACT TRUE)
-  endif()
-endif()

docs/Makefile

@@ -45,6 +45,7 @@ help:
 clean:
 	-rm -rf $(BUILDDIR)/*
 	-rm -rf source/pythonapi/generated/
+	-rm -rf doxygen/xml
 
 html:
 	$(SPHINXBUILD) -b html $(ALLSPHINXOPTS) $(BUILDDIR)/html

docs/doxygen/Doxyfile

@@ -0,0 +1,13 @@
+# Doxyfile 1.9.1
+
+# This file describes the settings to be used by the documentation system
+# doxygen (www.doxygen.org) for a project.
+
+# Difference with default Doxyfile 1.9.1
+PROJECT_NAME           = OpenMC
+QUIET                  = YES
+WARN_IF_UNDOCUMENTED   = NO
+INPUT                  = ../../include/openmc/capi.h
+GENERATE_HTML          = NO
+GENERATE_LATEX         = NO
+GENERATE_XML           = YES

docs/requirements-rtd.txt

@@ -1,13 +0,0 @@
-sphinx==5.0.2
-sphinx_rtd_theme==1.0.0
-sphinx-numfig
-jupyter
-sphinxcontrib-katex
-sphinxcontrib-svg2pdfconverter
-numpy
-scipy
-h5py
-pandas
-uncertainties
-matplotlib
-lxml

docs/source/capi/index.rst

@@ -46,41 +46,29 @@ Type Definitions
 Functions
 ---------
 
-.. c:function:: int openmc_calculate_volumes()
+..
+   Once documentation is complete in capi.h, use:
+   .. doxygenfile:: capi.h
+   to populate this documentation without using
+   .. doxygenfunction::
+   for every function.
 
-   Run a stochastic volume calculation
+.. doxygenfunction:: openmc_calculate_volumes
 
-   :return: Return status (negative if an error occurred)
-   :rtype: int
+.. doxygenfunction:: openmc_cell_get_fill
 
-.. c:function:: int openmc_cell_get_fill(int32_t index, int* type, int32_t** indices, int32_t* n)
+.. doxygenfunction:: openmc_cell_get_id
 
-   Get the fill for a cell
+.. doxygenfunction:: openmc_cell_get_temperature
+
+.. c:function:: int openmc_cell_get_density(int32_t index, const int32_t* instance, double* density)
+
+   Get the density of a cell
 
    :param int32_t index: Index in the cells array
-   :param int* type: Type of the fill
-   :param int32_t** indices: Array of material indices for cell
-   :param int32_t* n: Length of indices array
-   :return: Return status (negative if an error occurred)
-   :rtype: int
-
-.. c:function:: int openmc_cell_get_id(int32_t index, int32_t* id)
-
-   Get the ID of a cell
-
-   :param int32_t index: Index in the cells array
-   :param int32_t* id: ID of the cell
-   :return: Return status (negative if an error occurred)
-   :rtype: int
-
-.. c:function:: int openmc_cell_get_temperature(int32_t index, const int32_t* instance, double* T)
-
-   Get the temperature of a cell
-
-   :param int32_t index: Index in the cells array
-   :param int32_t* instance: Which instance of the cell. If a null pointer is passed, the temperature
-                             of the first instance is returned.
-   :param double* T: temperature of the cell
+   :param int32_t* instance: Which instance of the cell. If a null pointer is passed, the density
+                             multiplier of the first instance is returned.
+   :param double* density: Density of the cell in [g/cm3]
    :return: Return status (negative if an error occurred)
    :rtype: int
 
@@ -113,8 +101,22 @@ Functions
    :param double T: Temperature in Kelvin
    :param instance: Which instance of the cell. To set the temperature for all
                     instances, pass a null pointer.
-   :param set_contained: If the cell is not filled by a material, whether to set the temperatures
-                         of all filled cells
+   :param bool set_contained: If the cell is not filled by a material, whether
+                              to set the temperatures of all filled cells
+   :type instance: const int32_t*
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_cell_set_density(index index, double density, const int32_t* instance, bool set_contained)
+
+   Set the density of a cell.
+
+   :param int32_t index: Index in the cells array
+   :param double density: Density of the cell in [g/cm3]
+   :param instance: Which instance of the cell. To set the density multiplier for all
+                    instances, pass a null pointer.
+   :param bool set_contained: If the cell is not filled by a material, whether
+                              to set the density multiplier of all filled cells
    :type instance: const int32_t*
    :return: Return status (negative if an error occurred)
    :rtype: int
@@ -420,6 +422,16 @@ Functions
    :return: Return status (negative if an error occurred)
    :rtype: int
 
+.. c:function:: int openmc_mesh_filter_get_mesh(int32_t index, int32_t* index_mesh)
+
+   Get the mesh for a mesh filter
+
+   :param int32_t index: Index in the filters array
+   :param index_mesh: Index in the meshes array
+   :type index_mesh: int32_t*
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
 .. c:function:: int openmc_mesh_filter_set_mesh(int32_t index, int32_t index_mesh)
 
    Set the mesh for a mesh filter
@@ -429,6 +441,98 @@ Functions
    :return: Return status (negative if an error occurred)
    :rtype: int
 
+.. c:function:: int openmc_mesh_filter_get_translation(int32_t index, double translation[3])
+
+   Get the 3-D translation coordinates for a mesh filter
+
+   :param int32_t index: Index in the filters array
+   :param double[3] translation: 3-D translation coordinates
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_mesh_filter_set_translation(int32_t index, double translation[3])
+
+   Set the 3-D translation coordinates for a mesh filter
+
+   :param int32_t index: Index in the filters array
+   :param double[3] translation: 3-D translation coordinates
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_meshborn_filter_get_mesh(int32_t index, int32_t* index_mesh)
+
+   Get the mesh for a meshborn filter
+
+   :param int32_t index: Index in the filters array
+   :param index_mesh: Index in the meshes array
+   :type index_mesh: int32_t*
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_meshborn_filter_set_mesh(int32_t index, int32_t index_mesh)
+
+   Set the mesh for a meshborn filter
+
+   :param int32_t index: Index in the filters array
+   :param int32_t index_mesh: Index in the meshes array
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_meshborn_filter_get_translation(int32_t index, double translation[3])
+
+   Get the 3-D translation coordinates for a meshborn filter
+
+   :param int32_t index: Index in the filters array
+   :param double[3] translation: 3-D translation coordinates
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_meshborn_filter_set_translation(int32_t index, double translation[3])
+
+   Set the 3-D translation coordinates for a meshborn filter
+
+   :param int32_t index: Index in the filters array
+   :param double[3] translation: 3-D translation coordinates
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_meshsurface_filter_get_mesh(int32_t index, int32_t* index_mesh)
+
+   Get the mesh for a mesh surface filter
+
+   :param int32_t index: Index in the filters array
+   :param index_mesh: Index in the meshes array
+   :type index_mesh: int32_t*
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_meshsurface_filter_set_mesh(int32_t index, int32_t index_mesh)
+
+   Set the mesh for a mesh surface filter
+
+   :param int32_t index: Index in the filters array
+   :param int32_t index_mesh: Index in the meshes array
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_meshsurface_filter_get_translation(int32_t index, double translation[3])
+
+   Get the 3-D translation coordinates for a mesh surface filter
+
+   :param int32_t index: Index in the filters array
+   :param double[3] translation: 3-D translation coordinates
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_meshsurface_filter_set_translation(int32_t index, double translation[3])
+
+   Set the 3-D translation coordinates for a mesh surface filter
+
+   :param int32_t index: Index in the filters array
+   :param double[3] translation: 3-D translation coordinates
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
 .. c:function:: int openmc_next_batch()
 
    Simulate next batch of particles. Must be called after openmc_simulation_init().
@@ -453,6 +557,279 @@ Functions
    :return: Return status (negative if an error occurs)
    :rtype: int
 
+.. c:function:: int openmc_get_plot_index(int32_t id, int32_t* index)
+
+   Get the index in the plots array for a plot with a given ID.
+
+   :param int32_t id: Plot ID
+   :param int32_t* index: Index in the plots array
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_plot_get_id(int32_t index, int32_t* id)
+
+   Get the ID of a plot.
+
+   :param int32_t index: Index in the plots array
+   :param int32_t* id: Plot ID
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_plot_set_id(int32_t index, int32_t id)
+
+   Set the ID of a plot.
+
+   :param int32_t index: Index in the plots array
+   :param int32_t id: Plot ID
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: size_t openmc_plots_size()
+
+   Number of plots currently allocated.
+
+   :return: Number of plots in the plots array
+   :rtype: size_t
+
+.. c:function:: int openmc_solidraytrace_plot_create(int32_t* index)
+
+   Create a new solid raytrace plot.
+
+   :param int32_t* index: Index of the newly created plot
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_solidraytrace_plot_get_pixels(int32_t index, int32_t* width, int32_t* height)
+
+   Get output pixel dimensions for a solid raytrace plot.
+
+   :param int32_t index: Index in the plots array
+   :param int32_t* width: Image width in pixels
+   :param int32_t* height: Image height in pixels
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_solidraytrace_plot_set_pixels(int32_t index, int32_t width, int32_t height)
+
+   Set output pixel dimensions for a solid raytrace plot.
+
+   :param int32_t index: Index in the plots array
+   :param int32_t width: Image width in pixels
+   :param int32_t height: Image height in pixels
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_solidraytrace_plot_get_color_by(int32_t index, int32_t* color_by)
+
+   Get the domain type used for coloring (0=materials, 1=cells).
+
+   :param int32_t index: Index in the plots array
+   :param int32_t* color_by: Coloring mode
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_solidraytrace_plot_set_color_by(int32_t index, int32_t color_by)
+
+   Set the domain type used for coloring (0=materials, 1=cells).
+
+   :param int32_t index: Index in the plots array
+   :param int32_t color_by: Coloring mode
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_solidraytrace_plot_set_default_colors(int32_t index)
+
+   Set default random colors for the current ``color_by`` mode.
+
+   :param int32_t index: Index in the plots array
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_solidraytrace_plot_set_all_opaque(int32_t index)
+
+   Mark all domains in the current ``color_by`` mode as opaque.
+
+   :param int32_t index: Index in the plots array
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_solidraytrace_plot_set_opaque(int32_t index, int32_t id, bool visible)
+
+   Set whether a specific domain ID is opaque (visible) in the rendered image.
+
+   :param int32_t index: Index in the plots array
+   :param int32_t id: Cell/material ID (based on ``color_by``)
+   :param bool visible: Whether the domain is opaque
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_solidraytrace_plot_set_color(int32_t index, int32_t id, uint8_t r, uint8_t g, uint8_t b)
+
+   Set RGB color for a specific domain ID.
+
+   :param int32_t index: Index in the plots array
+   :param int32_t id: Cell/material ID (based on ``color_by``)
+   :param uint8_t r: Red channel
+   :param uint8_t g: Green channel
+   :param uint8_t b: Blue channel
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_solidraytrace_plot_get_color(int32_t index, int32_t id, uint8_t* r, uint8_t* g, uint8_t* b)
+
+   Get RGB color for a specific domain ID.
+
+   :param int32_t index: Index in the plots array
+   :param int32_t id: Cell/material ID (based on ``color_by``)
+   :param uint8_t* r: Red channel
+   :param uint8_t* g: Green channel
+   :param uint8_t* b: Blue channel
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_solidraytrace_plot_get_camera_position(int32_t index, double* x, double* y, double* z)
+
+   Get camera position.
+
+   :param int32_t index: Index in the plots array
+   :param double* x: X coordinate
+   :param double* y: Y coordinate
+   :param double* z: Z coordinate
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_solidraytrace_plot_set_camera_position(int32_t index, double x, double y, double z)
+
+   Set camera position.
+
+   :param int32_t index: Index in the plots array
+   :param double x: X coordinate
+   :param double y: Y coordinate
+   :param double z: Z coordinate
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_solidraytrace_plot_get_look_at(int32_t index, double* x, double* y, double* z)
+
+   Get camera target point.
+
+   :param int32_t index: Index in the plots array
+   :param double* x: X coordinate
+   :param double* y: Y coordinate
+   :param double* z: Z coordinate
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_solidraytrace_plot_set_look_at(int32_t index, double x, double y, double z)
+
+   Set camera target point.
+
+   :param int32_t index: Index in the plots array
+   :param double x: X coordinate
+   :param double y: Y coordinate
+   :param double z: Z coordinate
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_solidraytrace_plot_get_up(int32_t index, double* x, double* y, double* z)
+
+   Get the camera up vector.
+
+   :param int32_t index: Index in the plots array
+   :param double* x: X component
+   :param double* y: Y component
+   :param double* z: Z component
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_solidraytrace_plot_set_up(int32_t index, double x, double y, double z)
+
+   Set the camera up vector.
+
+   :param int32_t index: Index in the plots array
+   :param double x: X component
+   :param double y: Y component
+   :param double z: Z component
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_solidraytrace_plot_get_light_position(int32_t index, double* x, double* y, double* z)
+
+   Get light source position.
+
+   :param int32_t index: Index in the plots array
+   :param double* x: X coordinate
+   :param double* y: Y coordinate
+   :param double* z: Z coordinate
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_solidraytrace_plot_set_light_position(int32_t index, double x, double y, double z)
+
+   Set light source position.
+
+   :param int32_t index: Index in the plots array
+   :param double x: X coordinate
+   :param double y: Y coordinate
+   :param double z: Z coordinate
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_solidraytrace_plot_get_fov(int32_t index, double* fov)
+
+   Get horizontal field of view in degrees.
+
+   :param int32_t index: Index in the plots array
+   :param double* fov: Field of view in degrees
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_solidraytrace_plot_set_fov(int32_t index, double fov)
+
+   Set horizontal field of view in degrees.
+
+   :param int32_t index: Index in the plots array
+   :param double fov: Field of view in degrees
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_solidraytrace_plot_get_diffuse_fraction(int32_t index, double* diffuse_fraction)
+
+   Get diffuse-light fraction.
+
+   :param int32_t index: Index in the plots array
+   :param double* diffuse_fraction: Diffuse fraction in [0, 1]
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_solidraytrace_plot_set_diffuse_fraction(int32_t index, double diffuse_fraction)
+
+   Set diffuse-light fraction.
+
+   :param int32_t index: Index in the plots array
+   :param double diffuse_fraction: Diffuse fraction in [0, 1]
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_solidraytrace_plot_update_view(int32_t index)
+
+   Recompute internal camera/view transforms after camera changes.
+
+   :param int32_t index: Index in the plots array
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
+.. c:function:: int openmc_solidraytrace_plot_create_image(int32_t index, uint8_t* data_out, int32_t width, int32_t height)
+
+   Render the plot to an RGB image buffer.
+
+   :param int32_t index: Index in the plots array
+   :param uint8_t* data_out: Output buffer of shape ``height*width*3``
+   :param int32_t width: Image width in pixels
+   :param int32_t height: Image height in pixels
+   :return: Return status (negative if an error occurred)
+   :rtype: int
+
 .. c:function:: int openmc_reset()
 
    Resets all tally scores
@@ -474,6 +851,10 @@ Functions
    :return: Return status (negative if an error occurs)
    :rtype: int
 
+.. c:function:: void openmc_run_random_ray()
+
+   Run a random ray simulation
+
 .. c:function:: int openmc_set_n_batches(int32_t n_batches, bool set_max_batches, bool add_statepoint_batch)
 
    Set number of batches and number of max batches

docs/source/conf.py

@@ -11,7 +11,10 @@
 # All configuration values have a default; values that are commented out
 # serve to show the default.
 
-import sys, os
+import os
+from pathlib import Path
+import subprocess
+import sys
 
 # Determine if we're on Read the Docs server
 on_rtd = os.environ.get('READTHEDOCS', None) == 'True'
@@ -37,6 +40,7 @@ sys.path.insert(0, os.path.abspath('../..'))
 # Add any Sphinx extension module names here, as strings. They can be extensions
 # coming with Sphinx (named 'sphinx.ext.*') or your custom ones.
 extensions = [
+    'breathe',
     'sphinx.ext.autodoc',
     'sphinx.ext.napoleon',
     'sphinx.ext.autosummary',
@@ -47,12 +51,14 @@ extensions = [
 ]
 if not on_rtd:
     extensions.append('sphinxcontrib.rsvgconverter')
+    doxygen_dir = Path(__file__).parents[1] / 'doxygen'
+    subprocess.run(['doxygen'], cwd=doxygen_dir, check=True)
 
 # Add any paths that contain templates here, relative to this directory.
 templates_path = ['_templates']
 
 # The suffix of source filenames.
-source_suffix = '.rst'
+source_suffix = {'.rst': 'restructuredtext'}
 
 # The encoding of source files.
 #source_encoding = 'utf-8'
@@ -62,16 +68,17 @@ master_doc = 'index'
 
 # General information about the project.
 project = 'OpenMC'
-copyright = '2011-2023, Massachusetts Institute of Technology, UChicago Argonne LLC, and OpenMC contributors'
+copyright = '2011-2026, Massachusetts Institute of Technology, UChicago Argonne LLC, and OpenMC contributors'
 
 # The version info for the project you're documenting, acts as replacement for
 # |version| and |release|, also used in various other places throughout the
 # built documents.
 #
-# The short X.Y version.
-version = "0.14"
+
+import openmc
+
 # The full version, including alpha/beta/rc tags.
-release = "0.14.0"
+version = release = openmc.__version__
 
 # The language for content autogenerated by Sphinx. Refer to documentation
 # for a list of supported languages.
@@ -116,13 +123,17 @@ pygments_style = 'tango'
 # A list of ignored prefixes for module index sorting.
 #modindex_common_prefix = []
 
+# -- Options breathe + doxygen -------------------------------------------------
+
+breathe_projects = {"OpenMC": "../doxygen/xml"}
+breathe_default_project = "OpenMC"
+breathe_domain_by_file_pattern = {"*capi.h": "c"}
 
 # -- Options for HTML output ---------------------------------------------------
 
 # The theme to use for HTML and HTML Help pages
-import sphinx_rtd_theme
 html_theme = 'sphinx_rtd_theme'
-html_theme_path = [sphinx_rtd_theme.get_html_theme_path()]
+html_baseurl = "https://docs.openmc.org/en/stable/"
 
 html_logo = '_images/openmc_logo.png'
 

docs/source/devguide/agentic-tools.rst

@@ -0,0 +1,104 @@
+.. _devguide_agentic_tools:
+
+===========================
+Agentic Development Tools
+===========================
+
+OpenMC ships a set of tools designed for AI coding agents (such as
+`Claude Code`_) that agents can use to navigate and understand the codebase.
+
+.. _Claude Code: https://claude.ai/code
+
+Motivation
+----------
+
+Agentic tools like Claude Code are skilled at using grep to navigate and
+understand large code bases. However, grep can only find exact text matches —
+it cannot discover code that is *conceptually* related but uses different
+naming. Without a "global view" of the codebase that a human developer will
+build up over time, the agent is generally blind to any file it hasn't
+tokenized fully. While it can grep to see who else calls a function, it
+remains blind if other areas might be related but not share identical naming
+conventions.
+
+This problem is mitigated somewhat by using a model with a longer context
+window. OpenMC has somewhere around ~1 million tokens of C++ and ~1 million
+tokens of python. While Claude Code in early 2026 only has a context window
+of 200k tokens, beta versions have extended context windows of 1M tokens,
+and it's not unreasonable to assume that models may be available in the near
+future that greatly exceed these limits.
+
+However, even assuming the entire repository can be fit within a context
+window, there are several downsides to doing this.
+`Model performance degrades significantly as context size increases`_.
+Benchmark results are
+greatly improved if the model has less garbage to pick through. Additionally, API usage
+is typically billed as tokens in/out per turn. As the context file
+grows these costs become much larger. As such, there is still significant
+motivation to solving the above problem, so as to ensure only relevant
+information is drawn into context so as to maximize model performance and
+minimize costs.
+
+Setup
+-----
+
+The tools are registered as an `MCP (Model Context Protocol)`_ server in
+``.mcp.json`` at the repository root. AI agents that support MCP (such as
+Claude Code) discover them automatically on session start. The underlying
+Python scripts can also be run directly from the command line.
+
+All tools run entirely locally — no API keys or external service accounts are
+required. Python dependencies are installed automatically into an isolated
+virtual environment at ``.claude/cache/.venv/`` on first use.
+
+.. _Model performance degrades significantly as context size increases: https://www.anthropic.com/news/claude-opus-4-6
+.. _MCP (Model Context Protocol): https://modelcontextprotocol.io
+
+RAG Semantic Search
+-------------------
+
+The RAG (Retrieval-Augmented Generation) semantic search addresses this
+problem — it finds code by meaning, not just text match, surfacing related code
+across subsystems that ``grep`` would miss entirely. Two MCP tools are provided:
+
+- **openmc_rag_search** — Given a natural-language query, returns the most
+  relevant code chunks with file paths, line numbers, and a preview. Can search
+  code, documentation, or both. Can also find code related to a given file.
+- **openmc_rag_rebuild** — Rebuilds the search index. Should be called after
+  pulling new code or switching branches.
+
+How it works
+^^^^^^^^^^^^
+
+The search pipeline runs entirely on your local CPU:
+
+1. **Chunking.** All C++, Python, and RST files are split into overlapping
+   fixed-size windows (~1000 characters, 25% overlap). This ensures every line
+   of code appears in at least one chunk and most lines appear in two.
+
+2. **Embedding.** Each chunk is embedded into a 384-dimensional vector using
+   the `all-MiniLM-L6-v2`_ sentence-transformer model (22 million parameters).
+   This model runs on CPU with no GPU required. No API key is needed — the
+   model weights are downloaded once from Hugging Face and cached locally.
+
+3. **Indexing.** The vectors are stored in a local LanceDB_ database on disk.
+   Building the full index takes approximately 5 minutes on a machine with
+   10 CPU cores. The index is stored in ``.claude/cache/rag_index/`` and
+   persists across sessions.
+
+4. **Searching.** Your query is embedded using the same model, and the closest
+   chunks are retrieved by vector similarity. Results include the file path,
+   line range, file type, similarity distance, and a text preview.
+
+.. _all-MiniLM-L6-v2: https://huggingface.co/sentence-transformers/all-MiniLM-L6-v2
+.. _LanceDB: https://lancedb.com
+
+Requirements
+^^^^^^^^^^^^
+
+No system dependencies beyond **Python 3.12+** with ``pip``. An internet
+connection is required on first use to download the Python packages and
+embedding model weights; subsequent runs are fully offline. The Python packages
+(``sentence-transformers``, ``lancedb``) and their dependencies (including
+PyTorch, ~2GB) are installed automatically into an isolated virtual environment
+on first use.

docs/source/devguide/contributing.rst

@@ -109,9 +109,10 @@ Leadership Team
 The TC consists of the following individuals:
 
 - `Paul Romano <https://github.com/paulromano>`_
-- `Sterling Harper <https://github.com/smharper>`_
+- `Patrick Shriwise <https://github.com/pshriwise>`_
 - `Adam Nelson <https://github.com/nelsonag>`_
-- `Benoit Forget <https://github.com/bforget>`_
+- `Jonathan Shimwell <https://github.com/shimwell>`_
+- `John Tramm <https://github.com/jtramm>`_
 
 The Project Lead is Paul Romano.
 

docs/source/devguide/docbuild.rst

@@ -12,7 +12,12 @@ Python API. That is, from the root directory of the OpenMC repository:
 
 .. code-block:: sh
 
-    python -m pip install .[docs]
+    python -m pip install ".[docs]"
+
+The OpenMC documentation also uses Doxygen to automatically generate its
+C/C++ API documentation directly from the docstrings available in the source
+code. You will need to have a working installation of Doxygen to generate the
+documentation locally.
 
 -----------------------------------
 Building Documentation as a Webpage

docs/source/devguide/docker.rst

@@ -45,12 +45,11 @@ Now you can run the following to create a `Docker container`_ called
 This command will open an interactive shell running from within the
 Docker container where you have access to use OpenMC.
 
-.. note:: The ``docker run`` command supports many
-          `options <https://docs.docker.com/engine/reference/commandline/run/>`_
+.. note:: The ``docker run`` command supports many options_
           for spawning containers -- including `mounting volumes`_ from the
           host filesystem -- which many users will find useful.
 
-.. _Docker image: https://docs.docker.com/engine/reference/commandline/images/
+.. _Docker image: https://docs.docker.com/get-started/docker-concepts/the-basics/what-is-an-image/
 .. _Docker container: https://www.docker.com/resources/what-container
-.. _options: https://docs.docker.com/engine/reference/commandline/run/
-.. _mounting volumes: https://docs.docker.com/storage/volumes/
+.. _options: https://docs.docker.com/reference/cli/docker/container/run/
+.. _mounting volumes: https://docs.docker.com/engine/storage/volumes/

docs/source/devguide/index.rst

@@ -14,7 +14,9 @@ other related topics.
 
     contributing
     workflow
+    agentic-tools
     styleguide
+    policies
     tests
     user-input
     docbuild

docs/source/devguide/policies.rst

@@ -0,0 +1,35 @@
+.. _devguide_policies:
+
+========
+Policies
+========
+
+---------------------
+Python Version Policy
+---------------------
+
+OpenMC follows the Scientific Python Ecosystem Coordination guidelines `SPEC 0
+<https://scientific-python.org/specs/spec-0000/>`_ on minimum supported
+versions, which recommends that support for Python versions be dropped 3 years
+after their initial release.
+
+-------------------
+C++ Standard Policy
+-------------------
+
+C++ code in OpenMC must conform to the most recent C++ standard that is fully
+supported in the `version of the gcc compiler
+<https://gcc.gnu.org/projects/cxx-status.html>`_ that is distributed with the
+oldest version of Ubuntu that is still within its `standard support period
+<https://ubuntu.com/about/release-cycle>`_. Ubuntu 22.04 LTS will be supported
+through April 2027 and is distributed with gcc 11.4.0, which fully supports the
+C++17 standard.
+
+--------------------
+CMake Version Policy
+--------------------
+
+Similar to the C++ standard policy, the minimum supported version of CMake
+corresponds to whatever version is distributed with the oldest version of Ubuntu
+still within its standard support period. Ubuntu 22.04 LTS is distributed with
+CMake 3.22.

docs/source/devguide/styleguide.rst

@@ -30,7 +30,7 @@ whenever a file is saved. For example, `Visual Studio Code
 support for running clang-format.
 
 .. note::
-    OpenMC's CI uses `clang-format` version 15. A different version of `clang-format`
+    OpenMC's CI uses `clang-format` version 18. A different version of `clang-format`
     may produce different line changes and as a result fail the CI test.
 
 Miscellaneous
@@ -40,14 +40,15 @@ Follow the `C++ Core Guidelines`_ except when they conflict with another
 guideline listed here. For convenience, many important guidelines from that
 list are repeated here.
 
-Conform to the C++14 standard.
+Conform to the C++17 standard.
 
 Always use C++-style comments (``//``) as opposed to C-style (``/**/``). (It
 is more difficult to comment out a large section of code that uses C-style
 comments.)
 
 Do not use C-style casting. Always use the C++-style casts ``static_cast``,
-``const_cast``, or ``reinterpret_cast``. (See `ES.49 <http://isocpp.github.io/CppCoreGuidelines/CppCoreGuidelines#es49-if-you-must-use-a-cast-use-a-named-cast>`_)
+``const_cast``, or ``reinterpret_cast``. (See `ES.49
+<https://isocpp.github.io/CppCoreGuidelines/CppCoreGuidelines#es49-if-you-must-use-a-cast-use-a-named-cast>`_)
 
 Source Files
 ------------
@@ -55,7 +56,7 @@ Source Files
 Use a ``.cpp`` suffix for code files and ``.h`` for header files.
 
 Header files should always use include guards with the following style (See
-`SF.8 <http://isocpp.github.io/CppCoreGuidelines/CppCoreGuidelines#sf8-use-include-guards-for-all-h-files>`_):
+`SF.8 <https://isocpp.github.io/CppCoreGuidelines/CppCoreGuidelines#Rs-guards>`_):
 
 .. code-block:: C++
 
@@ -146,7 +147,7 @@ Style for Python code should follow PEP8_.
 
 Docstrings for functions and methods should follow numpydoc_ style.
 
-Python code should work with Python 3.7+.
+Python code should work with Python 3.8+.
 
 Use of third-party Python packages should be limited to numpy_, scipy_,
 matplotlib_, pandas_, and h5py_. Use of other third-party packages must be
@@ -156,11 +157,11 @@ Prefer pathlib_ when working with filesystem paths over functions in the os_
 module or other standard-library modules. Functions that accept arguments that
 represent a filesystem path should work with both strings and Path_ objects.
 
-.. _C++ Core Guidelines: http://isocpp.github.io/CppCoreGuidelines/CppCoreGuidelines
-.. _PEP8: https://www.python.org/dev/peps/pep-0008/
+.. _C++ Core Guidelines: https://isocpp.github.io/CppCoreGuidelines/CppCoreGuidelines
+.. _PEP8: https://peps.python.org/pep-0008/
 .. _numpydoc: https://numpydoc.readthedocs.io/en/latest/format.html
 .. _numpy: https://numpy.org/
-.. _scipy: https://www.scipy.org/
+.. _scipy: https://scipy.org/
 .. _matplotlib: https://matplotlib.org/
 .. _pandas: https://pandas.pydata.org/
 .. _h5py: https://www.h5py.org/

docs/source/devguide/tests.rst

@@ -37,6 +37,9 @@ Prerequisites
 - Some tests require `NJOY <https://www.njoy21.io/NJOY2016>`_ to preprocess
   cross section data. The test suite assumes that you have an ``njoy``
   executable available on your :envvar:`PATH`.
+- OpenMC should be compiled with ``-DOPENMC_ENABLE_STRICT_FP=on`` to ensure
+  reproducible floating-point results across platforms and optimization levels.
+  Without this flag, regression tests may not match reference values.
 
 Running Tests
 -------------
@@ -67,9 +70,11 @@ make sure you have satisfied all the prerequisites above. After you have done
 that, consider the following:
 
 - When building OpenMC, make sure you run CMake with
-  ``-DCMAKE_BUILD_TYPE=Debug``. Building with a release build will result in
-  some test failures due to differences in which compiler optimizations are
-  used.
+  ``-DOPENMC_ENABLE_STRICT_FP=on``. This prevents the compiler from applying
+  floating-point optimizations (such as replacing math library calls with
+  builtins or contracting multiply-add into FMA instructions) that can produce
+  bit-level differences across platforms and optimization levels. Any
+  ``CMAKE_BUILD_TYPE`` can be used.
 - Because tallies involve the sum of many floating point numbers, the
   non-associativity of floating point numbers can result in different answers
   especially when the number of threads is high (different order of operations).
@@ -84,6 +89,30 @@ that, consider the following:
   limit the number of threads that OpenBLAS uses internally; this can be done by
   setting the :envvar:`OPENBLAS_NUM_THREADS` environment variable to 1.
 
+Debugging Tests in CI
+---------------------
+
+Tests can be debugged in CI using a feature called 
+`tmate <https://github.com/mxschmitt/action-tmate?tab=readme-ov-file#debug-your-github-actions-by-using-tmate>`_. 
+CI debugging can be
+enabled by including "[gha-debug]" in the commit message. When the test fails, a
+link similar to the one shown below will be provided in the GitHub Actions
+output after failure occurs. Logging into the provided link will allow you to
+debug the test in the CI environment. The following is an example of the output
+shown in the CI log that provides the link to the tmate session:
+
+.. code-block:: text
+   :linenos:
+
+   Created new session successfully
+   ssh 2VcykjU7vNdvAzEjQcc839GM2@nyc1.tmate.io
+   https://tmate.io/t/2VcykjU7vNdvAzEjQcc839GM2
+   Entering main loop
+   Web shell: https://tmate.io/t/2VcykjU7vNdvAzEjQcc839GM2
+   SSH: ssh 2VcykjU7vNdvAzEjQcc839GM2@nyc1.tmate.io
+   ...
+
+
 Generating XML Inputs
 ---------------------
 
@@ -134,6 +163,12 @@ following files to your new test directory:
       compiler options during openmc configuration and build (e.g., no MPI, no
       debug/optimization).
 
+For tests using the Python API, both the **inputs_true.dat** and
+**results_true.dat** files can be generated automatically in the correct format
+via::
+
+    pytest --update <name-of-test>
+
 In addition to this description, please see the various types of tests that are
 already included in the test suite to see how to create them. If all is
 implemented correctly, the new test will automatically be discovered by pytest.

docs/source/devguide/user-input.rst

@@ -55,6 +55,6 @@ developer or send a message to the `developers mailing list`_.
 
 
 .. _property attribute: https://docs.python.org/3.6/library/functions.html#property
-.. _XML Schema Part 2: http://www.w3.org/TR/xmlschema-2/
-.. _boolean: http://www.w3.org/TR/xmlschema-2/#boolean
+.. _XML Schema Part 2: https://www.w3.org/TR/xmlschema-2/
+.. _boolean: https://www.w3.org/TR/xmlschema-2/#boolean
 .. _developers mailing list: https://groups.google.com/forum/?fromgroups=#!forum/openmc-dev

docs/source/devguide/workflow.rst

@@ -91,6 +91,30 @@ features and bug fixes. The general steps for contributing are as follows:
 6. After the pull request has been thoroughly vetted, it is merged back into the
    *develop* branch of openmc-dev/openmc.
 
+Setting Up Upstream Tracking (Required for Versioning)
+------------------------------------------------------
+
+By default, your fork **does not** include tags from the upstream OpenMC repository.  
+OpenMC relies on `git describe --tags` for versioning in source builds, and missing tags can lead  
+to incorrect version detection (i.e., ``0.0.0``). To ensure proper versioning, follow these steps:
+
+1. **Add the Upstream Repository**  
+   This allows you to fetch updates from the main OpenMC repository.
+
+   .. code-block:: sh
+
+       git remote add upstream https://github.com/openmc-dev/openmc.git
+
+2. **Fetch and Push Tags**  
+   Retrieve tags from the upstream repository and update your fork:
+
+   .. code-block:: sh
+
+       git fetch --tags upstream
+       git push --tags origin
+
+This ensures that both your **local** and **remote** fork have the correct versioning information.
+
 Private Development
 -------------------
 
@@ -126,10 +150,10 @@ reinstalling it). While the same effect can be achieved using the
 :envvar:`PYTHONPATH` environment variable, this is generally discouraged as it
 can interfere with virtual environments.
 
-.. _git: http://git-scm.com/
+.. _git: https://git-scm.com/
 .. _GitHub: https://github.com/
 .. _git flow: https://nvie.com/git-model
-.. _valgrind: https://www.valgrind.org/
+.. _valgrind: https://valgrind.org/
 .. _style guide: https://docs.openmc.org/en/latest/devguide/styleguide.html
 .. _pull request: https://docs.github.com/en/github/collaborating-with-issues-and-pull-requests/about-pull-requests
 .. _openmc-dev/openmc: https://github.com/openmc-dev/openmc

docs/source/index.rst

@@ -12,7 +12,7 @@ files produced by NJOY. Parallelism is enabled via a hybrid MPI and OpenMP
 programming model.
 
 OpenMC was originally developed by members of the `Computational Reactor Physics
-Group <http://crpg.mit.edu>`_ at the `Massachusetts Institute of Technology
+Group <https://crpg.mit.edu>`_ at the `Massachusetts Institute of Technology
 <https://web.mit.edu>`_ starting in 2011. Various universities, laboratories,
 and other organizations now contribute to the development of OpenMC. For more
 information on OpenMC, feel free to post a message on the `OpenMC Discourse

docs/source/io_formats/collision_track.rst

@@ -0,0 +1,46 @@
+.. _io_collision_track:
+
+===========================
+Collision Track File Format
+===========================
+
+When collision tracking is enabled with ``mcpl=false`` (the default), OpenMC
+writes binary data to an HDF5 file named ``collision_track.h5``. The same data
+may also be written after each batch when multiple files are requested
+(``collision_track.N.h5``) or when the run is performed in parallel. The file
+contains the information needed to reconstruct each recorded collision.
+
+The current revision of the collision track file format is 1.2.
+
+**/**
+
+:Attributes:
+            - **filetype** (*char[]*) -- String indicating the type of file.
+              For collision-track files the value is ``"collision_track"``.
+
+:Datasets:
+
+          - **collision_track_bank** (Compound type) -- Collision information
+            for each stored event. Each entry in the dataset corresponds to one
+            collision and contains the following fields:
+
+              - ``r`` (*double[3]*) -- Position of the collision in [cm].
+              - ``u`` (*double[3]*) -- Direction unit vector immediately after the collision.
+              - ``E`` (*double*) -- Incident particle energy before the collision in [eV].
+              - ``dE`` (*double*) -- Energy loss over the collision (:math:`E_\text{before} - E_\text{after}`) in [eV].
+              - ``time`` (*double*) -- Time of the collision in [s].
+              - ``wgt`` (*double*) -- Particle weight at the collision.
+              - ``event_mt`` (*int*) -- ENDF MT number identifying the reaction.
+              - ``delayed_group`` (*int*) -- Delayed neutron group index (non-zero for delayed events).
+              - ``cell_id`` (*int*) -- ID of the cell in which the collision occurred.
+              - ``nuclide_id`` (*int*) -- PDG number of the nuclide (100ZZZAAAM).
+              - ``material_id`` (*int*) -- ID of the material containing the collision site.
+              - ``universe_id`` (*int*) -- ID of the universe containing the collision site.
+              - ``n_collision`` (*int*) -- Collision counter for the particle history.
+              - ``particle`` (*int32_t*) -- Particle type (PDG number).
+              - ``parent_id`` (*int64_t*) -- Unique ID of the parent particle.
+              - ``progeny_id`` (*int64_t*) -- Progeny ID of the particle.
+
+In an MPI run, OpenMC writes the combined dataset by gathering collision-track
+entries from all ranks before flushing them to disk, so the final file appears
+as though it were produced serially.

docs/source/io_formats/depletion_chain.rst

@@ -56,6 +56,27 @@ attributes:
 
 .. _io_chain_reaction:
 
+--------------------
+``<source>`` Element
+--------------------
+
+The ``<source>`` element represents photon and electron sources associated with
+the decay of a nuclide and contains information to construct an
+:class:`openmc.stats.Univariate` object that represents this emission as an
+energy distribution. This element has the following attributes:
+
+  :type:
+    The type of :class:`openmc.stats.Univariate` source term.
+
+  :particle:
+    The type of particle emitted, e.g., 'photon' or 'electron'
+
+  :parameters:
+    The parameters of the source term, e.g., for a
+    :class:`openmc.stats.Discrete` source, the energies (in [eV]) at which the
+    particles are emitted and their relative intensities in [Bq/atom] (in other
+    words, decay constants).
+
 ----------------------
 ``<reaction>`` Element
 ----------------------

docs/source/io_formats/depletion_results.rst

@@ -4,7 +4,7 @@
 Depletion Results File Format
 =============================
 
-The current version of the depletion results file format is 1.1.
+The current version of the depletion results file format is 1.3.
 
 **/**
 
@@ -12,30 +12,31 @@ The current version of the depletion results file format is 1.1.
              - **version** (*int[2]*) -- Major and minor version of the
                statepoint file format.
 
-:Datasets: - **eigenvalues** (*double[][][2]*) -- k-eigenvalues at each
-             time/stage. This array has shape (number of timesteps, number of
-             stages, value). The last axis contains the eigenvalue and the
-             associated uncertainty
-           - **number** (*double[][][][]*) -- Total number of atoms. This array
-             has shape (number of timesteps, number of stages, number of
+:Datasets: - **eigenvalues** (*double[][2]*) -- k-eigenvalues at each timestep.
+             This array has shape (number of timesteps, 2). The second axis
+             contains the eigenvalue and its associated uncertainty.
+           - **number** (*double[][][]*) -- Total number of atoms at each
+             timestep. This array has shape (number of timesteps, number of
              materials, number of nuclides).
-           - **reaction rates** (*double[][][][][]*) -- Reaction rates used to
-             build depletion matrices. This array has shape (number of
-             timesteps, number of stages, number of materials, number of
-             nuclides, number of reactions).
+           - **reaction rates** (*double[][][][]*) -- Reaction rates at each
+             timestep. This array has shape (number of timesteps, number of
+             materials, number of nuclides, number of reactions). Only stored if
+             write_rates=True.
            - **time** (*double[][2]*) -- Time in [s] at beginning/end of each
              step.
-           - **source_rate** (*double[][]*) -- Power in [W] or source rate in
-             [neutron/sec]. This array has shape (number of timesteps, number
-             of stages).
+           - **source_rate** (*double[]*) -- Power in [W] or source rate in
+             [neutron/sec] for each timestep.
            - **depletion time** (*double[]*) -- Average process time in [s]
              spent depleting a material across all burnable materials and,
              if applicable, MPI processes.
+           - **keff_search_root** (*double[]*) -- Root of the keff search at the
+             end of the timestep, if applicable.
 
 **/materials/<id>/**
 
 :Attributes: - **index** (*int*) -- Index used in results for this material
              - **volume** (*double*) -- Volume of this material in [cm^3]
+             - **name** (*char[]*) -- Name of this material
 
 **/nuclides/<name>/**
 

docs/source/io_formats/geometry.rst

@@ -38,11 +38,9 @@ Each ``<surface>`` element can have the following attributes or sub-elements:
 
   :boundary:
      The boundary condition for the surface. This can be "transmission",
-     "vacuum", "reflective", or "periodic". Periodic boundary conditions can
-     only be applied to x-, y-, and z-planes. Only axis-aligned periodicity is
-     supported, i.e., x-planes can only be paired with x-planes. Specify which
-     planes are periodic and the code will automatically identify which planes
-     are paired together.
+     "vacuum", "reflective", or "periodic". Specify which planes are
+     periodic and the code will automatically identify which planes are
+     paired together.
 
     *Default*: "transmission"
 
@@ -318,9 +316,10 @@ the following attributes or sub-elements:
     *Default*: None
 
   :orientation:
-    The orientation of the hexagonal lattice. The string "x" indicates that two
-    sides of the lattice are parallel to the x-axis, whereas the string "y"
-    indicates that two sides are parallel to the y-axis.
+    The orientation of the hexagonal lattice. The string "x" indicates that each
+    lattice element has two faces that are perpendicular to the x-axis, whereas
+    the string "y" indicates that each lattice element has two faces that are
+    perpendicular to the y-axis.
 
     *Default*: "y"
 
@@ -407,13 +406,64 @@ Each ``<dagmc_universe>`` element can have the following attributes or sub-eleme
 
     *Default*: None
 
+  :cell:
+    Zero or more ``<cell>`` sub-elements may appear to override properties of
+    individual DAGMC volumes. Each ``<cell>`` element supports the following
+    attributes and sub-elements:
 
-  .. note:: A geometry.xml file containing only a DAGMC model for a file named `dagmc.h5m` (no CSG)
-            looks as follows
+    :id:
+      The integer cell ID in the DAGMC geometry to override. Required.
 
-            .. code-block:: xml
+    :name:
+      An optional string label for the cell.
 
-              <?xml version='1.0' encoding='utf-8'?>
-              <geometry>
-                <dagmc_universe filename="dagmc.h5m" id="1" />
-              </geometry>
+      *Default*: None
+
+    :material:
+      The material ID to assign to this cell. Use ``void`` for vacuum. Multiple
+      space-separated IDs may be given to specify a distribmat (distributed
+      material) assignment. Required.
+
+    :temperature:
+      Temperature(s) in [K] to assign to the cell. Must be greater than or equal
+      to 0. Multiple space-separated values may be given.
+
+      *Default*: None
+
+    :density:
+      Density in [g/cm³] to assign to the cell. Must be greater than 0. Requires a non-void
+      material fill. Multiple space-separated values may be given.
+
+      *Default*: None
+
+    :volume:
+      Volume of the cell in [cm³].
+
+      .. note:: DAGMC can compute cell volumes exactly from the triangulated
+                mesh surfaces. Specifying a manual volume risks inconsistency
+                with that capability.
+
+      *Default*: None
+
+    The following standard ``<cell>`` attributes are **not** supported inside
+    ``<dagmc_universe>`` and will raise an error if present: ``region``,
+    ``fill``, ``universe``, ``translation``, ``rotation``.
+
+    .. deprecated::
+      The ``<material_overrides>`` sub-element (containing ``<cell_override>``
+      children with ``<material_ids>``) is deprecated. A deprecation warning is
+      emitted and the overrides are converted to the ``<cell>`` format at parse
+      time. It is an error to specify both ``<material_overrides>`` and
+      ``<cell>`` sub-elements on the same ``<dagmc_universe>``.
+
+    *Default*: None
+
+.. note:: A geometry.xml file containing only a DAGMC model for a file named
+          `dagmc.h5m` (no CSG) looks as follows:
+
+    .. code-block:: xml
+
+        <?xml version='1.0' encoding='utf-8'?>
+        <geometry>
+          <dagmc_universe filename="dagmc.h5m" id="1" />
+        </geometry>

docs/source/io_formats/index.rst

@@ -44,6 +44,7 @@ Output Files
 
    statepoint
    source
+   collision_track
    summary
    properties
    depletion_results

docs/source/io_formats/mgxs_library.rst

@@ -133,6 +133,10 @@ Temperature-dependent data, provided for temperature <TTT>K.
              This dataset is optional. This is a 1-D vector if `representation`
              is "isotropic", or a 3-D vector if `representation` is "angle"
              with dimensions of [polar][azimuthal][groups].
+	     When this data is not available, an approximation using the
+	     group energy boundaries is used. For more information see
+	     the particle speed subsection in the multigroup-data section
+	     of the theory manual.
 
 **/<library name>/<TTT>K/scatter_data/**
 

docs/source/io_formats/particle_restart.rst

@@ -4,7 +4,7 @@
 Particle Restart File Format
 ============================
 
-The current version of the particle restart file format is 2.0.
+The current version of the particle restart file format is 2.1.
 
 **/**
 
@@ -26,8 +26,7 @@ The current version of the particle restart file format is 2.0.
            - **run_mode** (*char[]*) -- Run mode used, either 'fixed source',
              'eigenvalue', or 'particle restart'.
            - **id** (*int8_t*) -- Unique identifier of the particle.
-           - **type** (*int*) -- Particle type (0=neutron, 1=photon, 2=electron,
-             3=positron)
+           - **type** (*int32_t*) -- Particle type (PDG number)
            - **weight** (*double*) -- Weight of the particle.
            - **energy** (*double*) -- Energy of the particle in eV for
              continuous-energy mode, or the energy group of the particle for

docs/source/io_formats/plots.rst

@@ -7,13 +7,18 @@ Geometry Plotting Specification -- plots.xml
 Basic plotting capabilities are available in OpenMC by creating a plots.xml file
 and subsequently running with the ``--plot`` command-line flag. The root element
 of the plots.xml is simply ``<plots>`` and any number output plots can be
-defined with ``<plot>`` sub-elements.  Two plot types are currently implemented
+defined with ``<plot>`` sub-elements.  Four plot types are currently implemented
 in openMC:
 
 * ``slice``  2D pixel plot along one of the major axes. Produces a PNG image
   file.
 * ``voxel``  3D voxel data dump. Produces an HDF5 file containing voxel xyz
   position and cell or material id.
+* ``wireframe_raytrace``  2D pixel plot of a three-dimensional view of a
+  geometry using wireframes around cells or materials and coloring by depth
+  through each material.
+* ``solid_raytrace``  2D pixel plot of a three-dimensional view of a geometry
+  with solid colored surfaces of a set of cells or materials.
 
 
 ------------------
@@ -66,21 +71,22 @@ sub-elements:
     *Default*: None - Required entry
 
   :type:
-    Keyword for type of plot to be produced. Currently only "slice" and "voxel"
-    plots are implemented. The "slice" plot type creates 2D pixel maps saved in
-    the PNG file format. The "voxel" plot type produces a binary datafile
-    containing voxel grid positioning and the cell or material (specified by the
-    ``color`` tag) at the center of each voxel. Voxel plot files can be
-    processed into VTK files using the :ref:`scripts_voxel` script provided with
-    OpenMC and subsequently viewed with a 3D viewer such as VISIT or Paraview.
-    See the :ref:`io_voxel` for information about the datafile structure.
+    Keyword for type of plot to be produced. Currently "slice", "voxel",
+    "wireframe_raytrace", and "solid_raytrace" plots are implemented. The
+    "slice" plot type creates 2D pixel maps saved in the PNG file format. The
+    "voxel" plot type produces a binary datafile containing voxel grid
+    positioning and the cell or material (specified by the ``color`` tag) at the
+    center of each voxel. Voxel plot files can be processed into VTK files using
+    the :func:`openmc.voxel_to_vtk` function and subsequently viewed with a 3D
+    viewer such as VISIT or Paraview. See :ref:`io_voxel` for information about
+    the datafile structure.
 
     .. note:: High-resolution voxel files produced by OpenMC can be quite large,
               but the equivalent VTK files will be significantly smaller.
 
     *Default*: "slice"
 
-``<plot>`` elements of ``type`` "slice" and "voxel" must contain the ``pixels``
+All ``<plot>`` elements must contain the ``pixels``
 attribute or sub-element:
 
   :pixels:
@@ -96,7 +102,7 @@ attribute or sub-element:
                  ``width``/``pixels`` along that basis direction may not appear
                  in the plot.
 
-    *Default*: None - Required entry for "slice" and "voxel" plots
+    *Default*: None - Required entry for all plots
 
 ``<plot>`` elements of ``type`` "slice" can also contain the following
 attributes or sub-elements.  These are not used in "voxel" plots:
@@ -125,6 +131,11 @@ attributes or sub-elements.  These are not used in "voxel" plots:
       Specifies the custom color for the cell or material. Should be 3 integers
       separated by spaces.
 
+    :xs:
+      The attenuation coefficient for volume rendering of color in units of
+      inverse centimeters. Zero corresponds to transparency. Only for plot type
+      "wireframe_raytrace".
+
     As an example, if your plot is colored by material and you want material 23
     to be blue, the corresponding ``color`` element would look like:
 
@@ -151,6 +162,18 @@ attributes or sub-elements.  These are not used in "voxel" plots:
 
     *Default*: 255 255 255 (white)
 
+  :show_overlaps:
+    Indicates whether overlapping regions of different cells are shown.
+
+    *Default*: None
+
+  :overlap_color:
+    Specifies the RGB color of overlapping regions of different cells. Does not
+    do anything if ``show_overlaps`` is "false" or not specified. Should be 3
+    integers separated by spaces.
+
+    *Default*: 255 0 0 (red)
+
   :meshlines:
     The ``meshlines`` sub-element allows for plotting the boundaries of a
     regular mesh on top of a plot. Only one ``meshlines`` element is allowed per
@@ -179,3 +202,80 @@ attributes or sub-elements.  These are not used in "voxel" plots:
       *Default*: 0 0 0 (black)
 
     *Default*: None
+
+``<plot>`` elements of ``type`` "wireframe_raytrace" or "solid_raytrace" can contain the
+following attributes or sub-elements.
+
+  :camera_position:
+    Location in 3D Cartesian space the camera is at.
+
+
+    *Default*: None - Required for all ``wireframe_raytrace`` or
+    ``solid_raytrace`` plots
+
+  :look_at:
+    Location in 3D Cartesian space the camera is looking at.
+
+
+    *Default*: None - Required for all ``wireframe_raytrace`` or
+    ``solid_raytrace`` plots
+
+  :field_of_view:
+    The horizontal field of view in degrees. Defaults to roughly the same value
+    as for the human eye.
+
+    *Default*: 70
+
+  :orthographic_width:
+    If set to a nonzero value, an orthographic rather than perspective
+    projection for the camera is employed. An orthographic projection puts out
+    parallel rays from the camera of a width prescribed here in the horizontal
+    direction, with the width in the vertical direction decided by the pixel
+    aspect ratio.
+
+    *Default*: 0
+
+``<plot>`` elements of ``type`` "solid_raytrace" can contain the following attributes or
+sub-elements.
+
+  :opaque_ids:
+    List of integer IDs of cells or materials to be treated as visible in the
+    plot. Whether the integers are interpreted as cell or material IDs depends
+    on ``color_by``.
+
+    *Default*: None - Required for all phong plots
+
+  :light_position:
+    Location in 3D Cartesian space of the light.
+
+
+    *Default*: Same location as ``camera_position``
+
+  :diffuse_fraction:
+    Fraction of light originating from non-directional sources. If set to one,
+    the coloring is not influenced by surface curvature, and no shadows appear.
+    If set to zero, only regions illuminated by the light are not black.
+
+
+    *Default*: 0.1
+
+``<plot>`` elements of ``type`` "wireframe_raytrace" can contain the following
+attributes or sub-elements.
+
+  :wireframe_color:
+    RGB value of the wireframe's color
+
+    *Default*: 0, 0, 0 (black)
+
+  :wireframe_thickness:
+    Integer number of pixels that the wireframe takes up. The value is a radius
+    of the wireframe. Setting to zero removes any wireframing.
+
+    *Default*: 0
+
+  :wireframe_ids:
+    Integer IDs of cells or materials of regions to draw wireframes around.
+    Whether the integers are interpreted as cell or material IDs depends on
+    ``color_by``.
+
+    *Default*: None

docs/source/io_formats/properties.rst

@@ -4,7 +4,7 @@
 Properties File Format
 ======================
 
-The current version of the properties file format is 1.0.
+The current version of the properties file format is 1.1.
 
 **/**
 
@@ -25,6 +25,7 @@ The current version of the properties file format is 1.0.
 **/geometry/cells/cell <uid>/**
 
 :Datasets: - **temperature** (*double[]*) -- Temperature of the cell in [K].
+           - **density** (*double[]*) -- Density of the cell in [g/cm3].
 
 **/materials/**
 

docs/source/io_formats/settings.rst

@@ -7,6 +7,19 @@ Settings Specification -- settings.xml
 All simulation parameters and miscellaneous options are specified in the
 settings.xml file.
 
+-------------------------------
+``<atomic_relaxation>`` Element
+-------------------------------
+
+The ``<atomic_relaxation>`` element determines whether the atomic relaxation
+cascade, the X-ray fluorescence photons and Auger electrons emitted when an
+inner-shell vacancy is filled, is simulated following photoelectric and
+incoherent (Compton) scattering interactions. Disabling this can speed up
+photon transport calculations where the detailed secondary particle cascade is
+not of interest.
+
+  *Default*: true
+
 ---------------------
 ``<batches>`` Element
 ---------------------
@@ -20,6 +33,90 @@ source neutrons.
 
   *Default*: None
 
+-----------------------------
+``<collision_track>`` Element
+-----------------------------
+
+The ``<collision_track>`` element indicates to track information about particle
+collisions based on a set of criteria and store these events in a file named
+``collision_track.h5``. This file records details such as the position of the
+interaction, direction of the incoming particle, incident energy and deposited
+energy, weight, time of the interaction, and the delayed neutron group (0 for
+prompt neutrons). Additional information such as the cell ID, material ID,
+universe ID, nuclide ZAID, particle type, and event MT number are also stored.
+Users can specify one or more criterion to filter collisions. If no criteria are
+specified, it defaults to tracking all collisions across the model.
+
+.. warning::
+    Storing all collisions can be very memory intensive. For more targeted
+    tracking, users can employ a variety of parameters such as ``cell_ids``,
+    ``reactions``, ``universe_ids``, ``material_ids``, ``nuclides``, and
+    ``deposited_E_threshold`` to refine the selection of particle interactions
+    to be banked.
+
+This element can contain one or more of the following attributes or
+sub-elements:
+
+  :max_collisions:
+    An integer indicating the maximum number of collisions to be banked per file.
+
+    *Default*: 1000
+
+  :max_collision_track_files:
+    An integer indicating the number of collision_track files to be used.
+
+    *Default*: 1
+
+  :mcpl:
+    An optional boolean to enable MCPL_-format instead of the native HDF5-based
+    format. If activated, the output file name and type is changed to
+    ``collision_track.mcpl``.
+
+    *Default*: false
+
+    .. _MCPL: https://mctools.github.io/mcpl/mcpl.pdf
+
+  :cell_ids:
+    A list of integers representing cell IDs to define specific cells in which
+    collisions are to be banked.
+
+    *Default*: None
+
+  :universe_ids:
+    A list of integers representing the universe IDs to define specific
+    universes in which collisions are to be banked.
+
+    *Default*: None
+
+  :material_ids:
+    A list of integers representing the material IDs to define specific
+    materials in which collisions are to be banked.
+
+    *Default*: None
+
+  :nuclides:
+    A list of strings representing the nuclide, to define specific
+    define specific target nuclide collisions to be banked.
+
+    .. note::
+      Electron and positron collision-track events are not associated with
+      a specific nuclide. If a ``nuclides`` entry is specified, these events
+      are omitted.
+
+    *Default*: None
+
+  :reactions:
+    A list of integers representing the ENDF-6 format MT numbers or strings
+    (e.g. (n,fission)) to define specific reaction types to be banked.
+
+    *Default*: None
+
+  :deposited_E_threshold:
+    A float defining the minimum deposited energy per collision (in eV) to
+    trigger banking.
+
+    *Default*: 0.0
+
 ----------------------------------
 ``<confidence_intervals>`` Element
 ----------------------------------
@@ -81,6 +178,13 @@ time.
 
     *Default*: 1.0
 
+  :survival_normalization:
+    If this element is set to "true", this will enable the use of survival
+    biasing source normalization, whereby the weight parameters, weight and
+    weight_avg, are multiplied per history by the start weight of said history.
+
+    *Default*: false
+
   :energy_neutron:
     The energy under which neutrons will be killed.
 
@@ -171,6 +275,16 @@ history-based parallelism.
 
   *Default*: false
 
+--------------------------------
+``<free_gas_threshold>`` Element
+--------------------------------
+
+The ``<free_gas_threshold>`` element specifies the energy multiplier, expressed
+in units of :math:`kT`, that determines when the free gas scattering approach is
+used for elastic scattering. Values must be positive.
+
+  *Default*: 400.0
+
 -----------------------------------
 ``<generations_per_batch>`` Element
 -----------------------------------
@@ -181,6 +295,15 @@ ignored for all run modes other than "eigenvalue".
 
   *Default*: 1
 
+------------------------------
+``<ifp_n_generation>`` Element
+------------------------------
+
+The ``<ifp_n_generation>`` element indicates the number of generations to
+consider for the Iterated Fission Probability method.
+
+  *Default*: 10
+
 ----------------------
 ``<inactive>`` Element
 ----------------------
@@ -238,7 +361,7 @@ based on the recommended value in LA-UR-14-24530_.
 
   .. note:: This element is not used in the multi-group :ref:`energy_mode`.
 
-.. _LA-UR-14-24530: https://laws.lanl.gov/vhosts/mcnp.lanl.gov/pdf_files/la-ur-14-24530.pdf
+.. _LA-UR-14-24530: https://mcnp.lanl.gov/pdf_files/TechReport_2014_LANL_LA-UR-14-24530_Brown.pdf
 
 ---------------------------
 ``<material_cell_offsets>``
@@ -252,19 +375,45 @@ to false.
 
   *Default*: true
 
-----------------------------------------
-``<max_particles_in_flight>`` Element
-----------------------------------------
+--------------------------------
+``<max_lost_particles>`` Element
+--------------------------------
 
-This element indicates the number of neutrons to run in flight concurrently
+This element indicates the maximum number of lost particles.
+
+  *Default*: 10
+
+------------------------------------
+``<rel_max_lost_particles>`` Element
+------------------------------------
+
+
+This element indicates the maximum number of lost particles, relative to the
+total number of particles.
+
+  *Default*: 1.0e-6
+
+-------------------------------------
+``<max_particles_in_flight>`` Element
+-------------------------------------
+
+This element indicates the number of particles to run in flight concurrently
 when using event-based parallelism. A higher value uses more memory, but
 may be more efficient computationally.
 
   *Default*: 100000
 
----------------------------
+---------------------------------
+``<max_particle_events>`` Element
+---------------------------------
+
+This element indicates the maximum number of events a particle can undergo.
+
+  *Default*: 1000000
+
+-----------------------
 ``<max_order>`` Element
----------------------------
+-----------------------
 
 The ``<max_order>`` element allows the user to set a maximum scattering order
 to apply to every nuclide/material in the problem.  That is, if the data
@@ -276,6 +425,31 @@ then, OpenMC will only use up to the :math:`P_1` data.
   .. note:: This element is not used in the continuous-energy
     :ref:`energy_mode`.
 
+--------------------------------
+``<max_history_splits>`` Element
+--------------------------------
+
+The ``<max_history_splits>`` element indicates the number of times a particle
+can split during a history.
+
+  *Default*: 1000
+
+-----------------------------
+``<max_secondaries>`` Element
+-----------------------------
+
+The ``<max_secondaries>`` element indicates the maximum secondary bank size.
+
+  *Default*: 10000
+
+------------------------
+``<max_tracks>`` Element
+------------------------
+
+The ``<max_tracks>`` element indicates the maximum number of tracks written to a
+track file (per MPI process).
+
+  *Default*: 1000
 
 --------------------------------------
 ``<max_write_lost_particles>`` Element
@@ -386,6 +560,18 @@ generator during generation of colors in plots.
 
   *Default*: 1
 
+.. _properties_file:
+
+-----------------------------
+``<properties_file>`` Element
+-----------------------------
+
+  The ``properties_file`` element has no attributes and contains the path to a
+  properties HDF5 file to load cell temperatures/densities and material
+  densities.
+
+  *Default*: None
+
 ---------------------
 ``<ptables>`` Element
 ---------------------
@@ -398,6 +584,93 @@ or sub-elements and can be set to either "false" or "true".
 
   .. note:: This element is not used in the multi-group :ref:`energy_mode`.
 
+------------------------
+``<random_ray>`` Element
+------------------------
+
+The ``<random_ray>`` element enables random ray mode and contains a number of
+settings relevant to the solver. Tips for selecting these parameters can be
+found in the :ref:`random ray user guide <random_ray>`.
+
+  :distance_inactive:
+    The inactive ray length (dead zone length) in [cm].
+
+    *Default*: None
+
+  :distance_active:
+    The active ray length in [cm].
+
+    *Default*: None
+
+  :ray_source:
+    Specifies the starting ray distribution, and follows the format for
+    :ref:`source_element`. It must be uniform in space and angle and cover the
+    full domain. It does not represent a physical neutron or photon source -- it
+    is only used to sample integrating ray starting locations and directions.
+
+    *Default*: None
+
+  :adjoint_source:
+    Specifies an adjoint fixed source for adjoint transport simulations, and
+    follows the format for :ref:`source_element`. The distributions which make
+    up the adjoint source are subject to the same restrictions as forward
+    fixed sources in Random Ray mode.
+
+    *Default*: None
+
+  :adjoint:
+    Specifies whether to perform adjoint transport. The default is 'False',
+    corresponding to forward transport.
+
+    *Default*: None
+
+  :volume_estimator:
+    Specifies choice of volume estimator for the random ray solver. Options
+    are 'naive', 'simulation_averaged', or 'hybrid'. The default is 'hybrid'.
+
+    *Default*: None
+
+  :volume_normalized_flux_tallies:
+    Specifies whether to normalize flux tallies by volume (bool). The
+    default is 'False'. When enabled, flux tallies will be reported in units
+    of cm/cm^3. When disabled, flux tallies will be reported in units of cm
+    (i.e., total distance traveled by neutrons in the spatial tally
+    region).
+
+    *Default*: None
+
+  :sample_method:
+    Specifies the method for sampling the starting ray distribution. This
+    element can be set to "prng" or "halton".
+
+    *Default*: prng
+
+  :source_region_meshes:
+    Relates meshes to spatial domains for subdividing source regions with each domain.
+
+    :mesh:
+      Contains an ``id`` attribute and one or more ``<domain>`` sub-elements.
+
+      :id:
+        The unique identifier for the mesh.
+
+      :domain:
+        Each domain element has an ``id`` attribute and a ``type`` attribute.
+
+        :id:
+          The unique identifier for the domain.
+
+        :type:
+          The type of the domain. Can be ``material``, ``cell``, or ``universe``.
+
+  :diagonal_stabilization_rho:
+    The rho factor for use with diagonal stabilization. This technique is
+    applied when negative diagonal (in-group) elements are detected in
+    the scattering matrix of input MGXS data, which is a common feature
+    of transport corrected MGXS data.
+
+    *Default*: 1.0
+
 ----------------------------------
 ``<resonance_scattering>`` Element
 ----------------------------------
@@ -473,6 +746,19 @@ pseudo-random number generator.
 
   *Default*: 1
 
+-----------------------------------
+``<shared_secondary_bank>`` Element
+-----------------------------------
+
+  The ``shared_secondary_bank`` element indicates whether to use a shared
+  secondary particle bank. When enabled, secondary particles are collected into
+  a global bank, sorted for reproducibility, and load-balanced across MPI ranks
+  between generations. If not specified, the shared secondary bank is enabled
+  automatically for fixed-source simulations with weight windows active, and
+  disabled otherwise.
+
+.. _source_element:
+
 --------------------
 ``<source>`` Element
 --------------------
@@ -491,10 +777,15 @@ attributes/sub-elements:
     *Default*: 1.0
 
   :type:
-    Indicator of source type. One of ``independent``, ``file``, or ``compiled``.
+    Indicator of source type. One of ``independent``, ``file``, ``compiled``,
+    ``mesh``, or ``tokamak``. The type of the source will be determined by this
+    attribute if it is present.
 
   :particle:
-    The source particle type, either ``neutron`` or ``photon``.
+    The source particle type, specified as a PDG number or a string alias (e.g.,
+    ``neutron``/``n``, ``photon``/``gamma``, ``electron``, ``positron``,
+    ``proton``/``p``, ``deuteron``/``d``, ``triton``/``t``, ``alpha``, or GNDS
+    nuclide names like ``Fe57``).
 
     *Default*: neutron
 
@@ -534,24 +825,38 @@ attributes/sub-elements:
 
     :type:
       The type of spatial distribution. Valid options are "box", "fission",
-      "point", "cartesian", "cylindrical", and "spherical". A "box" spatial
-      distribution has coordinates sampled uniformly in a parallelepiped. A
-      "fission" spatial distribution samples locations from a "box"
+      "point", "cartesian", "cylindrical", "spherical", "mesh", and "cloud".
+
+      A "box" spatial distribution has coordinates sampled uniformly in a
+      parallelepiped.
+
+      A "fission" spatial distribution samples locations from a "box"
       distribution but only locations in fissionable materials are accepted.
+
       A "point" spatial distribution has coordinates specified by a triplet.
+
       A "cartesian" spatial distribution specifies independent distributions of
-      x-, y-, and z-coordinates. A "cylindrical" spatial distribution specifies
-      independent distributions of r-, phi-, and z-coordinates where phi is the
-      azimuthal angle and the origin for the cylindrical coordinate system is
-      specified by origin. A "spherical" spatial distribution specifies
-      independent distributions of r-, cos_theta-, and phi-coordinates where
-      cos_theta is the cosine of the angle with respect to the z-axis, phi is
-      the azimuthal angle, and the sphere is centered on the coordinate
-      (x0,y0,z0). A "mesh" spatial distribution samples source sites from a mesh element
+      x-, y-, and z-coordinates.
+
+      A "cylindrical" spatial distribution specifies independent distributions
+      of r-, phi-, and z-coordinates where phi is the azimuthal angle and the
+      origin for the cylindrical coordinate system is specified by origin.
+
+      A "spherical" spatial distribution specifies independent distributions of
+      r-, cos_theta-, and phi-coordinates where cos_theta is the cosine of the
+      angle with respect to the z-axis, phi is the azimuthal angle, and the
+      sphere is centered on the coordinate (x0,y0,z0).
+
+      A "mesh" spatial distribution samples source sites from a mesh element
       based on the relative strengths provided in the node. Source locations
       within an element are sampled isotropically. If no strengths are provided,
       the space within the mesh is uniformly sampled.
 
+      A "cloud" spatial distribution samples source sites from a list of spatial
+      positions provided in the node, based on the relative strengths provided
+      in the node. If no strengths are provided, the positions are uniformly
+      sampled.
+
       *Default*: None
 
     :parameters:
@@ -570,6 +875,7 @@ attributes/sub-elements:
 
       For a "cylindrical" distribution, no parameters are specified. Instead,
       the ``r``, ``phi``, ``z``, and ``origin`` elements must be specified.
+      Optionally, the ``r_dir`` and ``z_dir`` elements could be specified.
 
       For a "spherical" distribution, no parameters are specified. Instead,
       the ``r``, ``theta``, ``phi``, and ``origin`` elements must be specified.
@@ -601,6 +907,10 @@ attributes/sub-elements:
       of a univariate probability distribution (see the description in
       :ref:`univariate`).
 
+    :r_dir:
+      For "cylindrical" distributions, this element specifies the direction
+      of the cylinder r-axis at phi=0. Defaults to (1.0, 0.0, 0.0).
+
     :theta:
       For a "spherical" distribution, this element specifies the distribution
       of theta-coordinates. The necessary sub-elements/attributes are those of a
@@ -613,10 +923,39 @@ attributes/sub-elements:
       sub-elements/attributes are those of a univariate probability
       distribution (see the description in :ref:`univariate`).
 
+    :z_dir:
+      For "cylindrical" distributions, this element specifies the direction
+      of the cylinder z-axis. Defaults to (0.0, 0.0, 1.0).
+
     :origin:
       For "cylindrical and "spherical" distributions, this element specifies
       the coordinates for the origin of the coordinate system.
 
+    :mesh_id:
+      For "mesh" spatial distributions, this element specifies which mesh ID to
+      use for the geometric description of the mesh.
+
+    :coords:
+      For "cloud" distributions, this element specifies a list of coordinates
+      for each of the points in the cloud.
+
+    :strengths:
+      For "mesh" and "cloud" spatial distributions, this element specifies the
+      relative source strength of each mesh element or each point in the cloud.
+
+    :volume_normalized:
+      For "mesh" spatial distributions, this optional boolean element specifies
+      whether the vector of relative strengths should be multiplied by the mesh
+      element volume. This is most common if the strengths represent a source
+      per unit volume.
+
+      *Default*: false
+
+    :bias:
+      For "mesh" and "cloud" spatial distributions, this optional element
+      specifies floating point values corresponding to alternative probabilities
+      for each value/component to use for biased sampling.
+
   :angle:
     An element specifying the angular distribution of source sites. This element
     has the following attributes:
@@ -649,6 +988,10 @@ attributes/sub-elements:
       are those of a univariate probability distribution (see the description in
       :ref:`univariate`).
 
+    :bias:
+      For "isotropic" angular distributions, this optional element specifies a
+      "mu-phi" angular distribution used for biased sampling.
+
   :energy:
     An element specifying the energy distribution of source sites. The necessary
     sub-elements/attributes are those of a univariate probability distribution
@@ -664,6 +1007,125 @@ attributes/sub-elements:
 
     *Default*: false
 
+  :mesh:
+    For mesh sources, this indicates the ID of the corresponding mesh.
+
+  :source:
+    For mesh sources, this sub-element specifies the source for an individual
+    mesh element and follows the format for :ref:`source_element`. The number of
+    ``<source>`` sub-elements should correspond to the number of mesh elements.
+
+  For a source with ``type="tokamak"``, the spatial distribution is described by
+  a Miller-style flux-surface parameterization and the following sub-elements
+  are used instead of the ``space`` element:
+
+  :major_radius:
+    The major radius :math:`R_0` of the plasma in [cm].
+
+  :minor_radius:
+    The minor radius :math:`a` of the plasma in [cm]. Must be smaller than
+    ``major_radius``.
+
+  :elongation:
+    The plasma elongation :math:`\kappa` (must be > 0).
+
+  :triangularity:
+    The plasma triangularity :math:`\delta` (must be in [-1, 1]). Negative
+    values describe negative-triangularity plasmas.
+
+  :shafranov_shift:
+    The Shafranov shift :math:`\Delta` in [cm] (must be >= 0 and less than
+    ``minor_radius``/2).
+
+  :r_over_a:
+    A list of normalized minor-radius grid points :math:`r/a`. Must be strictly
+    increasing, start at 0, and end at 1.
+
+  :emission_density:
+    A list of neutron emission densities :math:`S(r)` evaluated at each
+    ``r_over_a`` grid point (arbitrary units, must be non-negative). Only the
+    shape matters, since the profile is normalized internally. Values are
+    interpolated linearly between grid points and the profile is refined on an
+    internal grid for radial sampling. Must have the same length as
+    ``r_over_a`` and contain at least one positive value.
+
+  :phi_start:
+    The starting toroidal angle in [rad].
+
+    *Default*: 0.0
+
+  :phi_extent:
+    The toroidal angle extent in [rad]. The source is sampled uniformly in
+    :math:`[\phi_\text{start},\ \phi_\text{start} + \phi_\text{extent}]`.
+
+    *Default*: :math:`2\pi`
+
+  :n_alpha:
+    The number of poloidal-angle grid points used to build the sampling CDFs
+    (must be > 2). Larger values reduce discretization bias; values below 51
+    produce a warning.
+
+    *Default*: 101
+
+  :vertical_shift:
+    A vertical shift of the plasma center in [cm].
+
+    *Default*: 0.0
+
+  :energy:
+    For a tokamak source, one or more ``energy`` sub-elements specify the
+    neutron energy distribution(s). Either a single distribution is given (used
+    at all radii) or exactly one distribution per ``r_over_a`` grid point is
+    given, in which case the energy is sampled from one of the two
+    distributions bracketing the sampled radius, selected stochastically with
+    probability proportional to the proximity of the radius to each grid point
+    (stochastic interpolation). Each follows the format of a univariate
+    probability distribution (see :ref:`univariate`).
+
+  :time:
+    An optional ``time`` sub-element specifying the time distribution of source
+    particles, following the format of a univariate probability distribution
+    (see :ref:`univariate`).
+
+    *Default*: particles are born at :math:`t=0`
+
+  .. note:: Biased sampling can be applied to the spatial and energy distributions
+            of a source by using the ``<bias>`` sub-element (see
+            :ref:`univariate` for details on how to specify bias distributions).
+
+  :constraints:
+    This sub-element indicates the presence of constraints on sampled source
+    sites (see :ref:`usersguide_source_constraints` for details). It may have
+    the following sub-elements:
+
+    :domain_ids:
+      The unique IDs of domains for which source sites must be within.
+
+      *Default*: None
+
+    :domain_type:
+      The type of each domain for source rejection ("cell", "material", or
+      "universe").
+
+      *Default*: None
+
+    :fissionable:
+      A boolean indicating whether source sites must be sampled within a
+      material that is fissionable in order to be accepted.
+
+    :time_bounds:
+      A pair of times in [s] indicating the lower and upper bound for a time
+      interval that source particles must be within.
+
+    :energy_bounds:
+      A pair of energies in [eV] indicating the lower and upper bound for an
+      energy interval that source particles must be within.
+
+    :rejection_strategy:
+      Either "resample", indicating that source sites should be resampled when
+      one is rejected, or "kill", indicating that a rejected source site is
+      assigned zero weight.
+
 .. _univariate:
 
 Univariate Probability Distributions
@@ -677,17 +1139,19 @@ variable and whose sub-elements/attributes are as follows:
 
 :type:
   The type of the distribution. Valid options are "uniform", "discrete",
-  "tabular", "maxwell", "watt", and "mixture". The "uniform" option produces
-  variates sampled from a uniform distribution over a finite interval. The
-  "discrete" option produces random variates that can assume a finite number
-  of values (i.e., a distribution characterized by a probability mass function).
-  The "tabular" option produces random variates sampled from a tabulated
-  distribution where the density function is either a histogram or
+  "tabular", "maxwell", "watt", "mixture", and "decay_spectrum". The "uniform"
+  option produces variates sampled from a uniform distribution over a finite
+  interval. The "discrete" option produces random variates that can assume a
+  finite number of values (i.e., a distribution characterized by a probability
+  mass function). The "tabular" option produces random variates sampled from a
+  tabulated distribution where the density function is either a histogram or
   linearly-interpolated between tabulated points. The "watt" option produces
   random variates is sampled from a Watt fission spectrum (only used for
   energies). The "maxwell" option produce variates sampled from a Maxwell
-  fission spectrum (only used for energies). The "mixture" option produces samples
-  from univariate sub-distributions with given probabilities.
+  fission spectrum (only used for energies). The "mixture" option produces
+  samples from univariate sub-distributions with given probabilities. The
+  "decay_spectrum" option produces photon energies sampled from decay photon
+  spectra in a depletion chain (only used for energies).
 
   *Default*: None
 
@@ -705,6 +1169,10 @@ variable and whose sub-elements/attributes are as follows:
   :math:`(x,p)` pairs defining the discrete/tabular distribution. All :math:`x`
   points are given first followed by corresponding :math:`p` points.
 
+  For a "decay_spectrum" distribution, ``parameters`` gives the atom densities
+  in [atom/b-cm] for the nuclides listed in the ``nuclides`` element, in the
+  same order.
+
   For a "watt" distribution, ``parameters`` should be given as two real numbers
   :math:`a` and :math:`b` that parameterize the distribution :math:`p(x) dx = c
   e^{-x/a} \sinh \sqrt{b \, x} dx`.
@@ -723,30 +1191,51 @@ variable and whose sub-elements/attributes are as follows:
   *Default*: histogram
 
 :pair:
-  For a "mixture" distribution, this element provides a distribution and its corresponding probability.
+  For a "mixture" distribution, this element provides a distribution and its
+  corresponding probability.
 
   :probability:
-    An attribute or ``pair`` that provides the probability of a univariate distribution within a "mixture" distribution.
+    An attribute or ``pair`` that provides the probability of a univariate
+    distribution within a "mixture" distribution.
 
   :dist:
-    This sub-element of a ``pair`` element provides information on the corresponding univariate distribution.
+    This sub-element of a ``pair`` element provides information on the
+    corresponding univariate distribution.
 
--------------------------
-``<state_point>`` Element
--------------------------
+:volume:
+  For a "decay_spectrum" distribution, this attribute specifies the source
+  region volume in cm\ :sup:`3`. It is used together with atom densities to
+  determine the absolute photon emission rate. When a source uses a
+  "decay_spectrum" energy distribution, the source strength is set from this
+  emission rate.
 
-The ``<state_point>`` element indicates at what batches a state point file
-should be written. A state point file can be used to restart a run or to get
-tally results at any batch. The default behavior when using this tag is to
-write out the source bank in the state_point file. This behavior can be
-customized by using the ``<source_point>`` element. This element has the
-following attributes/sub-elements:
+:nuclides:
+  For a "decay_spectrum" distribution, this element specifies a
+  whitespace-separated list of nuclide names contributing to the decay photon
+  source. The atom densities for these nuclides are given by the ``parameters``
+  element in the same order. Nuclides are resolved against the depletion chain,
+  and nuclides without decay photon spectra do not contribute to the
+  distribution.
 
-  :batches:
-    A list of integers separated by spaces indicating at what batches a state
-    point file should be written.
+:bias:
+  This optional element specifies a biased distribution for importance sampling.
+  For continuous distributions, the ``bias`` element should contain another
+  univariate distribution with the same support (interval) as the parent
+  distribution. For discrete distributions, the ``bias`` element should contain
+  floating point values corresponding to alternative probabilities for each
+  value/component to be used for biased sampling.
 
-    *Default*: Last batch only
+  *Default*: None
+
+---------------------------------------
+``<source_rejection_fraction>`` Element
+---------------------------------------
+
+The ``<source_rejection_fraction>`` element specifies the minimum fraction of
+external source sites that must be accepted when applying rejection sampling
+based on constraints.
+
+   *Default*: 0.05
 
 --------------------------
 ``<source_point>`` Element
@@ -797,6 +1286,32 @@ attributes/sub-elements:
 
     *Default*: false
 
+-------------------------
+``<state_point>`` Element
+-------------------------
+
+The ``<state_point>`` element indicates at what batches a state point file
+should be written. A state point file can be used to restart a run or to get
+tally results at any batch. The default behavior when using this tag is to
+write out the source bank in the state_point file. This behavior can be
+customized by using the ``<source_point>`` element. This element has the
+following attributes/sub-elements:
+
+  :batches:
+    A list of integers separated by spaces indicating at what batches a state
+    point file should be written.
+
+    *Default*: Last batch only
+
+--------------------
+``<stride>`` Element
+--------------------
+
+The ``stride`` element is used to specify how many random numbers are allocated
+for each source particle history.
+
+  *Default*: 152,917
+
 ------------------------------
 ``<surf_source_read>`` Element
 ------------------------------
@@ -816,7 +1331,12 @@ attributes/sub-elements:
 
 The ``<surf_source_write>`` element triggers OpenMC to bank particles crossing
 certain surfaces and write out the source bank in a separate file called
-``surface_source.h5``. This element has the following attributes/sub-elements:
+``surface_source.h5``. One or multiple surface IDs and one cell ID can be used
+to select the surfaces of interest. If no surface IDs are declared, every surface
+of the model is eligible to bank particles. In that case, a cell ID (using
+either the ``cell``, ``cellfrom`` or ``cellto`` attributes) can be used to select
+every surface of a specific cell. This element has the following
+attributes/sub-elements:
 
   :surface_ids:
     A list of integers separated by spaces indicating the unique IDs of surfaces
@@ -832,6 +1352,15 @@ certain surfaces and write out the source bank in a separate file called
 
     *Default*: None
 
+  :max_source_files:
+    An integer value indicating the number of surface source files to be written
+    containing the maximum number of particles each. The surface source bank
+    will be cleared in simulation memory each time a surface source file is
+    written. By default a ``surface_source.h5`` file will be created when the
+    maximum number of saved particles is reached.
+
+    *Default*: 1
+
   :mcpl:
     An optional boolean which indicates if the banked particles should be
     written to a file in the MCPL_-format instead of the native HDF5-based
@@ -842,6 +1371,51 @@ certain surfaces and write out the source bank in a separate file called
 
     .. _MCPL: https://mctools.github.io/mcpl/mcpl.pdf
 
+  :cell:
+    An integer representing the cell ID used to determine if particles crossing
+    identified surfaces are to be banked. Particles coming from or going to this
+    declared cell will be banked if they cross the identified surfaces.
+
+    *Default*: None
+
+  :cellfrom:
+    An integer representing the cell ID used to determine if particles crossing
+    identified surfaces are to be banked. Particles coming from this declared
+    cell will be banked if they cross the identified surfaces.
+
+    *Default*: None
+
+  :cellto:
+    An integer representing the cell ID used to determine if particles crossing
+    identified surfaces are to be banked. Particles going to this declared cell
+    will be banked if they cross the identified surfaces.
+
+    *Default*: None
+
+.. note:: The ``cell``, ``cellfrom`` and ``cellto`` attributes cannot be
+          used simultaneously.
+
+.. note:: Surfaces with boundary conditions that are not "transmission" or "vacuum"
+          are not eligible to store any particles when using ``cell``, ``cellfrom``
+          or ``cellto`` attributes. It is recommended to use surface IDs instead.
+
+------------------------------------
+``<surface_grazing_cutoff>`` Element
+------------------------------------
+
+The ``<surface_grazing_cutoff>`` element specifies the surface flux cosine cutoff.
+
+  *Default*: 0.001
+
+-----------------------------------
+``<surface_grazing_ratio>`` Element
+-----------------------------------
+
+The ``<surface_grazing_ratio>`` element specifies the surface flux cosine
+substitution ratio.
+
+  *Default*: 0.5
+
 ------------------------------
 ``<survival_biasing>`` Element
 ------------------------------
@@ -1025,6 +1599,15 @@ has the following attributes/sub-elements:
               for fixed source and small criticality calculations, but is very
               optimistic for highly coupled full-core reactor problems.
 
+-------------------------------------
+``<uniform_source_sampling>`` Element
+-------------------------------------
+
+The ``<uniform_source_sampling>`` element indicates whether to sample among
+multiple sources uniformly, applying their strengths as weights to sampled
+particles.
+
+  *Default*: False
 
 ------------------------
 ``<ufs_mesh>`` Element
@@ -1037,6 +1620,16 @@ Agency Monte Carlo Performance Benchmark Problem," Proceedings of *Physor 2012*,
 Knoxville, TN (2012). The mesh should cover all possible fissionable materials
 in the problem and is specified using a :ref:`mesh_element`.
 
+-------------------------------
+``<use_decay_photons>`` Element
+-------------------------------
+
+The ``<use_decay_photons>`` element indicates whether to produce decay photons
+from neutron reactions instead of prompt photons. This is used in conjunction
+with the direct 1-step method for shutdown dose rate calculations.
+
+  *Default*: False
+
 .. _verbosity:
 
 -----------------------
@@ -1138,7 +1731,8 @@ sub-elements/attributes:
     *Default*: None
 
   :particle_type:
-    The particle that the weight windows will apply to (e.g., 'neutron')
+    The particle that the weight windows will apply to, specified as a PDG
+    code or string (e.g., ``neutron``).
 
     *Default*: 'neutron'
 
@@ -1198,7 +1792,8 @@ mesh-based weight windows.
     *Default*: None
 
   :particle_type:
-    The particle that the weight windows will apply to (e.g., 'neutron')
+    The particle that the weight windows will apply to, specified as a PDG
+    code or string (e.g., ``neutron``).
 
     *Default*: neutron
 
@@ -1218,7 +1813,7 @@ mesh-based weight windows.
     *Default*: true
 
   :method:
-    Method used to update weight window values (currently only 'magic' is supported)
+    Method used to update weight window values (one of 'magic' or 'fw_cadis')
 
     *Default*: magic
 
@@ -1242,6 +1837,14 @@ mesh-based weight windows.
 
         *Default*: 5.0
 
+    For FW-CADIS:
+
+      :targets:
+        A sequence of IDs corresponding to the tallies which cover phase
+        space regions of interest for local variance reduction.
+
+        *Default*: None
+
 ---------------------------------------
 ``<weight_window_checkpoints>`` Element
 ---------------------------------------
@@ -1267,3 +1870,21 @@ following sub-elements/attributes:
 
   The ``weight_windows_file`` element has no attributes and contains the path to
   a weight windows HDF5 file to load during simulation initialization.
+
+-------------------------------
+``<weight_windows_on>`` Element
+-------------------------------
+
+  The ``weight_windows_on`` element indicates whether weight windows are
+  enabled.
+
+  *Default*: False
+
+----------------------------------
+``<write_initial_source>`` Element
+----------------------------------
+
+  The ``write_initial_source`` element indicates whether to write the initial
+  source distribution to file.
+
+  *Default*: False

docs/source/io_formats/source.rst

@@ -15,6 +15,8 @@ following the same format.
 **/**
 
 :Attributes: - **filetype** (*char[]*) -- String indicating the type of file.
+             - **version** (*int[2]*) -- Major and minor version of the source
+               file format.
 
 :Datasets:
 
@@ -22,5 +24,5 @@ following the same format.
              particle. The compound type has fields ``r``, ``u``, ``E``,
              ``time``, ``wgt``, ``delayed_group``, ``surf_id`` and ``particle``,
              which represent the position, direction, energy, time, weight,
-             delayed group, surface ID, and particle type (0=neutron, 1=photon,
-             2=electron, 3=positron), respectively.
+             delayed group, surface ID, and particle type (PDG number),
+             respectively.

docs/source/io_formats/statepoint.rst

@@ -4,7 +4,7 @@
 State Point File Format
 =======================
 
-The current version of the statepoint file format is 18.1.
+The current version of the statepoint file format is 18.2.
 
 **/**
 
@@ -23,6 +23,7 @@ The current version of the statepoint file format is 18.1.
                bank is present (1) or not (0).
 
 :Datasets: - **seed** (*int8_t*) -- Pseudo-random number generator seed.
+           - **stride** (*uint64_t*) -- Pseudo-random number generator stride.
            - **energy_mode** (*char[]*) -- Energy mode of the run, either
              'continuous-energy' or 'multi-group'.
            - **run_mode** (*char[]*) -- Run mode used, either 'eigenvalue' or
@@ -55,8 +56,8 @@ The current version of the statepoint file format is 18.1.
              ``time``, ``wgt``, ``delayed_group``, ``surf_id``, and
              ``particle``, which represent the position, direction, energy,
              time, weight, delayed group, surface ID, and particle type
-             (0=neutron, 1=photon, 2=electron, 3=positron), respectively. Only
-             present when `run_mode` is 'eigenvalue'.
+             (PDG number), respectively. Only present when `run_mode` is
+             'eigenvalue'.
 
 **/tallies/**
 
@@ -72,7 +73,10 @@ The current version of the statepoint file format is 18.1.
 
 **/tallies/meshes/mesh <uid>/**
 
-:Datasets: - **type** (*char[]*) -- Type of mesh.
+:Attributes: - **id** (*int*) -- ID of the mesh
+
+:Datasets: - **name** (*char[]*) -- Name of the mesh.
+           - **type** (*char[]*) -- Type of mesh.
            - **dimension** (*int*) -- Number of mesh cells in each dimension.
            - **Regular Mesh Only:**
               - **lower_left** (*double[]*) -- Coordinates of lower-left corner of
@@ -96,6 +100,8 @@ The current version of the statepoint file format is 18.1.
               - **library** (*char[]*) -- Mesh library used to represent the
                                           mesh ("moab" or "libmesh").
               - **length_multiplier** (*double*) Scaling factor applied to the mesh.
+              - **options** (*char[]*) -- Special options that control spatial
+                                          search data structures used.
               - **volumes** (*double[]*) -- Volume of each mesh cell.
               - **vertices** (*double[]*) -- x, y, z values of the mesh vertices.
               - **connectivity** (*int[]*) -- Connectivity array for the mesh
@@ -143,6 +149,8 @@ The current version of the statepoint file format is 18.1.
                tallies will have a value of 0 unless otherwise instructed.
              - **multiply_density** (*int*) -- Flag indicating whether reaction
                rates should be multiplied by atom density (1) or not (0).
+             - **higher_moments** (*int*) -- Flag indicating whether
+               higher-order tally moments are enabled (1) or not (0).
 
 :Datasets: - **n_realizations** (*int*) -- Number of realizations.
            - **n_filters** (*int*) -- Number of filters used.

docs/source/io_formats/summary.rst

@@ -4,7 +4,7 @@
 Summary File Format
 ===================
 
-The current version of the summary file format is 6.0.
+The current version of the summary file format is 6.1.
 
 **/**
 
@@ -38,6 +38,7 @@ The current version of the summary file format is 6.0.
              is an array if the cell uses distributed materials, otherwise it is
              a scalar.
            - **temperature** (*double[]*) -- Temperature of the cell in Kelvin.
+           - **density** (*double[]*) -- Density of the cell in [g/cm3].
            - **translation** (*double[3]*) -- Translation applied to the fill
              universe. This dataset is present only if fill_type is set to
              'universe'.

docs/source/io_formats/tallies.rst

@@ -100,6 +100,18 @@ The ``<tally>`` element accepts the following sub-elements:
 
      *Default*: None
 
+   :ignore_zeros:
+     Whether to allow zero tally bins to be ignored when assessing the
+     convergece of the precision trigger. If True, only nonzero tally scores
+     will be compared to the trigger's threshold.
+
+     .. note:: The ``ignore_zeros`` option can cause the tally trigger to fire
+               prematurely if there are no hits in any bins at the first
+               evalulation. It is the user's responsibility to specify enough
+               particles per batch to get a nonzero score in at least one bin.
+
+     *Default*: False
+
    :scores:
      The score(s) in this tally to which the trigger should be applied.
 
@@ -130,9 +142,9 @@ attributes/sub-elements:
 
   :type:
     The type of the filter. Accepted options are "cell", "cellfrom",
-    "cellborn", "surface", "material", "universe", "energy", "energyout", "mu",
-    "polar", "azimuthal", "mesh", "distribcell", "delayedgroup",
-    "energyfunction", and "particle".
+    "cellborn", "surface", "material", "universe", "energy", "energyout",
+    "mu", "polar", "azimuthal", "mesh", "distribcell", "delayedgroup",
+    "energyfunction", "particle", and "particleproduction".
 
   :bins:
      A description of the bins for each type of filter can be found in
@@ -306,8 +318,34 @@ should be set to:
   they use ``energy`` and ``y``.
 
 :particle:
-  A list of integers indicating the type of particles to tally ('neutron' = 1,
-  'photon' = 2, 'electron' = 3, 'positron' = 4).
+  A list of particle identifiers to tally, specified as strings (e.g.,
+  ``neutron``, ``photon``, ``He4``) or as integer PDG numbers.
+
+:particleproduction:
+  This filter tallies secondary particles produced in reactions, binned by
+  particle type and, optionally, by energy. Unlike other energy filters, the
+  weight applied is the weight of the secondary particle. To obtain secondary
+  particle production rates, use this filter with the ``events`` score.
+
+  The filter uses the following sub-elements instead of ``bins``:
+
+  :particles:
+    A space-separated list of secondary particle types to tally (e.g.,
+    ``photon``, ``neutron``, ``electron``).
+
+  :energies:
+    An optional monotonically increasing list of energy boundaries in [eV]
+    for binning the secondary particle energies. If omitted, total production
+    is tallied without energy binning.
+
+  For example, to tally photon and neutron production in three energy groups:
+
+  .. code-block:: xml
+
+      <filter id="1" type="particleproduction">
+        <particles>photon neutron</particles>
+        <energies>0.0 1.0e5 1.0e6 20.0e6</energies>
+      </filter>
 
 ------------------
 ``<mesh>`` Element
@@ -317,6 +355,11 @@ If a mesh is desired as a filter for a tally, it must be specified in a separate
 element with the tag name ``<mesh>``. This element has the following
 attributes/sub-elements:
 
+  :name:
+    An optional string name to identify the mesh in output files.
+
+    *Default*: ""
+
   :type:
     The type of mesh. This can be either "regular", "rectilinear",
     "cylindrical", "spherical", or "unstructured".
@@ -364,6 +407,10 @@ attributes/sub-elements:
     The mesh library used to represent an unstructured mesh. This can be either
     "moab" or "libmesh". (For unstructured mesh only.)
 
+  :options:
+    Special options that control spatial search data structures used. (For
+    unstructured mesh using MOAB only)
+
   :filename:
     The name of the mesh file to be loaded at runtime. (For unstructured mesh
     only.)

docs/source/io_formats/track.rst

@@ -4,7 +4,7 @@
 Track File Format
 =================
 
-The current revision of the particle track file format is 3.0.
+The current revision of the particle track file format is 3.1.
 
 **/**
 
@@ -32,6 +32,5 @@ The current revision of the particle track file format is 3.0.
                             the array for each primary/secondary particle. The
                             last offset should match the total size of the
                             array.
-                          - **particles** (*int[]*) -- Particle type for each
-                            primary/secondary particle (0=neutron, 1=photon,
-                            2=electron, 3=positron).
+                          - **particles** (*int32_t[]*) -- Particle type for
+                            each primary/secondary particle (PDG number).

docs/source/license.rst

@@ -4,7 +4,7 @@
 License Agreement
 =================
 
-Copyright © 2011-2023 Massachusetts Institute of Technology, UChicago Argonne
+Copyright © 2011-2026 Massachusetts Institute of Technology, UChicago Argonne
 LLC, and OpenMC contributors
 
 Permission is hereby granted, free of charge, to any person obtaining a copy of

docs/source/methods/charged_particles_physics.rst

@@ -0,0 +1,362 @@
+.. _methods_charged_particle_physics:
+
+========================
+Charged Particle Physics
+========================
+
+OpenMC neglects the spatial transport of charged particles (electrons and
+positrons), assuming they deposit all their energy locally and produce
+bremsstrahlung photons at their birth location. This approximation, called
+thick-target bremsstrahlung (TTB) approximation is justified by the fact that
+charged particles have much shorter stopping ranges compared to neutrons and
+photons, especially in high-density materials.
+
+-----------------------------
+Charged Particle Interactions
+-----------------------------
+
+Bremsstrahlung
+--------------
+
+When a charged particle is decelerated in the field of an atom, some of its
+kinetic energy is converted into electromagnetic radiation known as
+bremsstrahlung, or 'braking radiation'. In each event, an electron or positron
+with kinetic energy :math:`T` generates a photon with an energy :math:`E`
+between :math:`0` and :math:`T`. Bremsstrahlung is described by a cross section
+that is differential in photon energy, in the direction of the emitted photon,
+and in the final direction of the charged particle. However, in Monte Carlo
+simulations it is typical to integrate over the angular variables to obtain a
+single differential cross section with respect to photon energy, which is often
+expressed in the form
+
+.. math::
+    :label: bremsstrahlung-dcs
+
+    \frac{d\sigma_{\text{br}}}{dE} = \frac{Z^2}{\beta^2} \frac{1}{E}
+    \chi(Z, T, \kappa),
+
+where :math:`\kappa = E/T` is the reduced photon energy and :math:`\chi(Z, T,
+\kappa)` is the scaled bremsstrahlung cross section, which is experimentally
+measured.
+
+Because electrons are attracted to atomic nuclei whereas positrons are
+repulsed, the cross section for positrons is smaller, though it approaches that
+of electrons in the high energy limit. To obtain the positron cross section, we
+multiply :eq:`bremsstrahlung-dcs` by the :math:`\kappa`-independent factor used
+in Salvat_,
+
+.. math::
+    :label: positron-factor
+
+    \begin{aligned}
+    F_{\text{p}}(Z,T) =
+    & 1 - \text{exp}(-1.2359\times 10^{-1}t + 6.1274\times 10^{-2}t^2 - 3.1516\times 10^{-2}t^3 \\
+    & + 7.7446\times 10^{-3}t^4 - 1.0595\times 10^{-3}t^5 + 7.0568\times 10^{-5}t^6 \\
+    & - 1.8080\times 10^{-6}t^7),
+    \end{aligned}
+
+where
+
+.. math::
+    :label: positron-factor-t
+
+    t = \ln\left(1 + \frac{10^6}{Z^2}\frac{T}{\text{m}_\text{e}c^2} \right).
+
+:math:`F_{\text{p}}(Z,T)` is the ratio of the radiative stopping powers for
+positrons and electrons. Stopping power describes the average energy loss per
+unit path length of a charged particle as it passes through matter:
+
+.. math::
+    :label: stopping-power
+
+    -\frac{dT}{ds} = n \int E \frac{d\sigma}{dE} dE \equiv S(T),
+
+where :math:`n` is the number density of the material and :math:`d\sigma/dE` is
+the cross section differential in energy loss. The total stopping power
+:math:`S(T)` can be separated into two components: the radiative stopping
+power :math:`S_{\text{rad}}(T)`, which refers to energy loss due to
+bremsstrahlung, and the collision stopping power :math:`S_{\text{col}}(T)`,
+which refers to the energy loss due to inelastic collisions with bound
+electrons in the material that result in ionization and excitation. The
+radiative stopping power for electrons is given by
+
+.. math::
+    :label: radiative-stopping-power
+
+    S_{\text{rad}}(T) = n \frac{Z^2}{\beta^2} T \int_0^1 \chi(Z,T,\kappa)
+    d\kappa.
+
+
+To obtain the radiative stopping power for positrons,
+:eq:`radiative-stopping-power`  is multiplied by :eq:`positron-factor`.
+
+While the models for photon interactions with matter described above can safely
+assume interactions occur with free atoms, sampling the target atom based on
+the macroscopic cross sections, molecular effects cannot necessarily be
+disregarded for charged particle treatment. For compounds and mixtures, the
+bremsstrahlung cross section is calculated using Bragg's additivity rule as
+
+.. math::
+    :label: material-bremsstrahlung-dcs
+
+    \frac{d\sigma_{\text{br}}}{dE} = \frac{1}{\beta^2 E} \sum_i \gamma_i Z^2_i
+    \chi(Z_i, T, \kappa),
+
+where the sum is over the constituent elements and :math:`\gamma_i` is the
+atomic fraction of the :math:`i`-th element. Similarly, the radiative stopping
+power is calculated using Bragg's additivity rule as
+
+.. math::
+    :label: material-radiative-stopping-power
+
+    S_{\text{rad}}(T) = \sum_i w_i S_{\text{rad},i}(T),
+
+where :math:`w_i` is the mass fraction of the :math:`i`-th element and
+:math:`S_{\text{rad},i}(T)` is found for element :math:`i` using
+:eq:`radiative-stopping-power`. The collision stopping power, however, is a
+function of certain quantities such as the mean excitation energy :math:`I` and
+the density effect correction :math:`\delta_F` that depend on molecular
+properties. These quantities cannot simply be summed over constituent elements
+in a compound, but should instead be calculated for the material. The Bethe
+formula can be used to find the collision stopping power of the material:
+
+.. math::
+    :label: material-collision-stopping-power
+
+    S_{\text{col}}(T) = \frac{2 \pi r_e^2 m_e c^2}{\beta^2} N_A \frac{Z}{A_M}
+    [\ln(T^2/I^2) + \ln(1 + \tau/2) + F(\tau) - \delta_F(T)],
+
+where :math:`N_A` is Avogadro's number, :math:`A_M` is the molar mass,
+:math:`\tau = T/m_e`, and :math:`F(\tau)` depends on the particle type. For
+electrons,
+
+.. math::
+    :label: F-electron
+
+    F_{-}(\tau) = (1 - \beta^2)[1 + \tau^2/8 - (2\tau + 1) \ln2],
+
+while for positrons
+
+.. math::
+    :label: F-positron
+
+    F_{+}(\tau) = 2\ln2 - (\beta^2/12)[23 + 14/(\tau + 2) + 10/(\tau + 2)^2 +
+    4/(\tau + 2)^3].
+
+The density effect correction :math:`\delta_F` takes into account the reduction
+of the collision stopping power due to the polarization of the material the
+charged particle is passing through by the electric field of the particle.
+It can be evaluated using the method described by Sternheimer_, where the
+equation for :math:`\delta_F` is
+
+.. math::
+    :label: density-effect-correction
+
+    \delta_F(\beta) = \sum_{i=1}^n f_i \ln[(l_i^2 + l^2)/l_i^2] -
+    l^2(1-\beta^2).
+
+Here, :math:`f_i` is the oscillator strength of the :math:`i`-th transition,
+given by :math:`f_i = n_i/Z`, where :math:`n_i` is the number of electrons in
+the :math:`i`-th subshell. The frequency :math:`l` is the solution of the
+equation
+
+.. math::
+    :label: density-effect-l
+
+    \frac{1}{\beta^2} - 1 = \sum_{i=1}^{n} \frac{f_i}{\bar{\nu}_i^2 + l^2},
+
+where :math:`\bar{v}_i` is defined as
+
+.. math::
+    :label: density-effect-nubar
+
+    \bar{\nu}_i = h\nu_i \rho / h\nu_p.
+
+The plasma energy :math:`h\nu_p` of the medium is given by
+
+.. math::
+    :label: plasma-frequency
+
+    h\nu_p = \sqrt{\frac{(hc)^2 r_e \rho_m N_A Z}{\pi A}},
+
+where :math:`A` is the atomic weight and :math:`\rho_m` is the density of the
+material. In :eq:`density-effect-nubar`, :math:`h\nu_i` is the oscillator
+energy, and :math:`\rho` is an adjustment factor introduced to give agreement
+between the experimental values of the oscillator energies and the mean
+excitation energy. The :math:`l_i` in :eq:`density-effect-correction` are
+defined as
+
+.. math::
+    :label: density-effect-li
+
+    \begin{aligned}
+    l_i &= (\bar{\nu}_i^2 + 2/3f_i)^{1/2} ~~~~&\text{for}~~ \bar{\nu}_i > 0 \\
+    l_n &= f_n^{1/2} ~~~~&\text{for}~~ \bar{\nu}_n = 0,
+    \end{aligned}
+
+where the second case applies to conduction electrons. For a conductor,
+:math:`f_n` is given by :math:`n_c/Z`, where :math:`n_c` is the effective
+number of conduction electrons, and :math:`v_n = 0`. The adjustment factor
+:math:`\rho` is determined using the equation for the mean excitation energy:
+
+.. math::
+    :label: mean-excitation-energy
+
+    \ln I = \sum_{i=1}^{n-1} f_i \ln[(h\nu_i\rho)^2 + 2/3f_i(h\nu_p)^2]^{1/2} +
+    f_n \ln (h\nu_pf_n^{1/2}).
+
+.. _ttb:
+
+
+Thick-Target Bremsstrahlung Approximation
++++++++++++++++++++++++++++++++++++++++++
+
+Since charged particles lose their energy on a much shorter distance scale than
+neutral particles, not much error should be introduced by neglecting to
+transport electrons. However, the bremsstrahlung emitted from high energy
+electrons and positrons can travel far from the interaction site. Thus, even
+without a full electron transport mode it is necessary to model bremsstrahlung.
+We use a thick-target bremsstrahlung (TTB) approximation based on the models in
+Salvat_ and Kaltiaisenaho_ for generating bremsstrahlung photons, which assumes
+the charged particle loses all its energy in a single homogeneous material
+region.
+
+To model bremsstrahlung using the TTB approximation, we need to know the number
+of photons emitted by the charged particle and the energy distribution of the
+photons. These quantities can be calculated using the continuous slowing down
+approximation (CSDA). The CSDA assumes charged particles lose energy
+continuously along their trajectory with a rate of energy loss equal to the
+total stopping power, ignoring fluctuations in the energy loss. The
+approximation is useful for expressing average quantities that describe how
+charged particles slow down in matter. For example, the CSDA range approximates
+the average path length a charged particle travels as it slows to rest:
+
+.. math::
+    :label: csda-range
+
+    R(T) = \int^T_0 \frac{dT'}{S(T')}.
+
+Actual path lengths will fluctuate around :math:`R(T)`. The average number of
+photons emitted per unit path length is given by the inverse bremsstrahlung
+mean free path:
+
+.. math::
+    :label: inverse-bremsstrahlung-mfp
+
+    \lambda_{\text{br}}^{-1}(T,E_{\text{cut}})
+    = n\int_{E_{\text{cut}}}^T\frac{d\sigma_{\text{br}}}{dE}dE
+    = n\frac{Z^2}{\beta^2}\int_{\kappa_{\text{cut}}}^1\frac{1}{\kappa}
+    \chi(Z,T,\kappa)d\kappa.
+
+The lower limit of the integral in :eq:`inverse-bremsstrahlung-mfp` is non-zero
+because the bremsstrahlung differential cross section diverges for small photon
+energies but is finite for photon energies above some cutoff energy
+:math:`E_{\text{cut}}`. The mean free path
+:math:`\lambda_{\text{br}}^{-1}(T,E_{\text{cut}})` is used to calculate the
+photon number yield, defined as the average number of photons emitted with
+energy greater than :math:`E_{\text{cut}}` as the charged particle slows down
+from energy :math:`T` to :math:`E_{\text{cut}}`. The photon number yield is
+given by
+
+.. math::
+    :label: photon-number-yield
+
+    Y(T,E_{\text{cut}}) = \int^{R(T)}_{R(E_{\text{cut}})}
+    \lambda_{\text{br}}^{-1}(T',E_{\text{cut}})ds = \int_{E_{\text{cut}}}^T
+    \frac{\lambda_{\text{br}}^{-1}(T',E_{\text{cut}})}{S(T')}dT'.
+
+:math:`Y(T,E_{\text{cut}})` can be used to construct the energy spectrum of
+bremsstrahlung photons: the number of photons created with energy between
+:math:`E_1` and :math:`E_2` by a charged particle with initial kinetic energy
+:math:`T` as it comes to rest is given by :math:`Y(T,E_1) - Y(T,E_2)`.
+
+To simulate the emission of bremsstrahlung photons, the total stopping power
+and bremsstrahlung differential cross section for positrons and electrons must
+be calculated for a given material using :eq:`material-bremsstrahlung-dcs` and
+:eq:`material-radiative-stopping-power`. These quantities are used to build the
+tabulated bremsstrahlung energy PDF and CDF for that material for each incident
+energy :math:`T_k` on the energy grid. The following algorithm is then applied
+to sample the photon energies:
+
+1. For an incident charged particle with energy :math:`T`, sample the number of
+   emitted photons as
+
+   .. math::
+
+       N = \lfloor Y(T,E_{\text{cut}}) + \xi_1 \rfloor.
+
+2. Rather than interpolate the PDF between indices :math:`k` and :math:`k+1`
+   for which :math:`T_k < T < T_{k+1}`, which is computationally expensive, use
+   the composition method and sample from the PDF at either :math:`k` or
+   :math:`k+1`. Using linear interpolation on a logarithmic scale, the PDF can
+   be expressed as
+
+   .. math::
+
+       p_{\text{br}}(T,E) = \pi_k p_{\text{br}}(T_k,E) + \pi_{k+1}
+       p_{\text{br}}(T_{k+1},E),
+
+   where the interpolation weights are
+
+   .. math::
+
+       \pi_k = \frac{\ln T_{k+1} - \ln T}{\ln T_{k+1} - \ln T_k},~~~
+       \pi_{k+1} = \frac{\ln T - \ln T_k}{\ln T_{k+1} - \ln T_k}.
+
+   Sample either the index :math:`i = k` or :math:`i = k+1` according to the
+   point probabilities :math:`\pi_{k}` and :math:`\pi_{k+1}`.
+
+3. Determine the maximum value of the CDF :math:`P_{\text{br,max}}`.
+
+3. Sample the photon energies using the inverse transform method with the
+   tabulated CDF :math:`P_{\text{br}}(T_i, E)` i.e.,
+
+   .. math::
+
+       E = E_j \left[ (1 + a_j) \frac{\xi_2 P_{\text{br,max}} -
+       P_{\text{br}}(T_i, E_j)} {E_j p_{\text{br}}(T_i, E_j)} + 1
+       \right]^{\frac{1}{1 + a_j}}
+
+   where the interpolation factor :math:`a_j` is given by
+
+   .. math::
+
+       a_j = \frac{\ln p_{\text{br}}(T_i,E_{j+1}) - \ln p_{\text{br}}(T_i,E_j)}
+       {\ln E_{j+1} - \ln E_j}
+
+   and :math:`P_{\text{br}}(T_i, E_j) \le \xi_2 P_{\text{br,max}} \le
+   P_{\text{br}}(T_i, E_{j+1})`.
+
+We ignore the range of the electron or positron, i.e., the bremsstrahlung
+photons are produced in the same location that the charged particle was
+created. The direction of the photons is assumed to be the same as the
+direction of the incident charged particle, which is a reasonable approximation
+at higher energies when the bremsstrahlung radiation is emitted at small
+angles.
+
+
+Electron-Positron Annihilation
+------------------------------
+
+When a positron collides with an electron, both particles are annihilated and
+generally two photons with equal energy are created. If the kinetic energy of
+the positron is high enough, the two photons can have different energies, and
+the higher-energy photon is emitted preferentially in the direction of flight
+of the positron. It is also possible to produce a single photon if the
+interaction occurs with a bound electron, and in some cases three (or, rarely,
+even more) photons can be emitted. However, the annihilation cross section is
+largest for low-energy positrons, and as the positron energy decreases, the
+angular distribution of the emitted photons becomes isotropic.
+
+In OpenMC, we assume the most likely case in which a low-energy positron (which
+has already lost most of its energy to bremsstrahlung radiation) interacts with
+an electron which is free and at rest. Two photons with energy equal to the
+electron rest mass energy :math:`m_e c^2 = 0.511` MeV are emitted isotropically
+in opposite directions.
+
+
+.. _Kaltiaisenaho: https://aaltodoc.aalto.fi/bitstream/handle/123456789/21004/master_Kaltiaisenaho_Toni_2016.pdf
+
+.. _Salvat: https://doi.org/10.1787/32da5043-en
+
+.. _Sternheimer: https://doi.org/10.1103/PhysRevB.26.6067

docs/source/methods/cross_sections.rst

@@ -289,17 +289,59 @@ sections. This allows flexibility for the model to use highly anisotropic
 scattering information in the water while the fuel can be simulated with linear
 or even isotropic scattering.
 
+Particle Speed
+--------------
+
+When using a multigroup representation of cross sections, the particle speed has
+meaning only in an average sense. The particle speed is important when modeling
+dynamic behavior. OpenMC calculates the particle speed using the inverse
+velocity multigroup data if it is available. If such data is not available,
+OpenMC uses an approximate velocity using the group energy bounds in the
+following way:
+
+.. math::
+
+   \frac{1}{v_g} = \int_{E_{\text{min}}^g}^{E_{\text{max}}^g} \frac{1}{v(E)} \frac{\alpha}{E} dE
+
+Where :math:`E_{\text{min}}^g` and :math:`E_{\text{max}}^g` are the group energy
+boundaries for group :math:`g`. :math:`v(E)` is the neutron velocity calculated
+using relativistic kinematics, :math:`\alpha` is a normalization constant for the
+:math:`\frac{1}{E}` spectrum.
+
+This equation is valid when inside the group boundaries the neutron spectrum
+follows a typical :math:`\frac{1}{E}` slowing down spectrum. This assumption is
+widely used when generating fine group neutron cross section data libraries from
+continuous energy data.
+
+The solution to this equation is:
+
+.. math::
+
+   \frac{1}{v_g} = \frac{1}{c \log\left(\frac{E_{\text{max}}^g}{E_{\text{min}}^g}\right)}
+   \left[ 2(\operatorname{arctanh}(k_{\text{max}}^{-1}) - \operatorname{arctanh}(k_{\text{min}}^{-1}))
+   - (k_{\text{max}}-k_{\text{min}}) \right]
+
+where :math:`c` is the speed of light and :math:`k_{\text{max}}`,
+:math:`k_{\text{min}}` are defined by a change of variables:
+
+.. math::
+
+   k = \sqrt{1+\frac{2 m_n c^2}{E}}
+
+where :math:`E` is the particle kinetic energy and :math:`m_n` is the neutron
+rest mass.
+
 .. _logarithmic mapping technique:
-   https://laws.lanl.gov/vhosts/mcnp.lanl.gov/pdf_files/la-ur-14-24530.pdf
+   https://mcnp.lanl.gov/pdf_files/TechReport_2014_LANL_LA-UR-14-24530_Brown.pdf
 .. _Hwang: https://doi.org/10.13182/NSE87-A16381
 .. _Josey: https://doi.org/10.1016/j.jcp.2015.08.013
 .. _WMP Library: https://github.com/mit-crpg/WMP_Library
 .. _MCNP: https://mcnp.lanl.gov
-.. _Serpent: http://montecarlo.vtt.fi
-.. _NJOY: https://www.njoy21.io/NJOY21/
+.. _Serpent: https://serpent.vtt.fi
+.. _NJOY: https://www.njoy21.io/
 .. _ENDF/B data: https://www.nndc.bnl.gov/endf-b8.0/
 .. _Leppanen: https://doi.org/10.1016/j.anucene.2009.03.019
-.. _algorithms: http://ab-initio.mit.edu/wiki/index.php/Faddeeva_Package
+.. _algorithms: http://ab-initio.mit.edu/faddeeva/
 .. _NCrystal: https://github.com/mctools/ncrystal
 .. _NCrystal paper: https://doi.org/10.1016/j.cpc.2019.07.015
 .. _using plugins: https://doi.org/10.1016/j.cpc.2021.108082

docs/source/methods/depletion.rst

@@ -114,7 +114,7 @@ The predictor method only requires one evaluation and its error converges as
 twice as expensive as the predictor method, but achieves an error of
 :math:`\mathcal{O}(h^2)`. An exhaustive description of time integration methods
 and their merits can be found in the `thesis of Colin Josey
-<http://dspace.mit.edu/handle/1721.1/7582>`_.
+<https://dspace.mit.edu/handle/1721.1/7582>`_.
 
 OpenMC does not rely on a single time integration method but rather has several
 classes that implement different algorithms. For example, the

docs/source/methods/eigenvalue.rst

@@ -55,15 +55,17 @@ in :ref:`fission-bank-algorithms`.
 Source Convergence Issues
 -------------------------
 
+.. _methods-shannon-entropy:
+
 Diagnosing Convergence with Shannon Entropy
 -------------------------------------------
 
 As discussed earlier, it is necessary to converge both :math:`k_{eff}` and the
 source distribution before any tallies can begin. Moreover, the convergence rate
-of the source distribution is in general slower than that of
-:math:`k_{eff}`. One should thus examine not only the convergence of
-:math:`k_{eff}` but also the convergence of the source distribution in order to
-make decisions on when to start active batches.
+of the source distribution is in general slower than that of :math:`k_{eff}`.
+One should thus examine not only the convergence of :math:`k_{eff}` but also the
+convergence of the source distribution in order to make decisions on when to
+start active batches.
 
 However, the representation of the source distribution makes it a bit more
 difficult to analyze its convergence. Since :math:`k_{eff}` is a scalar
@@ -108,6 +110,13 @@ at plots of :math:`k_{eff}` and the Shannon entropy. A number of methods have
 been proposed (see e.g. [Romano]_, [Ueki]_), but each of these is not without
 problems.
 
+Shannon entropy is calculated differently for the random ray solver, as
+described :ref:`in the random ray theory section
+<methods-shannon-entropy-random-ray>`. Additionally, as the Shannon entropy only
+serves as a diagnostic tool for convergence of the fission source distribution,
+there is currently no diagnostic to determine if the scattering source
+distribution in random ray is converged.
+
 ---------------------------
 Uniform Fission Site Method
 ---------------------------
@@ -142,7 +151,7 @@ than unity. By ensuring that the expected number of fission sites in each mesh
 cell is constant, the collision density across all cells, and hence the variance
 of tallies, is more uniform than it would be otherwise.
 
-.. _Shannon entropy: https://laws.lanl.gov/vhosts/mcnp.lanl.gov/pdf_files/la-ur-06-3737.pdf
+.. _Shannon entropy: https://mcnp.lanl.gov/pdf_files/TechReport_2006_LANL_LA-UR-06-3737_Brown.pdf
 
 .. [Lieberoth] J. Lieberoth, "A Monte Carlo Technique to Solve the Static
    Eigenvalue Problem of the Boltzmann Transport Equation," *Nukleonik*, **11**,

docs/source/methods/energy_deposition.rst

@@ -25,19 +25,38 @@ KERMA (Kinetic Energy Release in Materials) [Mack97]_ coefficients for reaction
 :math:`\times` cross-section (e.g., eV-barn) and can be used much like a reaction
 cross section for the purpose of tallying energy deposition.
 
-KERMA coefficients can be computed using the energy-balance method with
-a nuclear data processing code like NJOY, which performs the following
-iteration over all reactions :math:`r` for all isotopes :math:`i`
-requested
+KERMA coefficients can be computed using the energy-balance method with a
+nuclear data processing code like NJOY, which estimates the KERMA coefficients
+using the following equation:
 
 .. math::
 
-    k_{i, r}(E) = \left(E + Q_{i, r} - \bar{E}_{i, r, n}
+    k_{i, r}(E) = \left(E + Q_{i, r} - \sum\limits_x \bar{E}_{i, r, x}
+    \right)\sigma_{i, r}(E),
+
+where the summation is over each secondary particle type :math:`x`. This
+equation states that the energy deposited is equal to the energy of the incident
+particle plus the reaction :math:`Q` value less the energy of secondary
+particles that are transported away from the reaction site. For neutron
+interactions, the energy-balance KERMA coefficient is
+
+.. math::
+
+    k_{i, r}(E) = \left(E + Q_{i, r} - \sum\limits_x \bar{E}_{i, r, n}
     - \bar{E}_{i, r, \gamma}\right)\sigma_{i, r}(E),
 
-removing the energy of neutral particles (neutrons and photons) that are
-transported away from the reaction site :math:`\bar{E}`, and the reaction
-:math:`Q` value.
+where :math:`\bar{E}_{i, r, n}` is the average energy of secondary neutrons and
+:math:`\bar{E}_{i, r, \gamma}` is the average energy of secondary photons. For
+photon and charged particle interactions the KERMA coefficient is
+
+.. math::
+    :label: energy-balance-photon
+
+    k_{i, r}(E) = \left(E + Q_{i, r} - \sum\limits_x \bar{E}_{i, r, x}
+    \right)\sigma_{i, r}(E).
+
+where the :math:`Q` value is zero for all interactions except for pair
+production and positron annihilation.
 
 -------
 Fission
@@ -120,7 +139,7 @@ run with :math:`N918` reflecting fission heating computed from NJOY.
 This modified heating data is stored as the MT=901 reaction and will be scored
 if ``heating-local`` is included in :attr:`openmc.Tally.scores`.
 
-Coupled neutron-photon transport
+Coupled Neutron-Photon Transport
 --------------------------------
 
 Here, OpenMC instructs ``heatr`` to assume that energy from photons is not
@@ -138,6 +157,50 @@ Let :math:`N301` represent the total heating number returned from this
 This modified heating data is stored as the MT=301 reaction and will be scored
 if ``heating`` is included in :attr:`openmc.Tally.scores`.
 
+Photons and Charged Particles
+-----------------------------
+
+In OpenMC, energy deposition from photons or charged particles is scored using
+the energy balance method based on Equation :eq:`energy-balance-photon`. Special
+consideration is given to electrons and positrons as described below.
+
++++++++++++++++++
+Charged Particles
++++++++++++++++++
+
+OpenMC tracks photons interaction by interaction so the energy deposited in each
+collision is easily attributed back to the nuclide and reaction for which the
+photon interacted with. Charged particles (electrons and photons) aren't tracked
+in the same way. For charged particles, OpenMC assumes that all their energy
+(less the energy of bremsstrahlung radiation) is deposited in the material in
+which they were born. In this way it is harder to trace how much energy should
+be attributed in each nuclide.
+
+According to the CSDA approximation (see :ref:`ttb`) the energy deposited by a
+charged particle with kinetic energy :math:`T` in the :math:`i`-th element can
+be calculated as:
+
+.. math::
+
+    E_{i} = \int_{0}^{R(T)} w_{i}S_{\text{col,i}} dx
+
+where :math:`R(T)` is the CSDA range of the charged particle,
+:math:`S_{\text{col},i}` is the collision stopping power of the charged particle
+in the :math:`i`-th element and :math:`w_i` is the mass fraction of the
+:math:`i`-th element. According to the Bethe formula the collision stopping
+power of the :math:`i`-th element is proportional to :math:`Z_i/A_i`, so the
+fractional collision stopping power from the :math:`i`-th element is:
+
+.. math::
+
+    \frac{w_{i}S_{\text{col},i}(T)}{S_{\text{col}}(T)} =
+    \frac{\frac{w_{i}Z_{i}}{A_{i}}}{\sum_{i}\frac{w_{i}Z_{i}}{A_{i}}} =
+    \frac{\gamma_i Z_{i}}{\sum_{i}\gamma_i Z_{i}}.
+
+where :math:`\gamma_i` is the atomic fraction of the :math:`i`-th element.
+Therefore, the energy deposited by charged particles should be attributed to
+a given element according to its fractional charge density.
+
 ----------
 References
 ----------

docs/source/methods/geometry.rst

@@ -1066,5 +1066,5 @@ surface is known as in :ref:`reflection`.
 .. _constructive solid geometry: https://en.wikipedia.org/wiki/Constructive_solid_geometry
 .. _surfaces: https://en.wikipedia.org/wiki/Surface
 .. _MCNP: https://mcnp.lanl.gov
-.. _Serpent: http://montecarlo.vtt.fi
+.. _Serpent: https://serpent.vtt.fi
 .. _Monte Carlo Performance benchmark: https://github.com/mit-crpg/benchmarks/tree/master/mc-performance/openmc

docs/source/methods/index.rst

@@ -14,9 +14,12 @@ Theory and Methodology
     random_numbers
     neutron_physics
     photon_physics
+    charged_particles_physics
     tallies
     eigenvalue
     depletion
     energy_deposition
     parallelization
     cmfd
+    variance_reduction
+    random_ray

docs/source/methods/neutron_physics.rst

@@ -290,7 +290,10 @@ create and store fission sites for the following generation. First, the average
 number of prompt and delayed neutrons must be determined to decide whether the
 secondary neutrons will be prompt or delayed. This is important because delayed
 neutrons have a markedly different spectrum from prompt neutrons, one that has a
-lower average energy of emission. The total number of neutrons emitted
+lower average energy of emission. Furthermore, in simulations where tracking
+time of neutrons is important, we need to consider the emission time delay of
+the secondary neutrons, which is dependent on the decay constant of the
+delayed neutron precursor. The total number of neutrons emitted
 :math:`\nu_t` is given as a function of incident energy in the ENDF format. Two
 representations exist for :math:`\nu_t`. The first is a polynomial of order
 :math:`N` with coefficients :math:`c_0,c_1,\dots,c_N`. If :math:`\nu_t` has this
@@ -306,8 +309,8 @@ interpolation law. The number of prompt neutrons released per fission event
 :math:`\nu_p` is also given as a function of incident energy and can be
 specified in a polynomial or tabular format. The number of delayed neutrons
 released per fission event :math:`\nu_d` can only be specified in a tabular
-format. In practice, we only need to determine :math:`nu_t` and
-:math:`nu_d`. Once these have been determined, we can calculated the delayed
+format. In practice, we only need to determine :math:`\nu_t` and
+:math:`\nu_d`. Once these have been determined, we can calculate the delayed
 neutron fraction
 
 .. math::
@@ -335,8 +338,14 @@ neutrons. Otherwise, we produce :math:`\lfloor \nu \rfloor + 1` neutrons. Then,
 for each fission site produced, we sample the outgoing angle and energy
 according to the algorithms given in :ref:`sample-angle` and
 :ref:`sample-energy` respectively. If the neutron is to be born delayed, then
-there is an extra step of sampling a delayed neutron precursor group since they
-each have an associated secondary energy distribution.
+there is an extra step of sampling a delayed neutron precursor group to get the
+associated secondary energy distribution and the decay constant
+:math:`\lambda`, which is needed to sample the emission delay time :math:`t_d`:
+
+.. math::
+    :label: sample-delay-time
+
+    t_d = -\frac{\ln \xi}{\lambda}.
 
 The sampled outgoing angle and energy of fission neutrons along with the
 position of the collision site are stored in an array called the fission
@@ -1743,19 +1752,19 @@ types.
 
 .. _Watt fission spectrum: https://doi.org/10.1103/PhysRev.87.1037
 
-.. _Foderaro: http://hdl.handle.net/1721.1/1716
+.. _Foderaro: https://dspace.mit.edu/handle/1721.1/1716
 
 .. _OECD: https://www.oecd-nea.org/tools/abstract/detail/NEA-1792
 
 .. _NJOY: https://www.njoy21.io/NJOY2016/
 
-.. _PREPRO: https://www-nds.iaea.org/ndspub/endf/prepro/
+.. _PREPRO: https://www-nds.iaea.org/public/endf/prepro/
 
 .. _ENDF-6 Format: https://www.oecd-nea.org/dbdata/data/manual-endf/endf102.pdf
 
-.. _Monte Carlo Sampler: https://permalink.lanl.gov/object/tr?what=info:lanl-repo/lareport/LA-09721-MS
+.. _Monte Carlo Sampler: https://mcnp.lanl.gov/pdf_files/TechReport_1983_LANL_LA-9721-MS_EverettCashwell.pdf
 
-.. _LA-UR-14-27694: https://permalink.lanl.gov/object/tr?what=info:lanl-repo/lareport/LA-UR-14-27694
+.. _LA-UR-14-27694: https://www.osti.gov/biblio/1159204
 
 .. _MC21: https://www.osti.gov/biblio/903083
 
@@ -1763,6 +1772,4 @@ types.
 
 .. _Sutton and Brown: https://www.osti.gov/biblio/307911
 
-.. _lectures: https://laws.lanl.gov/vhosts/mcnp.lanl.gov/pdf_files/la-ur-05-4983.pdf
-
-.. _MCNP Manual: https://laws.lanl.gov/vhosts/mcnp.lanl.gov/pdf_files/la-ur-03-1987.pdf
+.. _lectures: https://mcnp.lanl.gov/pdf_files/TechReport_2005_LANL_LA-UR-05-4983_Brown.pdf

docs/source/methods/parallelization.rst

@@ -609,17 +609,17 @@ is actually independent of the number of nodes:
 
 .. _first paper: https://doi.org/10.2307/2280232
 
-.. _work of Forrest Brown: http://hdl.handle.net/2027.42/24996
+.. _work of Forrest Brown: https://deepblue.lib.umich.edu/handle/2027.42/24996
 
 .. _Brissenden and Garlick: https://doi.org/10.1016/0306-4549(86)90095-2
 
-.. _MPICH: http://www.mpich.org
+.. _MPICH: https://www.mpich.org
 
 .. _binomial tree: https://www.mcs.anl.gov/~thakur/papers/ijhpca-coll.pdf
 
 .. _Geary: https://doi.org/10.2307/2342070
 
-.. _Barnett: http://citeseerx.ist.psu.edu/viewdoc/summary?doi=10.1.1.51.7772
+.. _Barnett: https://citeseerx.ist.psu.edu/viewdoc/summary?doi=10.1.1.51.7772
 
 .. _single-instruction multiple-data: https://en.wikipedia.org/wiki/SIMD
 

docs/source/methods/photon_physics.rst

@@ -667,342 +667,6 @@ and Auger electrons:
 
 5. Repeat from step 1 for vacancy left by the transition electron.
 
-Electron-Positron Annihilation
-------------------------------
-
-When a positron collides with an electron, both particles are annihilated and
-generally two photons with equal energy are created. If the kinetic energy of
-the positron is high enough, the two photons can have different energies, and
-the higher-energy photon is emitted preferentially in the direction of flight
-of the positron. It is also possible to produce a single photon if the
-interaction occurs with a bound electron, and in some cases three (or, rarely,
-even more) photons can be emitted. However, the annihilation cross section is
-largest for low-energy positrons, and as the positron energy decreases, the
-angular distribution of the emitted photons becomes isotropic.
-
-In OpenMC, we assume the most likely case in which a low-energy positron (which
-has already lost most of its energy to bremsstrahlung radiation) interacts with
-an electron which is free and at rest. Two photons with energy equal to the
-electron rest mass energy :math:`m_e c^2 = 0.511` MeV are emitted isotropically
-in opposite directions.
-
-Bremsstrahlung
---------------
-
-When a charged particle is decelerated in the field of an atom, some of its
-kinetic energy is converted into electromagnetic radiation known as
-bremsstrahlung, or 'braking radiation'. In each event, an electron or positron
-with kinetic energy :math:`T` generates a photon with an energy :math:`E`
-between :math:`0` and :math:`T`. Bremsstrahlung is described by a cross section
-that is differential in photon energy, in the direction of the emitted photon,
-and in the final direction of the charged particle. However, in Monte Carlo
-simulations it is typical to integrate over the angular variables to obtain a
-single differential cross section with respect to photon energy, which is often
-expressed in the form
-
-.. math::
-    :label: bremsstrahlung-dcs
-
-    \frac{d\sigma_{\text{br}}}{dE} = \frac{Z^2}{\beta^2} \frac{1}{E}
-    \chi(Z, T, \kappa),
-
-where :math:`\kappa = E/T` is the reduced photon energy and :math:`\chi(Z, T,
-\kappa)` is the scaled bremsstrahlung cross section, which is experimentally
-measured.
-
-Because electrons are attracted to atomic nuclei whereas positrons are
-repulsed, the cross section for positrons is smaller, though it approaches that
-of electrons in the high energy limit. To obtain the positron cross section, we
-multiply :eq:`bremsstrahlung-dcs` by the :math:`\kappa`-independent factor used
-in Salvat_,
-
-.. math::
-    :label: positron-factor
-
-    \begin{aligned}
-    F_{\text{p}}(Z,T) =
-    & 1 - \text{exp}(-1.2359\times 10^{-1}t + 6.1274\times 10^{-2}t^2 - 3.1516\times 10^{-2}t^3 \\
-    & + 7.7446\times 10^{-3}t^4 - 1.0595\times 10^{-3}t^5 + 7.0568\times 10^{-5}t^6 \\
-    & - 1.8080\times 10^{-6}t^7),
-    \end{aligned}
-
-where
-
-.. math::
-    :label: positron-factor-t
-
-    t = \ln\left(1 + \frac{10^6}{Z^2}\frac{T}{\text{m}_\text{e}c^2} \right).
-
-:math:`F_{\text{p}}(Z,T)` is the ratio of the radiative stopping powers for
-positrons and electrons. Stopping power describes the average energy loss per
-unit path length of a charged particle as it passes through matter:
-
-.. math::
-    :label: stopping-power
-
-    -\frac{dT}{ds} = n \int E \frac{d\sigma}{dE} dE \equiv S(T),
-
-where :math:`n` is the number density of the material and :math:`d\sigma/dE` is
-the cross section differential in energy loss. The total stopping power
-:math:`S(T)` can be separated into two components: the radiative stopping
-power :math:`S_{\text{rad}}(T)`, which refers to energy loss due to
-bremsstrahlung, and the collision stopping power :math:`S_{\text{col}}(T)`,
-which refers to the energy loss due to inelastic collisions with bound
-electrons in the material that result in ionization and excitation. The
-radiative stopping power for electrons is given by
-
-.. math::
-    :label: radiative-stopping-power
-
-    S_{\text{rad}}(T) = n \frac{Z^2}{\beta^2} T \int_0^1 \chi(Z,T,\kappa)
-    d\kappa.
-
-
-To obtain the radiative stopping power for positrons,
-:eq:`radiative-stopping-power`  is multiplied by :eq:`positron-factor`.
-
-While the models for photon interactions with matter described above can safely
-assume interactions occur with free atoms, sampling the target atom based on
-the macroscopic cross sections, molecular effects cannot necessarily be
-disregarded for charged particle treatment. For compounds and mixtures, the
-bremsstrahlung cross section is calculated using Bragg's additivity rule as
-
-.. math::
-    :label: material-bremsstrahlung-dcs
-
-    \frac{d\sigma_{\text{br}}}{dE} = \frac{1}{\beta^2 E} \sum_i \gamma_i Z^2_i
-    \chi(Z_i, T, \kappa),
-
-where the sum is over the constituent elements and :math:`\gamma_i` is the
-atomic fraction of the :math:`i`-th element. Similarly, the radiative stopping
-power is calculated using Bragg's additivity rule as
-
-.. math::
-    :label: material-radiative-stopping-power
-
-    S_{\text{rad}}(T) = \sum_i w_i S_{\text{rad},i}(T),
-
-where :math:`w_i` is the mass fraction of the :math:`i`-th element and
-:math:`S_{\text{rad},i}(T)` is found for element :math:`i` using
-:eq:`radiative-stopping-power`. The collision stopping power, however, is a
-function of certain quantities such as the mean excitation energy :math:`I` and
-the density effect correction :math:`\delta_F` that depend on molecular
-properties. These quantities cannot simply be summed over constituent elements
-in a compound, but should instead be calculated for the material. The Bethe
-formula can be used to find the collision stopping power of the material:
-
-.. math::
-    :label: material-collision-stopping-power
-
-    S_{\text{col}}(T) = \frac{2 \pi r_e^2 m_e c^2}{\beta^2} N_A \frac{Z}{A_M}
-    [\ln(T^2/I^2) + \ln(1 + \tau/2) + F(\tau) - \delta_F(T)],
-
-where :math:`N_A` is Avogadro's number, :math:`A_M` is the molar mass,
-:math:`\tau = T/m_e`, and :math:`F(\tau)` depends on the particle type. For
-electrons,
-
-.. math::
-    :label: F-electron
-
-    F_{-}(\tau) = (1 - \beta^2)[1 + \tau^2/8 - (2\tau + 1) \ln2],
-
-while for positrons
-
-.. math::
-    :label: F-positron
-
-    F_{+}(\tau) = 2\ln2 - (\beta^2/12)[23 + 14/(\tau + 2) + 10/(\tau + 2)^2 +
-    4/(\tau + 2)^3].
-
-The density effect correction :math:`\delta_F` takes into account the reduction
-of the collision stopping power due to the polarization of the material the
-charged particle is passing through by the electric field of the particle.
-It can be evaluated using the method described by Sternheimer_, where the
-equation for :math:`\delta_F` is
-
-.. math::
-    :label: density-effect-correction
-
-    \delta_F(\beta) = \sum_{i=1}^n f_i \ln[(l_i^2 + l^2)/l_i^2] -
-    l^2(1-\beta^2).
-
-Here, :math:`f_i` is the oscillator strength of the :math:`i`-th transition,
-given by :math:`f_i = n_i/Z`, where :math:`n_i` is the number of electrons in
-the :math:`i`-th subshell. The frequency :math:`l` is the solution of the
-equation
-
-.. math::
-    :label: density-effect-l
-
-    \frac{1}{\beta^2} - 1 = \sum_{i=1}^{n} \frac{f_i}{\bar{\nu}_i^2 + l^2},
-
-where :math:`\bar{v}_i` is defined as
-
-.. math::
-    :label: density-effect-nubar
-
-    \bar{\nu}_i = h\nu_i \rho / h\nu_p.
-
-The plasma energy :math:`h\nu_p` of the medium is given by
-
-.. math::
-    :label: plasma-frequency
-
-    h\nu_p = \sqrt{\frac{(hc)^2 r_e \rho_m N_A Z}{\pi A}},
-
-where :math:`A` is the atomic weight and :math:`\rho_m` is the density of the
-material. In :eq:`density-effect-nubar`, :math:`h\nu_i` is the oscillator
-energy, and :math:`\rho` is an adjustment factor introduced to give agreement
-between the experimental values of the oscillator energies and the mean
-excitation energy. The :math:`l_i` in :eq:`density-effect-correction` are
-defined as
-
-.. math::
-    :label: density-effect-li
-
-    \begin{aligned}
-    l_i &= (\bar{\nu}_i^2 + 2/3f_i)^{1/2} ~~~~&\text{for}~~ \bar{\nu}_i > 0 \\
-    l_n &= f_n^{1/2} ~~~~&\text{for}~~ \bar{\nu}_n = 0,
-    \end{aligned}
-
-where the second case applies to conduction electrons. For a conductor,
-:math:`f_n` is given by :math:`n_c/Z`, where :math:`n_c` is the effective
-number of conduction electrons, and :math:`v_n = 0`. The adjustment factor
-:math:`\rho` is determined using the equation for the mean excitation energy:
-
-.. math::
-    :label: mean-excitation-energy
-
-    \ln I = \sum_{i=1}^{n-1} f_i \ln[(h\nu_i\rho)^2 + 2/3f_i(h\nu_p)^2]^{1/2} +
-    f_n \ln (h\nu_pf_n^{1/2}).
-
-.. _ttb:
-
-Thick-Target Bremsstrahlung Approximation
-+++++++++++++++++++++++++++++++++++++++++
-
-Since charged particles lose their energy on a much shorter distance scale than
-neutral particles, not much error should be introduced by neglecting to
-transport electrons. However, the bremsstrahlung emitted from high energy
-electrons and positrons can travel far from the interaction site. Thus, even
-without a full electron transport mode it is necessary to model bremsstrahlung.
-We use a thick-target bremsstrahlung (TTB) approximation based on the models in
-Salvat_ and Kaltiaisenaho_ for generating bremsstrahlung photons, which assumes
-the charged particle loses all its energy in a single homogeneous material
-region.
-
-To model bremsstrahlung using the TTB approximation, we need to know the number
-of photons emitted by the charged particle and the energy distribution of the
-photons. These quantities can be calculated using the continuous slowing down
-approximation (CSDA). The CSDA assumes charged particles lose energy
-continuously along their trajectory with a rate of energy loss equal to the
-total stopping power, ignoring fluctuations in the energy loss. The
-approximation is useful for expressing average quantities that describe how
-charged particles slow down in matter. For example, the CSDA range approximates
-the average path length a charged particle travels as it slows to rest:
-
-.. math::
-    :label: csda-range
-
-    R(T) = \int^T_0 \frac{dT'}{S(T')}.
-
-Actual path lengths will fluctuate around :math:`R(T)`. The average number of
-photons emitted per unit path length is given by the inverse bremsstrahlung
-mean free path:
-
-.. math::
-    :label: inverse-bremsstrahlung-mfp
-
-    \lambda_{\text{br}}^{-1}(T,E_{\text{cut}})
-    = n\int_{E_{\text{cut}}}^T\frac{d\sigma_{\text{br}}}{dE}dE
-    = n\frac{Z^2}{\beta^2}\int_{\kappa_{\text{cut}}}^1\frac{1}{\kappa}
-    \chi(Z,T,\kappa)d\kappa.
-
-The lower limit of the integral in :eq:`inverse-bremsstrahlung-mfp` is non-zero
-because the bremsstrahlung differential cross section diverges for small photon
-energies but is finite for photon energies above some cutoff energy
-:math:`E_{\text{cut}}`. The mean free path
-:math:`\lambda_{\text{br}}^{-1}(T,E_{\text{cut}})` is used to calculate the
-photon number yield, defined as the average number of photons emitted with
-energy greater than :math:`E_{\text{cut}}` as the charged particle slows down
-from energy :math:`T` to :math:`E_{\text{cut}}`. The photon number yield is
-given by
-
-.. math::
-    :label: photon-number-yield
-
-    Y(T,E_{\text{cut}}) = \int^{R(T)}_{R(E_{\text{cut}})}
-    \lambda_{\text{br}}^{-1}(T',E_{\text{cut}})ds = \int_{E_{\text{cut}}}^T
-    \frac{\lambda_{\text{br}}^{-1}(T',E_{\text{cut}})}{S(T')}dT'.
-
-:math:`Y(T,E_{\text{cut}})` can be used to construct the energy spectrum of
-bremsstrahlung photons: the number of photons created with energy between
-:math:`E_1` and :math:`E_2` by a charged particle with initial kinetic energy
-:math:`T` as it comes to rest is given by :math:`Y(T,E_1) - Y(T,E_2)`.
-
-To simulate the emission of bremsstrahlung photons, the total stopping power
-and bremsstrahlung differential cross section for positrons and electrons must
-be calculated for a given material using :eq:`material-bremsstrahlung-dcs` and
-:eq:`material-radiative-stopping-power`. These quantities are used to build the
-tabulated bremsstrahlung energy PDF and CDF for that material for each incident
-energy :math:`T_k` on the energy grid. The following algorithm is then applied
-to sample the photon energies:
-
-1. For an incident charged particle with energy :math:`T`, sample the number of
-   emitted photons as
-
-   .. math::
-
-       N = \lfloor Y(T,E_{\text{cut}}) + \xi_1 \rfloor.
-
-2. Rather than interpolate the PDF between indices :math:`k` and :math:`k+1`
-   for which :math:`T_k < T < T_{k+1}`, which is computationally expensive, use
-   the composition method and sample from the PDF at either :math:`k` or
-   :math:`k+1`. Using linear interpolation on a logarithmic scale, the PDF can
-   be expressed as
-
-   .. math::
-
-       p_{\text{br}}(T,E) = \pi_k p_{\text{br}}(T_k,E) + \pi_{k+1}
-       p_{\text{br}}(T_{k+1},E),
-
-   where the interpolation weights are
-
-   .. math::
-
-       \pi_k = \frac{\ln T_{k+1} - \ln T}{\ln T_{k+1} - \ln T_k},~~~
-       \pi_{k+1} = \frac{\ln T - \ln T_k}{\ln T_{k+1} - \ln T_k}.
-
-   Sample either the index :math:`i = k` or :math:`i = k+1` according to the
-   point probabilities :math:`\pi_{k}` and :math:`\pi_{k+1}`.
-
-3. Determine the maximum value of the CDF :math:`P_{\text{br,max}}`.
-
-3. Sample the photon energies using the inverse transform method with the
-   tabulated CDF :math:`P_{\text{br}}(T_i, E)` i.e.,
-
-   .. math::
-
-       E = E_j \left[ (1 + a_j) \frac{\xi_2 P_{\text{br,max}} -
-       P_{\text{br}}(T_i, E_j)} {E_j p_{\text{br}}(T_i, E_j)} + 1
-       \right]^{\frac{1}{1 + a_j}}
-
-   where the interpolation factor :math:`a_j` is given by
-
-   .. math::
-
-       a_j = \frac{\ln p_{\text{br}}(T_i,E_{j+1}) - \ln p_{\text{br}}(T_i,E_j)}
-       {\ln E_{j+1} - \ln E_j}
-
-   and :math:`P_{\text{br}}(T_i, E_j) \le \xi_2 P_{\text{br,max}} \le
-   P_{\text{br}}(T_i, E_{j+1})`.
-
-We ignore the range of the electron or positron, i.e., the bremsstrahlung
-photons are produced in the same location that the charged particle was
-created. The direction of the photons is assumed to be the same as the
-direction of the incident charged particle, which is a reasonable approximation
-at higher energies when the bremsstrahlung radiation is emitted at small
-angles.
 
 .. _photon_production:
 
@@ -1059,16 +723,14 @@ emitted photon.
 
 .. _anomalous scattering: http://pd.chem.ucl.ac.uk/pdnn/diff1/anomscat.htm
 
-.. _Kahn's rejection method: https://laws.lanl.gov/vhosts/mcnp.lanl.gov/pdf_files/aecu-3259_kahn.pdf
+.. _Kahn's rejection method: https://doi.org/10.2172/4353680
 
 .. _Klein-Nishina: https://en.wikipedia.org/wiki/Klein%E2%80%93Nishina_formula
 
-.. _LA-UR-04-0487: https://laws.lanl.gov/vhosts/mcnp.lanl.gov/pdf_files/la-ur-04-0487.pdf
+.. _LA-UR-04-0487: https://mcnp.lanl.gov/pdf_files/TechReport_2004_LANL_LA-UR-04-0487_Sood.pdf
 
-.. _LA-UR-04-0488: https://laws.lanl.gov/vhosts/mcnp.lanl.gov/pdf_files/la-ur-04-0488.pdf
+.. _LA-UR-04-0488: https://mcnp.lanl.gov/pdf_files/TechReport_2004_LANL_LA-UR-04-0488_SoodWhite.pdf
 
 .. _Kaltiaisenaho: https://aaltodoc.aalto.fi/bitstream/handle/123456789/21004/master_Kaltiaisenaho_Toni_2016.pdf
 
-.. _Salvat: https://www.oecd-nea.org/globalsearch/download.php?doc=77434
-
-.. _Sternheimer: https://doi.org/10.1103/PhysRevB.26.6067
+.. _Salvat: https://doi.org/10.1787/32da5043-en

docs/source/methods/random_numbers.rst

@@ -7,7 +7,7 @@ Random Number Generation
 In order to sample probability distributions, one must be able to produce random
 numbers. The standard technique to do this is to generate numbers on the
 interval :math:`[0,1)` from a deterministic sequence that has properties that
-make it appear to be random, e.g. being uniformly distributed and not exhibiting
+make it appear to be random, e.g., being uniformly distributed and not exhibiting
 correlation between successive terms. Since the numbers produced this way are
 not truly "random" in a strict sense, they are typically referred to as
 pseudorandom numbers, and the techniques used to generate them are pseudorandom
@@ -15,6 +15,11 @@ number generators (PRNGs). Numbers sampled on the unit interval can then be
 transformed for the purpose of sampling other continuous or discrete probability
 distributions.
 
+There are many different algorithms for pseudorandom number generation. OpenMC
+currently uses `permuted congruential generator`_ (PCG), which builds on top of
+the simpler linear congruential generator (LCG). Both algorithms are described
+below.
+
 ------------------------------
 Linear Congruential Generators
 ------------------------------
@@ -37,8 +42,8 @@ be generated with a method chosen at random. Some theory should be used."
 Typically, :math:`M` is chosen to be a power of two as this enables :math:`x
 \mod M` to be performed using the bitwise AND operator with a bit mask. The
 constants for the linear congruential generator used by default in OpenMC are
-:math:`g = 2806196910506780709`, :math:`c = 1`, and :math:`M = 2^{63}` (see
-`L'Ecuyer`_).
+:math:`g = 2806196910506780709`, :math:`c = 1`, and :math:`M = 2^{63}` (from
+`L'Ecuyer <https://doi.org/10.1090/S0025-5718-99-00996-5>`_).
 
 Skip-ahead Capability
 ---------------------
@@ -50,7 +55,8 @@ want to skip ahead :math:`N` random numbers and :math:`N` is large, the cost of
 sampling :math:`N` random numbers to get to that position may be prohibitively
 expensive. Fortunately, algorithms have been developed that allow us to skip
 ahead in :math:`O(\log_2 N)` operations instead of :math:`O(N)`. One algorithm
-to do so is described in a paper by Brown_. This algorithm relies on the following
+to do so is described in a `paper by Brown
+<https://www.osti.gov/biblio/976209>`_. This algorithm relies on the following
 relationship:
 
 .. math::
@@ -58,15 +64,26 @@ relationship:
 
     \xi_{i+k} = g^k \xi_i + c \frac{g^k - 1}{g - 1} \mod M
 
-Note that equation :eq:`lcg-skipahead` has the same general form as equation :eq:`lcg`, so
-the idea is to determine the new multiplicative and additive constants in
-:math:`O(\log_2 N)` operations.
-
-.. only:: html
-
-   .. rubric:: References
+Note that equation :eq:`lcg-skipahead` has the same general form as equation
+:eq:`lcg`, so the idea is to determine the new multiplicative and additive
+constants in :math:`O(\log_2 N)` operations.
 
 
-.. _L'Ecuyer: https://doi.org/10.1090/S0025-5718-99-00996-5
-.. _Brown: https://laws.lanl.gov/vhosts/mcnp.lanl.gov/pdf_files/anl-rn-arb-stride.pdf
+--------------------------------
+Permuted Congruential Generators
+--------------------------------
+
+The `permuted congruential generator`_ (PCG) algorithm aims to improve upon the
+LCG algorithm by permuting the output. The algorithm works on the basic
+principle of first advancing the generator state using the LCG algorithm and
+then applying a permutation function on the LCG state to obtain the output. This
+results in increased statistical quality as measured by common statistical tests
+while exhibiting a very small performance overhead relative to the LCG algorithm
+and an equivalent memory footprint. For further details, see the original
+technical report by `O'Neill
+<https://www.pcg-random.org/pdf/hmc-cs-2014-0905.pdf>`_.  OpenMC uses the
+PCG-RXS-M-XS variant with a 64-bit state and 64-bit output.
+
 .. _linear congruential generator: https://en.wikipedia.org/wiki/Linear_congruential_generator
+
+.. _permuted congruential generator: https://en.wikipedia.org/wiki/Permuted_congruential_generator

docs/source/methods/tallies.rst

@@ -4,9 +4,9 @@
 Tallies
 =======
 
-Note that the methods discussed in this section are written specifically for
-continuous-energy mode but equivalent apply to the multi-group mode if the
-particle's energy is replaced with the particle's group
+The methods discussed in this section are written specifically for continuous-
+energy mode. However, they can also apply to the multi-group mode if the
+particle's energy is instead interpreted as the particle's group.
 
 ------------------
 Filters and Scores
@@ -205,7 +205,73 @@ had a collision at every event. Thus, for tallies with outgoing-energy filters
 or for tallies of scattering moments (which require the scattering cosine of
 the change-in-angle), we must use an analog estimator.
 
-.. TODO: Add description of surface current tallies
+-----------------------------------
+Surface-Integrated Flux and Current
+-----------------------------------
+
+Surface tallies allow you to measure particle behavior as they cross specific
+boundaries in your geometry. Unlike volume tallies, which integrate over a
+volumetric region, surface tallies capture the current or flux passing through a
+surface. Surface tallies are estimated using an analog estimator.
+
+Current Score
+-------------
+
+When tallying the current across a surface, we simply count the weight of
+particles that cross the surface of interest:
+
+
+.. math::
+    :label: analog-current-estimator
+
+    J = \frac{1}{W} \sum_{i \in S} w_i.
+
+where :math:`J` is the area-integrated current passing through surface
+:math:`S`, :math:`W` is the total starting weight of the particles, and
+:math:`w_i` is the weight of the particle as it crosses the surface :math:`S`.
+
+Flux Score
+----------
+
+When tallying flux over a surface, we use the relationship between current and
+flux:
+
+
+.. math::
+    :label: surface-flux-estimator
+
+    \phi_S = \frac{1}{W} \sum_{i \in S} \frac{w_i}{|\mu|}.
+
+where :math:`\phi_S` is the area-integrated flux over surface :math:`S`,
+:math:`W` is the total starting weight of the particles, :math:`w_i` is the
+weight of the particle as it crosses the surface :math:`S` and :math:`\mu` is
+the cosine of angle between the particle direction and the surface normal.
+
+This equation diverges when the particle crossing the surface is nearly parallel
+to it (that is, as :math:`\mu` approaches zero). To remove this divergence,
+OpenMC scores:
+
+.. math::
+    :label: modified-surface-flux-estimator
+
+    \phi_S = \frac{1}{W} \sum_{i \in S} w_i f(\mu).
+
+and the function :math:`f` is defined by:
+
+.. math::
+     f(\mu) = \begin{cases}
+      \frac{1}{|\mu|} & |\mu| > \mu_\text{cut} \\
+      \frac{1}{c\mu_\text{cut}} & |\mu| \le \mu_\text{cut}
+     \end{cases}
+
+where :math:`\mu_\text{cut}` is the grazing cosine cutoff and :math:`c` is the
+cosine substitution ratio. The parameters :math:`\mu_\text{cut}` and :math:`c`
+can be set by the user via the :attr:`openmc.Settings.surface_grazing_cutoff`
+and :attr:`openmc.Settings.surface_grazing_ratio` attributes, respectively. The
+default values for these parameters are 0.001 and 0.5 as recommended by
+`Favorite, Thomas, and Booth <https://doi.org/10.13182/NSE09-72>`_.
+
+.. _tallies_statistics:
 
 ----------
 Statistics
@@ -268,6 +334,14 @@ normal, log-normal, Weibull, etc. The central limit theorem states that as
 Estimating Statistics of a Random Variable
 ------------------------------------------
 
+After running OpenMC, each tallied quantity has a reported mean and standard
+deviation. The below sections explain how these quantities are computed. Note
+that OpenMC uses **batch statistics**, meaning that each observation for a tally
+random variable corresponds to the aggregation of tally contributions from
+multiple source particles that are grouped together into a single batch. See
+:ref:`usersguide_particles` for more information on how the number of source
+particles and statistical batches are specified.
+
 Mean
 ++++
 
@@ -377,6 +451,130 @@ of this is that the longer you run a simulation, the better you know your
 results. Therefore, by running a simulation long enough, it is possible to
 reduce the stochastic uncertainty to arbitrarily low levels.
 
+Skewness
+++++++++
+
+The `skewness`_ of a population quantifies the asymmetry of the probability
+distribution around its mean. Positive and negative skewness indicate a
+longer/heavier right and left tail respectively. Let :math:`x_1,\ldots,x_n` be
+the per-realization values for a bin, with sample mean :math:`\bar{x}` and
+sample central moments:
+
+.. math::
+
+   m_k \;=\; \frac{1}{n}\sum_{i=1}^{n}\bigl(x_i-\bar{x}\bigr)^k.
+
+OpenMC reports the *adjusted Fisher-Pearson skewness* (defined for :math:`n \ge
+3`), which is commonly used in many statistical packages:
+
+.. math::
+
+   G_1 \;=\; \frac{\sqrt{n \cdot (n-1)}}{\,n-2\,}\cdot\frac{m_3}{m_2^{3/2}}.
+
+where :math:`m_2` and :math:`m_3` correspond to the biased sample second and
+third central moment respectively.
+
+Kurtosis
+++++++++
+
+The `kurtosis`_ of a population quantifies tail weight (also called tailedness)
+of the probability distribution relative to a normal distribution. Positive
+excess kurtosis indicates *heavier tails* whereas negative excess kurtosis
+indicates *lighter tails*. Kurtosis is especially useful for identifying bins
+where occasional extreme scores dominate uncertainty. OpenMC reports the
+*adjusted excess kurtosis* (defined for :math:`n \ge 4`):
+
+.. math::
+
+   G_2 \;=\; \frac{(n-1)}{(n-2)(n-3)}
+             \left[(n+1)\,\frac{m_4}{m_2^{2}} \;-\; 3(n-1)\right].
+
+where :math:`m_2` and :math:`m_4` correspond to the biased sample second and
+fourth central moment respectively. For a perfectly normal distribution, the
+excess kurtosis is :math:`0`.
+
+Variance of Variance
+++++++++++++++++++++
+
+The variance of the variance (also known as the coefficient of variation
+squared) measures *stability of the sample variance* :math:`s^2` and, by
+extension, the reliability of reported relative errors. High VOV means that
+error bars themselves are noisy—often due to heavy tails, skewness, or too few
+realizations.
+
+.. math::
+
+    VOV = \frac{s^2(s_{\bar{X}}^2)}{s_{\bar{X}}^4 } = \frac{m_4}{m_2^2} - \frac{1}{n}
+
+where :math:`s_{\bar{X}}^2` is the estimated variance of the mean and
+:math:`s^2(s_{\bar{X}}^2)` is the estimated variance in :math:`s_{\bar{X}}^2`.
+The MCNP manual suggests a hard threshold such that :math:`VOV < 0.1` to improve
+the probability of forming a reliable confidence interval. However, OpenMC does
+not enforce an universal cut-off because the suitability of any single threshold
+depends strongly on problem specifics (estimator choice, variance-reduction
+settings, tally binning, or even effective sample size).
+
+
+Normality Tests (D'Agostino-Pearson)
+++++++++++++++++++++++++++++++++++++
+
+These normality test verify the hypothesis that fluctuations are *approximately
+normal*, a working assumption behind many Monte Carlo diagnostics and
+`confidence-interval heuristics`_. Tests are provided for: (i) skewness-only,
+(ii) kurtosis-only, and (iii) the *omnibus* combination. OpenMC uses the
+finite-sample-adjusted skewness :math:`G_1` and excess kurtosis :math:`G_2`
+above to construct standardized normal scores :math:`Z_1` (from :math:`G_1`) and
+:math:`Z_2` (from :math:`G_2`) via the D'Agostino-Pearson transformations. The
+omnibus statistic is
+
+.. math::
+
+   K^2 \;=\; Z_1^{\,2} \;+\; Z_2^{\,2}
+   \;\sim\; \chi^2_{(2)} \quad \text{under } H_0:\ \text{normality}.
+
+OpenMC reports :math:`Z_1`, :math:`Z_2`, :math:`K^2`, and their p-values when
+prerequisites are met (skewness for :math:`n\ge 3`, kurtosis and omnibus for
+:math:`n\ge 4`). Given a user-chosen significance level :math:`\alpha` (default
+is :math:`0.05`), reject :math:`H_0` if :math:`\text{p-value}<\alpha`; otherwise
+fail to reject. OpenMC leaves the interpretation to the user, who should
+consider VOV together with skewness, kurtosis, and normality tests results when
+judging whether reported confidence intervals are credible for their application
+[#norm-tests]_.
+
+.. [#norm-tests]
+   Higher-moments accumulation must be enabled with ``higher_moments = True``
+   for running these diagnostics including the skewness, kurtosis, and normality
+   tests.
+
+Figure of Merit
++++++++++++++++
+
+The figure of merit (FOM) is an indicator that accounts for both the statistical
+uncertainty and the execution time and represents how much information is
+obtained per unit time in the simulation. The FOM is defined as
+
+.. math::
+    :label: figure_of_merit
+
+    FOM = \frac{1}{r^2 t},
+
+where :math:`t` is the total execution time and :math:`r` is the relative error
+defined as
+
+.. math::
+    :label: relative_error
+
+    r = \frac{s_{\bar{X}}}{\bar{x}}.
+
+Based on this definition, one can see that a higher FOM is desirable. The FOM is
+useful as a comparative tool. For example, if a variance reduction technique is
+being applied to a simulation, the FOM with variance reduction can be compared
+to the FOM without variance reduction to ascertain whether the reduction in
+variance outweighs the potential increase in execution time (e.g., due to
+particle splitting). It is important to note that MCNP reports the FOM using CPU
+time (wall-clock time multiplied by the number of threads/cores), whereas OpenMC
+reports the FOM using only the wall-clock time :math:`t`.
+
 Confidence Intervals
 ++++++++++++++++++++
 
@@ -484,6 +682,8 @@ improve the estimate of the percentile.
 
    .. rubric:: References
 
+.. _confidence-interval heuristics: https://doi.org/10.1080/00031305.1990.10475751
+
 .. _following approximation: https://doi.org/10.1080/03610918708812641
 
 .. _Bessel's correction: https://en.wikipedia.org/wiki/Bessel's_correction
@@ -504,6 +704,10 @@ improve the estimate of the percentile.
 
 .. _converges in distribution: https://en.wikipedia.org/wiki/Convergence_of_random_variables#Convergence_in_distribution
 
+.. _skewness: https://en.wikipedia.org/wiki/Skewness
+
+.. _kurtosis: https://en.wikipedia.org/wiki/Kurtosis
+
 .. _confidence intervals: https://en.wikipedia.org/wiki/Confidence_interval
 
 .. _Student's t-distribution: https://en.wikipedia.org/wiki/Student%27s_t-distribution
@@ -512,4 +716,4 @@ improve the estimate of the percentile.
 
 .. _unpublished rational approximation: https://stackedboxes.org/2017/05/01/acklams-normal-quantile-function/
 
-.. _MC21: http://www.osti.gov/bridge/servlets/purl/903083-HT5p1o/903083.pdf
+.. _MC21: https://www.osti.gov/servlets/purl/903083

docs/source/methods/variance_reduction.rst

@@ -0,0 +1,216 @@
+.. _methods_variance_reduction:
+
+==================
+Variance Reduction
+==================
+
+.. _methods_variance_reduction_intro:
+
+------------
+Introduction
+------------
+
+Transport problems can sometimes involve a significant degree of attenuation
+between the source and a detector (tally) region, which can result in a flux
+differential of ten orders of magnitude (or more) throughout the simulation
+domain. As Monte Carlo uncertainties tend to be inversely proportional to the
+physical flux density, it can be extremely difficult to accurately resolve
+tallies in locations that are optically far from the source. This issue is
+particularly common in fixed source simulations, where some tally locations may
+not experience a single scoring event, even after billions of analog histories.
+
+Variance reduction techniques aim to either flatten the global uncertainty
+distribution, such that all regions of phase space have a fairly similar
+uncertainty, or to reduce the uncertainty in specific locations (such as a
+detector). There are three strategies available in OpenMC for variance
+reduction: weight windows generated via the MAGIC method or the FW-CADIS method,
+and source biasing. Both weight windowing strategies work by developing a mesh
+that can be utilized by subsequent Monte Carlo solves to split particles heading
+towards areas of lower flux densities while terminating particles in higher flux
+regions. In contrast, source biasing modifies source site sampling behavior to
+preferentially track particles more likely to reach phase space regions of
+interest.
+
+------------
+MAGIC Method
+------------
+
+The Method of Automatic Generation of Importances by Calculation, or `MAGIC
+method <https://doi.org/10.1016/j.fusengdes.2011.01.059>`_, is an iterative
+technique that uses spatial flux information :math:`\phi(r)` obtained from a
+normal Monte Carlo solve to produce weight windows :math:`w(r)` that can be
+utilized by a subsequent iteration of Monte Carlo. While the first generation of
+weight windows produced may only help to reduce variance slightly, use of these
+weights to generate another set of weight windows results in a progressively
+improving iterative scheme.
+
+Equation :eq:`magic` defines how the lower bound of weight windows
+:math:`w_{\ell}(r)` are generated with MAGIC using forward flux information.
+Here, we can see that the flux at location :math:`r` is normalized by the
+maximum flux in any group at that location. We can also see that the weights are
+divided by a factor of two, which accounts for the typical :math:`5\times`
+factor separating the lower and upper weight window bounds in OpenMC.
+
+.. math::
+    :label: magic
+
+    w_{\ell}(r) = \frac{\phi(r)}{2\,\text{max}(\phi(r))}
+
+A major advantage of this technique is that it does not require any special
+transport machinery; it simply uses multiple Monte Carlo simulations to
+iteratively improve a set of weight windows (which are typically defined on a
+mesh covering the simulation domain). The downside to this method is that as the
+flux differential increases between areas near and far from the source, it
+requires more outer Monte Carlo iterations, each of which can be expensive in
+itself. Additionally, computation of weight windows based on regular (forward)
+neutron flux tally information does not produce the most numerically effective
+set of weight windows. Nonetheless, MAGIC remains a simple and effective
+technique for generating weight windows.
+
+--------
+FW-CADIS
+--------
+
+As discussed in the previous section, computation of weight windows based on
+regular (forward) neutron flux tally information does not produce the most
+numerically efficient set of weight windows. It is highly preferable to generate
+weight windows based on spatial adjoint flux :math:`\phi^{\dag}(r)`
+information. The adjoint flux is essentially the "reverse" simulation problem,
+where we sample a random point and assume this is where a particle was absorbed,
+and then trace it backwards (upscattering in energy), until we sample the point
+where it was born from.
+
+The Forward-Weighted Consistent Adjoint Driven Importance Sampling method, or
+`FW-CADIS method <https://doi.org/10.13182/NSE12-33>`_, produces weight windows
+for global or local variance reduction given adjoint flux information throughout 
+the entire domain. The weight window lower bound is defined in Equation
+:eq:`fw_cadis`, and also involves a normalization step not shown here.
+
+.. math::
+    :label: fw_cadis
+
+    w_{\ell}(r) = \frac{1}{2\phi^{\dag}(r)}
+
+While the algorithm itself is quite simple, it requires estimates of the global
+adjoint flux distribution, which is difficult to generate directly with Monte
+Carlo transport. Thus, FW-CADIS typically uses an alternative solver (often
+deterministic) that can be more readily adapted for generating adjoint flux
+information, and which is often much cheaper than Monte Carlo given that a rough
+solution is often sufficient for weight window generation.
+
+The FW-CADIS implementation in OpenMC utilizes its own internal random ray
+multigroup transport solver to generate the adjoint source distribution. No
+coupling to any external transport is solver is necessary. The random ray solver
+operates on the same geometry as the Monte Carlo solver, so no redefinition of
+the simulation geometry is required. More details on how the adjoint flux is
+computed are given in the :ref:`adjoint methods section <adjoint>`.
+
+More information on the workflow is available in the :ref:`user guide
+<variance_reduction>`, but generally production of weight windows with FW-CADIS
+involves several stages (some of which are highly automated). These tasks
+include generation of approximate multigroup cross section data for use by the
+random ray solver, running of the random ray solver in normal (forward flux)
+mode to generate a source for the adjoint solver, running of the random ray
+solver in adjoint mode to generate adjoint flux tallies, and finally the
+production of weight windows via the FW-CADIS method. As is discussed in the
+user guide, most of these steps are automated together, making the additional
+burden on the user fairly small.
+
+The major advantage of this technique is that it typically produces much more
+numerically efficient weight windows as compared to those generated with MAGIC,
+sometimes with an order-of-magnitude improvement in the figure of merit
+(Equation :eq:`variance_fom`), which accounts for both the variance and the
+execution time. Another major advantage is that the cost of the random ray
+solver is typically negligible compared to the cost of the subsequent Monte
+Carlo solve itself, making it a very cheap method to deploy. The downside to
+this method is that it introduces a second transport method into the mix (random
+ray), such that there are more free input parameters for the user to know about
+and adjust, potentially making the method more complex to use. However, as many
+of the parameters have natural choices, much of this parameterization can be
+handled automatically behind the scenes without the need for the user to be
+aware of this.
+
+.. math::
+    :label: variance_fom
+
+    \text{FOM} = \frac{1}{\text{Time} \times \sigma^2}
+
+Finally, one unique capability of the FW-CADIS weight window generator is to 
+produce weight windows for local variance reduction, given a list of the 
+responses of interest. This is controlled by optionally specifying target 
+tallies from the :class:`openmc.model.Model` to the 
+:class:`openmc.WeightWindowGenerator`, as illustrated in the
+:ref:`user guide<variance_reduction>`. If target tallies for local variance 
+reduction are supplied, then the adjoint sources are only populated after the 
+initial forward simulation in the source regions associated with those tallies. 
+In other regions, the adjoint source term is instead set to zero. The Random  
+Ray solver then determines the adjoint flux map used to generate FW-CADIS 
+weight windows following the usual technique.
+
+.. _methods_source_biasing:
+
+--------------
+Source Biasing
+--------------
+
+In contrast to the previous two methods that introduce population controls
+during transport, source biasing modifies the sampling of the external source
+distribution. The basic premise of the technique is that for each spatial,
+angular, energy, or time distribution of a source, an additional distribution
+can be specified provided that the two share a common support (set of points
+where the distribution is nonzero). Samples are then drawn from this "bias"
+distribution, which can be chosen to preferentially direct particles towards
+phase space regions of interest. In order to avoid biasing the tally results,
+however, a weight adjustment is applied to each sampled site as described below.
+
+Assume that the unbiased probability density function of a random variable
+:math:`X:x \rightarrow \mathbb{R}` is given by :math:`f(x)`, but that using the
+biased distribution :math:`g(x)` will result in a greater number of particle
+trajectories reaching some phase space region of interest. Then a sample
+:math:`x_0` may be drawn from :math:`g(x)` while maintaining a fair game,
+provided that its weight is adjusted as:
+
+.. math::
+    :label: source_bias
+
+    w = w_0 \times \frac{f(x_0)}{g(x_0)}
+
+where :math:`w_0` is the weight of an unbiased sample from :math:`f(x)`,
+typically unity.
+
+Returning now to Equation :eq:`source_bias`, the requirement for common support
+becomes evident. If :math:`\mathrm{supp} (g)` fully contains but is not
+identical to :math:`\mathrm{supp} (f)`, then some samples from :math:`g(x)` will
+correspond to points where :math:`f(x) = 0`. Thus these source sites would be
+assigned a starting weight of 0, meaning the particles would be killed
+immediately upon transport, effectively wasting computation time. Conversely, if
+:math:`\mathrm{supp} (g)` is fully contained by but not identical to
+:math:`\mathrm{supp} (f)`, the contributions of some regions outside
+:math:`\mathrm{supp} (g)` will not be counted towards the integral, potentially
+biasing the tally. The weight assigned to such points would be undefined since
+:math:`g(x) = \mathbf{0}` at these points.
+
+When an independent source is sampled in OpenMC, the particle's coordinate in
+each variable of phase space :math:`(\mathbf{r},\mathbf{\Omega},E,t)` is
+successively drawn from an independent probability distribution. Multiple
+variables can be biased, in which case the resultant weight :math:`w` applied to
+the particle is the product of the weights assigned from all sampled
+distributions: space, angle, energy, and time, as shown in Equation
+:eq:`tot_wgt`.
+
+.. math::
+    :label: tot_wgt
+
+    w = w_r \times w_{\Omega} \times w_E \times w_t
+
+Finally, source biasing and weight windows serve different purposes. Source
+biasing changes how particles are born, allowing the initial source sites to be
+sampled preferentially from important regions of phase space (space, angle,
+energy, and time) with an accompanying weight adjustment. Weight windows, by
+contrast, apply population control during transport (splitting and Russian
+roulette) to help particles reach and contribute in important regions as they
+move through the system. Because particle transport proceeds as usual after a
+biased source is sampled, particle attenuation in optically thick regions
+outside the source volume will not be affected by source biasing; in such
+scenarios, transport biasing techniques such as weight windows are often more
+effective.

docs/source/publications.rst

@@ -138,8 +138,8 @@ Geometry and Visualization
   *Trans. Am. Nucl. Soc.*, **114**, 391-394 (2016).
 
 - Derek M. Lax, "`Memory efficient indexing algorithm for physical properties in
-  OpenMC <http://hdl.handle.net/1721.1/97862>`_," S. M. Thesis, Massachusetts
-  Institute of Technology (2015).
+  OpenMC <https://dspace.mit.edu/handle/1721.1/97862>`_," S. M. Thesis,
+  Massachusetts Institute of Technology (2015).
 
 - Derek Lax, William Boyd, Nicholas Horelik, Benoit Forget, and Kord Smith, "A
   memory efficient algorithm for classifying unique regions in constructive
@@ -399,7 +399,8 @@ Doppler Broadening
 - Jonathan A. Walsh, Benoit Forget, Kord S. Smith, Brian C. Kiedrowski, and
   Forrest B. Brown, "`Direct, on-the-fly calculation of unresolved resonance
   region cross sections in Monte Carlo simulations
-  <http://hdl.handle.net/1721.1/108644>`_," *Proc. Joint Int. Conf. M&C+SNA+MC*,
+  <https://dspace.mit.edu/handle/1721.1/108644>`_," *Proc. Joint Int. Conf.
+  M&C+SNA+MC*,
   Nashville, Tennessee, Apr. 19--23 (2015).
 
 - Colin Josey, Benoit Forget, and Kord Smith, "`Windowed multipole sensitivity
@@ -596,7 +597,8 @@ Depletion
 
 - Matthew S. Ellis, Colin Josey, Benoit Forget, and Kord Smith, "`Spatially
   Continuous Depletion Algorithm for Monte Carlo Simulations
-  <http://hdl.handle.net/1721.1/107880>`_," *Trans. Am. Nucl. Soc.*, **115**,
+  <https://dspace.mit.edu/handle/1721.1/107880>`_," *Trans. Am. Nucl. Soc.*,
+  **115**,
   1221-1224 (2016).
 
 - Anas Gul, K. S. Chaudri, R. Khan, and M. Azeen, "`Development and verification

docs/source/pythonapi/base.rst

@@ -25,6 +25,8 @@ Simulation Settings
    openmc.IndependentSource
    openmc.FileSource
    openmc.CompiledSource
+   openmc.MeshSource
+   openmc.TokamakSource
    openmc.SourceParticle
    openmc.VolumeCalculation
    openmc.Settings
@@ -34,8 +36,8 @@ Simulation Settings
    :nosignatures:
    :template: myfunction.rst
 
+   openmc.read_source_file
    openmc.write_source_file
-   openmc.wwinp_to_wws
 
 Material Specification
 ----------------------
@@ -118,6 +120,7 @@ Constructing Tallies
    openmc.Filter
    openmc.UniverseFilter
    openmc.MaterialFilter
+   openmc.MaterialFromFilter
    openmc.CellFilter
    openmc.CellFromFilter
    openmc.CellBornFilter
@@ -125,10 +128,14 @@ Constructing Tallies
    openmc.CollisionFilter
    openmc.SurfaceFilter
    openmc.MeshFilter
+   openmc.MeshBornFilter
+   openmc.MeshMaterialFilter
    openmc.MeshSurfaceFilter
    openmc.EnergyFilter
    openmc.EnergyoutFilter
+   openmc.ParticleProductionFilter
    openmc.MuFilter
+   openmc.MuSurfaceFilter
    openmc.PolarFilter
    openmc.AzimuthalFilter
    openmc.DistribcellFilter
@@ -138,18 +145,31 @@ Constructing Tallies
    openmc.SpatialLegendreFilter
    openmc.SphericalHarmonicsFilter
    openmc.TimeFilter
+   openmc.WeightFilter
    openmc.ZernikeFilter
    openmc.ZernikeRadialFilter
+   openmc.ParentNuclideFilter
    openmc.ParticleFilter
+   openmc.ReactionFilter
+   openmc.MeshMaterialVolumes
+   openmc.Trigger
+   openmc.TallyDerivative
+   openmc.Tally
+   openmc.Tallies
+
+Meshes
+------
+
+.. autosummary::
+   :toctree: generated
+   :nosignatures:
+   :template: myclassinherit.rst
+
    openmc.RegularMesh
    openmc.RectilinearMesh
    openmc.CylindricalMesh
    openmc.SphericalMesh
    openmc.UnstructuredMesh
-   openmc.Trigger
-   openmc.TallyDerivative
-   openmc.Tally
-   openmc.Tallies
 
 Geometry Plotting
 -----------------
@@ -159,8 +179,10 @@ Geometry Plotting
    :nosignatures:
    :template: myclass.rst
 
-   openmc.Plot
-   openmc.ProjectionPlot
+   openmc.SlicePlot
+   openmc.VoxelPlot
+   openmc.WireframeRayTracePlot
+   openmc.SolidRayTracePlot
    openmc.Plots
 
 Running OpenMC
@@ -186,6 +208,7 @@ Post-processing
    :template: myclass.rst
 
    openmc.Particle
+   openmc.ParticleList
    openmc.ParticleTrack
    openmc.StatePoint
    openmc.Summary
@@ -197,6 +220,9 @@ Post-processing
    :nosignatures:
    :template: myfunction.rst
 
+   openmc.read_collision_track_file
+   openmc.read_collision_track_hdf5
+   openmc.read_collision_track_mcpl
    openmc.voxel_to_vtk
 
 The following classes and functions are used for functional expansion reconstruction.
@@ -239,8 +265,16 @@ Variance Reduction
    :template: myclass
 
    openmc.WeightWindows
+   openmc.WeightWindowsList
    openmc.WeightWindowGenerator
+
+.. autosummary::
+   :toctree: generated
+   :nosignatures:
+   :template: myfunction.rst
+
    openmc.hdf5_to_wws
+   openmc.wwinp_to_wws
 
 
 Coarse Mesh Finite Difference Acceleration

docs/source/pythonapi/capi.rst

@@ -16,15 +16,16 @@ Functions
    current_batch
    export_properties
    export_weight_windows
-   import_weight_windows
    finalize
    find_cell
    find_material
+   dagmc_universe_cell_ids
    global_bounding_box
    global_tallies
    hard_reset
    id_map
    import_properties
+   import_weight_windows
    init
    is_statepoint_batch
    iter_batches
@@ -39,10 +40,12 @@ Functions
    reset_timers
    run
    run_in_memory
+   run_random_ray
    sample_external_source
-   simulation_init
    simulation_finalize
+   simulation_init
    source_bank
+   statepoint_load
    statepoint_write
 
 Classes
@@ -53,12 +56,98 @@ Classes
    :nosignatures:
    :template: myclass.rst
 
+   AzimuthalFilter
    Cell
+   CellFilter
+   CellInstanceFilter
+   CellbornFilter
+   CellfromFilter
+   CollisionFilter
+   CylindricalMesh
+   DelayedGroupFilter
+   DistribcellFilter
    EnergyFilter
-   MaterialFilter
+   EnergyFunctionFilter
+   EnergyoutFilter
+   Filter
+   LegendreFilter
    Material
+   MaterialFilter
+   MaterialFromFilter
+   Mesh
    MeshFilter
+   MeshBornFilter
    MeshSurfaceFilter
+   MuFilter
    Nuclide
+   ParentNuclideFilter
+   ParticleFilter
+   ParticleProductionFilter
+   PolarFilter
+   ReactionFilter
+   RectilinearMesh
    RegularMesh
+   SpatialLegendreFilter
+   SphericalHarmonicsFilter
+   SphericalMesh
+   SolidRayTracePlot
+   SurfaceFilter
    Tally
+   TemporarySession
+   UniverseFilter
+   UnstructuredMesh
+   WeightFilter
+   WeightWindows
+   ZernikeFilter
+   ZernikeRadialFilter
+
+Data
+----
+
+.. data:: cells
+
+   Mapping of cell ID to :class:`openmc.lib.Cell` instances.
+
+   :type: dict
+
+.. data:: filters
+
+   Mapping of filter ID to :class:`openmc.lib.Filter` instances.
+
+   :type: dict
+
+.. data:: materials
+
+   Mapping of material ID to :class:`openmc.lib.Material` instances.
+
+   :type: dict
+
+.. data:: meshes
+
+   Mapping of mesh ID to :class:`openmc.lib.Mesh` instances.
+
+   :type: dict
+
+.. data:: plots
+
+   Mapping of plot ID to :class:`openmc.lib.SolidRayTracePlot` instances.
+
+   :type: dict
+
+.. data:: nuclides
+
+   Mapping of nuclide name to :class:`openmc.lib.Nuclide` instances.
+
+   :type: dict
+
+.. data:: tallies
+
+   Mapping of tally ID to :class:`openmc.lib.Tally` instances.
+
+   :type: dict
+
+.. data:: weight_windows
+
+   Mapping of weight window ID to :class:`openmc.lib.WeightWindows` instances.
+
+   :type: dict