(design_internal_page)= # Design notes (Internal) This area contains assortment of detailed design notes on various topics for the internal development team. (design_CE_page)= ## Constraint Embedding (design_aggsg_page)= ### Aggregation Subtree Decomposition #### Overview Many multibody systems have an underlying rooted-tree topology together with constraints that refer to multiple tree nodes. In this document, each such constraint is called a *cutjoint*. A cutjoint is represented by an unordered list of node IDs rather than by one additional graph edge. The Aggregation Subtree Decomposition (ASD) identifies tree regions that can be treated as compound or *supernodes*. The concrete `Karana::Core::BasicTreeImpl::aggregationTreeASD()` helper performs the tree-related part of this decomposition for one target cutjoint. It does not construct the complete quotient graph. The central construction is: 1. Find the minimal spanning tree of a cutjoint's nodes. 1. Treat that spanning tree's root as the cutjoint's lowest common ancestor (LCA). 1. Exclude the LCA from the candidate supernode, even when the LCA was explicitly listed in the cutjoint. 1. Resolve proper partial overlaps with candidates from the complete cutjoint list. 1. Return an ordinary directed tree whose root is the final LCA and whose remaining nodes define the supernode. #### Scope This document describes the behavior of: ```cpp // Find the ASD subtree associated with one target cutjoint. std::unique_ptr> aggregationTreeASD( const std::vector& target_cutjoint, const std::vector>& cutjoints) const; ``` The method belongs only to the concrete `BasicTreeImpl` class. It is not part of the abstract `Tree` interface. The method returns one target-associated aggregation tree. A higher-level ASD procedure may call this operation while constructing all supernodes and their nesting relationships. #### Inputs ##### Original tree Let $$ T = (V, E_T) $$ be a rooted directed tree. Tree edges point from parent to child. ##### Cutjoint A cutjoint is a finite node set $$ C_i \subseteq V $$ containing at least two distinct nodes. The C++ API represents this set as a `std::vector`. The vector representation has set semantics: - Node order is irrelevant. - Repeated node IDs are ignored. - Every referenced node must belong to the original tree. - At least two distinct nodes must remain after duplicates are removed. For example, `{5, 7, 10}`, `{10, 5, 7}`, and `{5, 7, 10, 5}` represent the same cutjoint. The second method argument contains the complete cutjoint collection $$ \mathcal{C} = \{C_1, C_2, \ldots, C_n\}. $$ The target cutjoint must belong to this collection under node-set equality. #### Candidate Construction For a cutjoint $C_i$, let $$ H_i = \operatorname{span}_T(C_i) $$ be the minimal subtree of $T$ that contains every node in $C_i$. The root of this spanning tree is the common LCA of the cutjoint nodes: $$ a_i = \operatorname{root}(H_i) = \operatorname{LCA}_T(C_i). $$ The candidate supernode is $$ S_i = V(H_i) \setminus \{a_i\}. $$ This rule is unconditional: $a_i$ is excluded even if $a_i \in C_i$. Consequently, the explicitly listed cutjoint nodes are not required to be a subset of the candidate supernode. The precise spanning-tree requirement is therefore: $$ S_i = V(H_i) \setminus \{a_i\}. $$ The figures below use red for an excluded LCA, blue for an explicitly listed node that belongs to the candidate, and green for an included connector node. Dashed green boundaries identify complete candidate supernodes. ##### Figure 1: Multi-node cutjoint Consider the cutjoint $C = \{5, 10, 11\}$. Nodes labeled `listed` were explicitly listed in the cutjoint. ```{graphviz} digraph multi_node_cutjoint { // Draw the original parent-to-child direction from top to bottom. graph [bgcolor="transparent", rankdir=TB, nodesep=0.45, ranksep=0.55]; node [shape=circle, style=filled, fontname="Helvetica"]; edge [color="#4a4a4a"]; lca [label="3\nLCA\nexcluded", fillcolor="#f4cccc"]; // The dashed cluster is the LCA-excluded candidate supernode. subgraph cluster_candidate { label="candidate supernode {5, 7, 10, 11}"; color="#38761d"; fontcolor="#274e13"; style="dashed,rounded"; n5 [label="5\nlisted", fillcolor="#cfe2f3"]; n7 [label="7\nconnector", fillcolor="#d9ead3"]; n10 [label="10\nlisted", fillcolor="#cfe2f3"]; n11 [label="11\nlisted", fillcolor="#cfe2f3"]; } lca -> n5; lca -> n7; lca -> n11; n7 -> n10; } ``` Node `7` is included because it connects listed node `10` to the spanning-tree root. Node `3` is excluded because it is the LCA and becomes the returned-tree root. ##### Figure 2: A listed cutjoint node is the LCA Now consider $C = \{3, 10\}$, where `3` is an ancestor of `10`. ```{graphviz} digraph listed_lca { // Node 3 is both listed and the spanning-tree root, so it is excluded. graph [bgcolor="transparent", rankdir=TB, ranksep=0.6]; node [shape=circle, style=filled, fontname="Helvetica"]; edge [color="#4a4a4a"]; lca [label="3\nlisted LCA\nexcluded", fillcolor="#f4cccc"]; subgraph cluster_candidate { label="candidate supernode {7, 10}"; color="#38761d"; fontcolor="#274e13"; style="dashed,rounded"; n7 [label="7\nconnector", fillcolor="#d9ead3"]; n10 [label="10\nlisted", fillcolor="#cfe2f3"]; } lca -> n7; n7 -> n10; } ``` Although node `3` was explicitly listed, it is still excluded. This is not an exception or an error: cutjoint membership identifies the spanning tree, while the LCA-exclusion rule defines the candidate supernode. The same rule applies when the LCA is the original tree root. The original root is allowed as the returned aggregation-tree root and remains outside the supernode. #### Candidate Representation Conceptually, each candidate stores two pieces of information: $$ K_i = (S_i, a_i), $$ where: - $S_i$ is the LCA-excluded candidate node set. - $a_i$ is the excluded candidate LCA. Keeping the LCA as metadata is necessary because absorbed candidates may have different LCAs. Their LCAs determine the root of the final returned tree. #### Selecting the Target Candidate The first argument identifies one target cutjoint $C_t$. Before candidate construction, the method: 1. Normalizes every cutjoint by removing repeated nodes. 1. Validates all referenced nodes. 1. Verifies that every normalized cutjoint contains at least two nodes. 1. Finds a set-equal copy of $C_t$ in the complete collection $\mathcal{C}$. If the target does not occur in the complete collection, the method rejects the request. The initial aggregation region is the target candidate: $$ S \leftarrow S_t. $$ Its initial LCA collection is $$ A \leftarrow \{a_t\}. $$ #### Overlap Resolution ##### Proper partial overlap A pool candidate $S_i$ merges with the current region $S$ only when the two sets overlap and neither contains the other: $$ S \cap S_i \neq \varnothing, $$ $$ S \nsubseteq S_i, $$ and $$ S_i \nsubseteq S. $$ This relationship is called *proper partial overlap*. When it occurs, the method updates $$ S \leftarrow S \cup S_i $$ and records $a_i$ in the component LCA collection $A$. ##### Containment represents nesting Containment alone does not cause a merge. If $$ S_i \subseteq S $$ or $$ S \subseteq S_i, $$ the candidates can represent nested supernodes and remain separate. Equal candidates also require no merge. ##### Shared LCAs do not imply overlap Candidate LCAs are excluded before overlap testing. Two candidates that share only an LCA therefore remain separate. For example, cutjoints `{3, 10}` and `{3, 5}` have the same LCA `3`, but their candidates `{7, 10}` and `{5}` are disjoint. ##### Figure 3: Partial overlap versus containment Proper partial overlap merges candidate $A = \{10, 12\}$ and candidate $B = \{10, 14\}$ at node `10`: ```{graphviz} digraph partial_overlap { // Orange marks the node shared by both proper-partial-overlap candidates. graph [ bgcolor="transparent", rankdir=TB, ranksep=0.6, label="proper partial overlap: A = {10, 12}, B = {10, 14}", labelloc=t ]; node [shape=circle, style=filled, fontname="Helvetica"]; edge [color="#4a4a4a"]; lca [label="7\nfinal LCA\nexcluded", fillcolor="#f4cccc"]; subgraph cluster_merged { label="merged supernode {10, 12, 14}"; color="#38761d"; fontcolor="#274e13"; style="dashed,rounded"; n10 [label="10\nA and B", fillcolor="#fce5cd"]; n12 [label="12\nA only", fillcolor="#cfe2f3"]; n14 [label="14\nB only", fillcolor="#d9ead3"]; } lca -> n10; lca -> n12; n10 -> n14; } ``` Containment remains nested: ```{graphviz} digraph candidate_containment { // The inner candidate remains a distinct nested region inside the outer one. graph [bgcolor="transparent", rankdir=TB, ranksep=0.65]; node [shape=circle, style=filled, fontname="Helvetica"]; edge [color="#4a4a4a"]; outer_lca [label="3\nouter LCA\nexcluded", fillcolor="#f4cccc"]; subgraph cluster_outer { label="outer candidate {5, 7, 10}"; color="#3d85c6"; fontcolor="#073763"; style="dashed,rounded"; n5 [label="5\nouter only", fillcolor="#cfe2f3"]; n7 [label="7\ninner LCA\nouter member", fillcolor="#cfe2f3"]; subgraph cluster_inner { label="inner candidate {10}"; color="#674ea7"; fontcolor="#351c75"; style="dashed,rounded"; n10 [label="10\ninner member", fillcolor="#d9d2e9"]; } } outer_lca -> n5; outer_lca -> n7; n7 -> n10; } ``` Selecting the inner candidate returns supernode `{10}`. Selecting the outer candidate returns supernode `{5, 7, 10}`. Containment alone does not combine them. ##### Transitive closure Overlap resolution runs to a fixed point. A merge may add nodes that create a proper partial overlap with a candidate examined earlier. The method therefore rescans the candidate pool until a complete pass makes no change. #### Final LCA Suppose the target candidate and all candidates absorbed by proper partial overlap have LCA metadata $$ A = \{a_{i_1}, a_{i_2}, \ldots, a_{i_k}\}. $$ The final returned-tree root is $$ a = \operatorname{LCA}_T(A). $$ A merge may therefore move the returned LCA upward. The LCA is not fixed to the target candidate's original LCA. The final LCA may be the original tree root. This is valid and does not require special rejection. An individual candidate LCA is excluded from that candidate. After merging, it may nevertheless occur inside the final supernode if another candidate's spanning tree contains it. Only the final returned-tree root $a$ is guaranteed to be external to the returned supernode. #### Returned Tree The method returns an ordinary directed tree containing: $$ V_R = \{a\} \cup S, $$ where $a$ is the final LCA and $S$ is the merged supernode node set. The interpretation is: ```text returned_tree.root() final LCA; not in the supernode all other returned-tree vertices nodes in the supernode ``` Original tree edges are reused to connect the returned vertices. Every supernode vertex must have its original parent either inside the supernode or equal to the returned root. Otherwise, the selected union cannot be represented as one valid returned subtree and the method rejects it. ##### Figure 4: Returned-tree interpretation ```{graphviz} digraph returned_tree { // The returned Tree contains both the external LCA and supernode members. graph [bgcolor="transparent", rankdir=TB, nodesep=0.45, ranksep=0.6]; node [shape=circle, style=filled, fontname="Helvetica"]; edge [color="#4a4a4a"]; root [label="a\nreturned root\nexternal LCA", fillcolor="#f4cccc"]; subgraph cluster_supernode { label="semantic supernode S = {x, y, z, w}"; color="#38761d"; fontcolor="#274e13"; style="dashed,rounded"; x [label="x", fillcolor="#d9ead3"]; y [label="y", fillcolor="#d9ead3"]; z [label="z", fillcolor="#d9ead3"]; w [label="w", fillcolor="#d9ead3"]; } root -> x; root -> y; root -> z; y -> w; } ``` Only the nodes inside the dashed boundary define the supernode. The returned root `a` is present in the ordinary tree but excluded from the supernode. Cutjoint relationships are not inserted as edges in the returned tree. The result stores only original directed-tree topology. #### Algorithm ```text function aggregationTreeASD(target_cutjoint, cutjoints): normalized_pool = [] for cutjoint in cutjoints: nodes = validate_and_remove_duplicates(cutjoint) normalized_pool.append(nodes) target = validate_and_remove_duplicates(target_cutjoint) if no node-set-equal target exists in normalized_pool: reject the request selected = candidate_from_spanning_tree(target) candidates = candidate_from_spanning_tree(c) for c in normalized_pool aggregation_nodes = selected.nodes component_lcas = [selected.lca] absorbed = [false for each candidate] changed = true while changed: changed = false for candidate in candidates: if candidate was already absorbed: continue if candidate contains aggregation_nodes: continue // Preserve nesting. if aggregation_nodes contains candidate.nodes: continue // Preserve nesting or equality. if candidate.nodes intersects aggregation_nodes: aggregation_nodes = aggregation_nodes union candidate.nodes component_lcas.append(candidate.lca) mark candidate absorbed changed = true final_lca = common_lca(component_lcas) return original_tree_subtree( root = final_lca, remaining_nodes = aggregation_nodes) ``` #### Required Behaviors ##### LCA exclusion For every independently generated candidate, $$ a_i \notin S_i. $$ This remains true when $a_i$ was explicitly listed in $C_i$. ##### Spanning-tree coverage Every non-root node of a cutjoint's minimal spanning tree belongs to its candidate: $$ V(H_i) \setminus \{a_i\} = S_i. $$ This replaces the obsolete requirement that every explicitly listed cutjoint node must belong to the candidate. ##### Original-root acceptance The original tree root is allowed to be a candidate LCA or the final merged LCA. It is returned as the aggregation-tree root and excluded from the supernode in the usual way. ##### Target membership The target cutjoint must occur in the complete cutjoint collection under node-set equality. ##### Nested containment Candidate equality or containment does not enlarge the target region. This permits an inner supernode to be represented independently inside a later outer supernode. ##### Proper-overlap closure All candidates connected to the selected region through transitive proper partial overlap are absorbed. ##### Deterministic topology Input order and duplicate cutjoint nodes do not affect candidate membership or the returned topology. Returned vertices are assembled in the original tree's topological order. #### Invalid Inputs and Ill-posed Results The method rejects: - A cutjoint containing fewer than two distinct nodes. - A cutjoint referring to a node absent from the original tree. - A target cutjoint absent from the complete cutjoint collection. - A merged node set that cannot be connected through the final LCA using the original parent edges included in the returned tree. The method does not reject a result merely because its LCA is the original tree root. #### Relationship to Full ASD `aggregationTreeASD()` performs the graph-related work needed to identify the candidate associated with one target cutjoint and its proper-overlap component. It does not itself: - Replace original nodes with compound-body nodes. - Construct every nested supernode in one call. - Build the complete quotient tree. - Transfer physical bodies, joints, or constraints into compound bodies. Those operations belong to the higher-level decomposition and model-building procedure. In particular, the phrase *cutjoint absorption* must not be interpreted as a requirement that every explicitly listed cutjoint node lie inside the supernode. ASD uses the cutjoint node set to locate a spanning tree, then excludes that spanning tree's LCA by definition. #### Summary For each cutjoint node set $C_i$: $$ H_i = \operatorname{span}_T(C_i), \qquad a_i = \operatorname{root}(H_i), \qquad S_i = V(H_i) \setminus \{a_i\}. $$ The target cutjoint must belong to the complete cutjoint collection. Starting from its candidate, the method absorbs candidates connected through transitive proper partial overlap while leaving equal or containing candidates separate for nesting. The common LCA of the absorbed candidates becomes the root of the returned tree. That root is always excluded from the represented supernode, including when it is an explicitly listed cutjoint node or the original tree root. (design_SDF_page)= ## SDF importers/exporters ### Converting kdFlex multibodies to SDF This leg projects the multibody data in {py:class}`Karana.Dynamics.SOADyn_types.SubGraphDS` into an SDF 1.9 model. A {py:class}`Karana.KUtils.BasicPrefab.BasicPrefabDS` input is accepted only as a container: the converter selects its `params.subtree`. Prefab and KModel data are outside the conversion boundary. ```{figure} images/kdflex_to_sdf_flow.svg :alt: Forward conversion products. :align: center The forward conversion products and their responsibilities. ``` #### Process The converter loads the source through the concrete kdFlex `fromFile()` API and walks the typed body tree. Each body becomes a link with its mass, center of mass, and inertia tensor. Tree hinges become ordinary SDF joints. Each {py:class}`Karana.Dynamics.SOADyn_types.LoopConstraintCutJointDS` becomes an additional joint marked `kdflex_loop_joint="true"`. kdFlex stores independent parent- and child-body transforms for a connection. SDF ordinarily describes an assembled link pose plus a joint pose. To retain both kdFlex transforms, the exporter creates one frame on each endpoint and chooses the child-link pose so the endpoints coincide: $$ {}^{P}T_C = {}^{P}T_J ({}^{C}T_J)^{-1}. $$ Sensor, force, and constraint nodes also become body-attached SDF frames at the same translation and orientation. Consequently the reverse leg can recover a node without treating a body origin as a joint location. #### Deterministic names Arbitrary kdFlex names are encoded with URL-safe Base64 without padding. The encoding is reversible and stable, including for spaces or punctuation. The four forms are: ```text kdflex__node__ROLES__BODY__NODE kdflex__joint_parent__JOINT__ROLES__NODE kdflex__joint_child__JOINT__ROLES__NODE kdflex__collision__KIND__BODY__NAME__LAYERS__CONVEX ``` Node names encode roles, body, and node. Endpoint names encode the joint and endpoint node. Collision names encode source kind, body, object name, hexadecimal layer mask, and convexity. `ROLES` is a sorted compact token: `c` means constraint, `f` force, `s` sensor, and `n` no node role. `KIND` is `part` or `file`. The exporter also writes `MODEL.sdf.kdflex.toml`. Its `cut_joints`, `joint_endpoints`, and `node_frames` entries allow the reverse leg to interpret an SDF whose deterministic names have been edited. The generated SDF remains self-describing without this optional file. #### Collision geometry Only kdFlex scene objects whose layer mask contains `Karana.Scene.LAYER_COLLISION` become SDF collision elements. The exact supported subset is: | kdFlex source | SDF target | Important restriction | |---|---|---| | `BoxGeometry` | `box` | dimensions include the three-axis scale | | `SphereGeometry` | `sphere` | scale must be uniform | | `CylinderGeometry` | `cylinder` | the two radial scales must match | | `CapsuleGeometry` | `capsule` | the two radial scales must match | | `StaticMeshGeometry` | `mesh` | a filename is required | | `SceneFileObjectSpec` | `mesh` | its scalar scale is repeated on three axes | Object pose is copied. The deterministic name carries the original layer mask, object name, source kind, and mesh `is_convex` flag because SDF collision geometry has no equivalent kdFlex metadata fields. Material and contact behavior are not transferred. #### Issue log Every invocation creates JSON at `OUTPUT.errors.log`, including successful conversions with an empty `entries` array. `--error-log PATH` selects another location. Each entry has `severity`, a stable `code`, a human-readable `message`, and object-specific `context`. Unsupported hinges are exported as fixed and logged. Unsupported constraints and inexact collision shapes are logged; collision shapes that cannot be represented exactly are skipped. Review this file before accepting a conversion. A zero process exit status means the conversion completed, not necessarily that the issue count is zero. #### Gaps and caveats - Unsupported kdFlex hinge families become fixed SDF joints and are logged. - Only position bounds are transferred; joint dynamics, coupling, prescribed state, and coordinate constraints are not. - General constraint frames, `inb_nodal_matrix`, and virtual-root scene data are not represented. - Collision contact, friction, filtering, material, and self-collision settings are not transferred. - Visual mesh URIs are copied as written and may not remain portable. - SI values replace the source unit spelling and numeric formatting. - Base-body attachment to the kdFlex virtual root is implicit in SDF. See [the detailed forward gap register](../kdflex2sdf_GAPS.md) and the [reverse guide](sdf_to_kdflex.md). #### Command ```console python3 kdflex_sdf_datastruct_converter.py to-sdf \ input.yaml model.sdf --error-log model.forward.errors.log ``` ### Converting SDF graphs to kdFlex multibodies This leg converts one direct SDF model into {py:class}`Karana.Dynamics.SOADyn_types.SubGraphDS`. The complete link-and-joint graph is considered; the input is not required to be a tree. ```{figure} images/sdf_to_kdflex_flow.svg :alt: SDF graph decomposition into tree and leftover edges. :align: center Graph decomposition and DataStruct construction in the reverse leg. ``` #### Process Links are graph vertices and joints are undirected edges during selection. The converter first honors any joints listed in the optional TOML `cut_joints` array. It then applies deterministic document-order union-find to the remaining edges. Accepted edges form a spanning forest; cycle-producing and policy-cut edges are leftovers. Each forest is oriented into the single-parent hierarchy required by kdFlex. Its root becomes a base body with a `FULL6DOF` virtual-root hinge. Tree edges become body hinges, while every leftover edge becomes a {py:class}`Karana.Dynamics.SOADyn_types.LoopConstraintCutJointDS`. Edge reversal during orientation also swaps endpoint roles, preserving local transforms. The converter constructs a real `SubGraphDS`, saves it with `toFile()`, reloads it, and validates the concrete type. It does not parse or synthesize kdFlex YAML dictionaries. #### Naming and optional TOML Files written by the forward converter carry reversible `kdflex__node__...`, `kdflex__joint_parent__...`, `kdflex__joint_child__...`, and `kdflex__collision__...` names. The reverse leg decodes these names directly, so a sidecar is not essential. If `MODEL.sdf.kdflex.toml` exists, it is loaded automatically. Use `--edge-sidecar FILE` to override that location. Version 2 uses this structure: ```toml format_version = 2 cut_joints = ["steering_loop"] [[joint_endpoints]] joint = "steering_loop" parent_frame = "parent_marker" child_frame = "child_marker" [[node_frames]] frame = "camera_marker" body = "chassis" name = "camera_mount" roles = ["sensor"] ``` This is particularly useful for externally authored or renamed SDF. It cannot create pose information that is absent or resolve arbitrary semantic frame chains. #### Collision geometry Direct link collisions containing box, sphere, cylinder, or capsule geometry become collision-layer `ScenePartSpec` objects. A uniformly scaled mesh becomes a `SceneFileObjectSpec`. Converter-generated names restore the original kdFlex name, layers, source kind, and mesh convexity; ordinary SDF collisions receive collision-only defaults. Plane, heightmap, ellipsoid, polyline, and other unsupported shapes are skipped and logged. A mesh with unequal scale components is skipped because the target DataStruct has only scalar scale. Collision `surface` data is not represented by these scene DataStructs and generates an omission entry. #### Issue log Every invocation creates `OUTPUT.errors.log`; `--error-log PATH` overrides it. Entries have stable codes and context for the affected joint, constraint, frame, body, or collision. Fatal exceptions are also recorded as `CONVERSION_FAILED` before being re-raised. Treat the log as part of the conversion result and archive it with the SDF, optional TOML, and kdFlex file. #### Gaps and caveats - Only fixed, revolute/continuous, prismatic, universal/revolute2, and ball joints have mappings. Unsupported types fail conversion and are logged. - Arbitrary SDF pose-frame graphs, nested models, worlds, includes, and joints to `world` are not resolved. - Without endpoint frames, a conventional SDF joint provides only one local attachment pose; the other kdFlex endpoint uses identity. - Tree choice can alter hierarchy and coordinate signs. TOML can force cuts but cannot force a root or force an edge into the tree. - Axis `expressed_in`, dynamics, most limits, mimic/gearbox behavior, plugins, and initial state are omitted. - A rotated SDF inertial frame is not applied to the imported inertia tensor. - Visuals and SDF sensor elements do not become kdFlex scene objects or nodes. - General constraint frames, virtual-root scene data, and `inb_nodal_matrix` are not constructed. The regression checks all supported body and node poses and compares the full 6×6 overall multibody spatial-inertia matrix through the kdFlex runtime API. Those checks establish preservation for the covered model; they do not make an unsupported logged feature lossless. See [the detailed reverse gap register](../sdf2kdflex_GAPS.md) and the [forward guide](kdflex_to_sdf.md). #### Command ```console python3 kdflex_sdf_datastruct_converter.py to-kdflex \ model.sdf output.json --edge-sidecar model.sdf.kdflex.toml \ --error-log output.reverse.errors.log ``` ## Using Atmosphere and Aerodynamics Models Together Written by Codex for jain on 2026-08-06. ### Overall goal An atmospheric-flight simulation must answer two different questions at every integration step: 1. **What air surrounds the vehicle here and now?** The atmosphere model supplies density, wind, and other environmental properties. 1. **What loads does that air create on this vehicle?** The aerodynamics model combines the atmospheric state with vehicle motion, geometry, and aerodynamic coefficients to produce force and moment. Keeping these responsibilities separate lets the same vehicle aerodynamics run on Earth, Mars, or another atmospheric body. It also lets a developer replace a simple standard atmosphere with a high-fidelity model without rewriting the vehicle equations. ```mermaid flowchart LR S[Simulation time and vehicle state] F[Planet-fixed atmosphere frame] P[AtmosphereProvider adapter] C[Atmosphere modeling core] A[AerodynamicsDamped KModel] L[Force and moment on vehicle nodes] S --> A F --> A A -->|time and Cartesian position| P P -->|planet-specific coordinates| C C -->|density, wind, temperature,
pressure, speed of sound| P P -->|AtmosphereState| A A --> L ``` The important boundary is the provider adapter. The aerodynamic model knows which frame contains the atmosphere, but it does not know the planet's radius, ellipsoid, geodetic convention, rotation rate, or atmosphere product. ### The aerodynamics model The C++ `Karana::Models::AerodynamicsDamped` KModel is implemented in the neighboring [header](include/AerodynamicsDamped.h). Its standalone `AerodynamicsDampedCore` computes quasi-steady lift, drag, side force, static moments, and angular-rate damping. The thin KModel wrapper obtains live kdFlex frame kinematics and applies a `Karana::Math::SpatialForce` to each configured force node. For each aerodynamic surface, the simulation provides: - A kdFlex force node whose local axes are `+X` forward, `+Y` right, and `+Z` down. - Reference area, chord, span, and aerodynamic coefficients. - A planet- or atmosphere-fixed frame. - Either a unified atmosphere provider or simpler density and wind callbacks. The KModel determines the surface velocity relative to the atmosphere frame. Consequently, a rotating body-fixed frame automatically accounts for the planet's rotation. The simulation must not add the surface velocity caused by planetary rotation a second time. The principal load relationship is $$ \bar q = \frac{1}{2}\rho\lVert\mathbf V_{air}\rVert^2, $$ where $\rho$ comes from the atmosphere provider and $\mathbf V_{air}$ is the surface velocity relative to the air. The model multiplies dynamic pressure by the configured reference geometry and coefficients to obtain force and moment. See [README.md](README.md) for all coefficient equations, axis conventions, physics references, build instructions, and tests. ### The general-purpose atmosphere provider `AerodynamicsDamped::AtmosphereProvider` is the preferred integration point for a coupled atmosphere model. It has the conceptual C++ signature: ```cpp AtmosphereState provider(double time_s, const km::Vec3 &position_m); ``` The arguments have deliberately planet-neutral meanings: - `time_s` is elapsed simulation time in seconds. - `position_m` is the surface position in the axes of the atmosphere frame supplied when the KModel was created. The returned `AtmosphereState` contains: | Field | Meaning | Units used by the C++ interface | |---|---|---| | `density_kg_m3` | Atmospheric mass density | kg/m³ | | `wind_velocity_in_atmosphere` | Wind resolved in atmosphere-frame axes | m/s | | `temperature_k` | Temperature, when available | K | | `pressure_pa` | Static pressure, when available | Pa | | `speed_of_sound_m_s` | Local speed of sound, when available | m/s | Aerodynamics currently consumes density and wind. Temperature, pressure, and speed of sound are carried for diagnostics and future compressibility models. Unavailable optional scalar fields default to `NaN`. The provider is called once per surface during each derivative evaluation. It therefore produces one internally consistent sample rather than asking a stateful atmosphere separately for density and wind. #### Simpler density and wind providers For an analytic or constant environment, separate callbacks remain available: ```text density_provider(time_s, position_m) -> density_kg_m3 wind_provider(time_s, position_m) -> wind_velocity_in_atmosphere_m_s ``` Do not supply these callbacks together with `AtmosphereProvider`; the KModel rejects that ambiguous configuration. With no density callback, the model uses `params.air_density_kg_m3`. With no wind callback, the air is stationary in the configured atmosphere frame. ```mermaid sequenceDiagram participant I as kdFlex integrator participant K as AerodynamicsDamped participant P as AtmosphereProvider participant C as Planet atmosphere core participant N as Vehicle force node I->>K: preDeriv(time, state) K->>K: Read surface pose and relative velocity K->>P: sample(time_s, position in atmosphere frame) P->>C: Convert coordinates and evaluate model C-->>P: Native atmospheric state P-->>K: AtmosphereState K->>K: Compute dynamic pressure and aerodynamic loads K->>N: Accumulate SpatialForce ``` ### The atmospheric modeling core An atmospheric core contains the actual environmental physics or empirical data model. It should be usable without constructing a kdFlex multibody or KModel. Typical inputs include: - A model epoch and elapsed time. - Planet-specific longitude, latitude, and altitude. - Solar and geomagnetic activity, season, or weather data. - Model configuration and data-file locations. Typical outputs are density, temperature, pressure, speed of sound, and local wind. A core may use geodetic coordinates, planetocentric coordinates, areoid height, or another native convention. Those details belong behind the adapter. This separation is useful for testing. The atmosphere core can be checked against its publisher's reference cases, while the provider can be tested for coordinate order, frame direction, units, and wind-vector transformation. ### Adapting the Earth atmosphere The standalone `Karana::EarthGram::Atmosphere` core in `/workspace/code/EarthGRAM` wraps NASA Earth-GRAM. It accepts geodetic longitude, latitude, altitude, and elapsed time. It does not know the vehicle's kdFlex frames. An Earth adapter therefore performs these steps: 1. Receive elapsed time and Cartesian position in an Earth-fixed frame such as ITRF93. 1. Use `TerrainGrid.createEarthWGS84()` and `TerrainGrid.fromCartesianVertex()` to obtain **longitude, latitude, and altitude, in that order**. 1. Populate an Earth-GRAM `GeodeticPoint` and sample the standalone core. 1. Copy density, temperature, pressure, and speed of sound into `AtmosphereState`. 1. If wind is enabled, transform Earth-GRAM's local east-north-up wind into the Cartesian axes of the Earth-fixed atmosphere frame before returning it. ```mermaid flowchart TD E[ITRF93 Cartesian position] W[WGS-84 TerrainGrid] G[Longitude, latitude, altitude] EG[Earth-GRAM core] ENU[Earth-GRAM local ENU wind] EC[Earth-fixed Cartesian wind] AS[AtmosphereState] E --> W --> G --> EG EG --> AS EG --> ENU -->|local-to-Earth-fixed transform| EC --> AS ``` The NESC orbital implementation provides a working density-and-temperature example in [orbital_case.py](../../sims/NESC_2015/orbital/NESC2015checkcases_common/orbital_case.py). Its Earth-specific WGS-84 conversion is intentionally outside AerodynamicsDamped. ### Earth-GRAM runtime and deployment dependencies An atmosphere adapter can be correct in code and still fail at startup if its data products or kernels are not deployed. Treat these files as part of the simulation configuration, record their versions, and check them before a long run. Earth-GRAM's `Config.data_root` names the directory that contains the `modeldata` and `NCEPdata` subdirectories. The workspace default is: ```text /workspace/DATA/GramSuite/Earth/data ``` With that default, the current NCEP/MET configuration expects this layout: ```text /workspace/DATA/GramSuite/Earth/data/ ├── modeldata/ │ ├── sdata.txt │ ├── topo.txt │ └── zdata.txt └── NCEPdata/ └── FixedBin/ ├── Nb971501.bin ├── ... └── Nb971512.bin ``` The three `modeldata` files are runtime dependencies, not examples: | File | Earth-GRAM purpose | |---|---| | `sdata.txt` | Stationary perturbation data | | `topo.txt` | Topographical data | | `zdata.txt` | Zonal-mean data | The twelve `Nb9715MM.bin` files cover the calendar months for the configured 1997–2015 NCEP climatology. A run opens the binary matching its configured month. Keeping all twelve allows the same installation to run any month. `Config.spice_root` separately names a directory with three required subdirectories. The workspace default is `/workspace/code/EarthGRAM/resources/spice` and resolves to: ```text spice_root/ ├── lsk/naif0012.tls ├── pck/pck00011.tpc └── spice/de442s.bsp ``` The current files may be symbolic links, but every link target must remain readable in the runtime environment. The leap-second kernel establishes time conversion, the text PCK supplies planetary constants, and the BSP supplies the Earth/Sun ephemeris used by Earth-GRAM. | Dependency | Build time | Runtime | Notes | |---|:---:|:---:|---| | GramSuite headers | Yes | No | Compile the standalone wrapper | | Earth/GRAMCommon static libraries | Yes | No | Linked into the wrapper | | CSPICE headers and static library | Yes | No | Linked into the wrapper | | `modeldata` files | No | Yes | All three files are required | | NCEP `FixedBin` data | No | Yes | At least the run's month is required | | Earth-GRAM SPICE kernels | No | Yes | LSK, PCK, and BSP are required | | Earth texture | No | No | Visualization only | The GramSuite `sample_inputs` directory is not a runtime dependency. It is useful only for vendor examples and independent regression comparisons. #### Relocating the runtime data The files do not have to remain under `/workspace`. Set both roots before constructing the stateful Earth-GRAM atmosphere: ```python config.data_root = "/opt/atmosphere/Earth/data" config.spice_root = "/opt/atmosphere/Earth/spice" earth_atmosphere = keg.Atmosphere(config) ``` Earth-GRAM will then expect `modeldata/sdata.txt`, `modeldata/topo.txt`, `modeldata/zdata.txt`, and `NCEPdata/FixedBin` beneath the new `data_root`. It will expect `lsk`, `pck`, and `spice` beneath the new `spice_root`. For reproducible deployment, validate the directories, individual files, permissions, symbolic-link targets, and model versions during simulation startup. Checking only that the root directory exists is not sufficient to prove that every later atmosphere query can load its data. ### Switching to another planetary atmosphere Changing planets should require a new frame and provider adapter, not a new aerodynamics model. For Mars, for example: 1. Supply a Mars body-fixed atmosphere frame to AerodynamicsDamped. 1. Implement a Mars adapter with the same `(time_s, position_m)` contract. 1. Convert Mars-fixed Cartesian position using the ellipsoid, areoid, or planetocentric convention required by the selected Mars atmosphere core. 1. Configure and call Mars-GRAM or another Mars atmosphere implementation. 1. Transform any local north-east-down or east-north-up wind into Mars-fixed Cartesian axes. 1. Return the common `AtmosphereState`. No WGS-84, ITRF93, Earth rotation rate, or Earth-GRAM type should appear in the Mars adapter. Conversely, no Mars datum or Mars-GRAM type should enter AerodynamicsDamped. | Concern | Earth implementation | Mars implementation | Shared component | |---|---|---|---| | Atmosphere frame | ITRF93 | Mars-fixed | Frame passed to KModel | | Coordinates | WGS-84 | Mars datum | Provider adapter pattern | | Atmosphere physics | Earth-GRAM | Mars-GRAM or equivalent | `AtmosphereState` | | Vehicle coefficients | Unchanged | Unchanged | AerodynamicsDamped | | Load application | Unchanged | Unchanged | kdFlex force nodes | ### Basic integration example The following Python excerpt assumes that the simulation has already created a state propagator, a vehicle force node, and an Earth-fixed frame. It shows the important provider boundary without obscuring it with complete vehicle setup. ```python import numpy as np import Karana.Math as km from Karana.Math.Kquantities import ureg from Karana.Scene._Terrain_Py import TerrainGrid import _AerodynamicsDamped_Py as kad import _EarthGRAM_Py as keg terrain = TerrainGrid.createEarthWGS84() earth_atmosphere = keg.Atmosphere(earth_gram_config) def atmosphereProvider(time_s: float, position_m: np.ndarray) -> kad.AtmosphereState: """Convert an Earth-fixed point and return one coherent atmosphere sample.""" # TerrainGrid returns longitude first, followed by latitude and altitude. longitude, latitude, altitude = terrain.fromCartesianVertex(position_m) point = keg.GeodeticPoint() point.elapsed_time = time_s point.longitude = float(longitude) point.latitude = float(latitude) point.altitude = float(altitude) earth_state = earth_atmosphere.sample(point) state = kad.AtmosphereState() state.density_kg_m3 = earth_state.density state.temperature_k = earth_state.temperature state.pressure_pa = earth_state.pressure state.speed_of_sound_m_s = earth_state.speed_of_sound # This minimal example chooses still air in the rotating Earth-fixed frame. # A wind-enabled adapter must transform Earth-GRAM ENU wind to Earth-fixed # Cartesian axes before assigning this field. state.wind_velocity_in_atmosphere = np.zeros(3) * ureg.meter / ureg.second return state aerodynamics = kad.AerodynamicsDamped.create( "vehicle_aerodynamics", propagator, [aerodynamic_force_node], earth_fixed_frame, None, # No separate density provider. None, # No separate wind provider. None, # Use the built-in coefficient model, not a custom load provider. atmosphereProvider, ) surface = kad.AeroSurfaceConfig() surface.name = "vehicle_reference_surface" surface.area_m2 = 12.0 surface.chord_m = 2.0 surface.span_m = 8.0 surface.cl0 = 0.2 surface.cl_alpha = 4.5 surface.cd0 = 0.03 surface.induced_drag_factor = 0.05 aerodynamics.params.surfaces = [surface] # Ask kdFlex to validate the completed model graph after configuration. propagator.ensureHealthy() ``` The C++ KModel automatically converts compatible Python quantities at its quantity-aware properties. The plain numeric fields in `AtmosphereState` use the SI units listed earlier. ### Summary recipe Use this checklist when adding atmospheric aerodynamics to a vehicle: 1. **Build vehicle dynamics.** Create the multibody, propagator, physical body, and one force node for each aerodynamic surface. 1. **Choose the atmosphere frame.** Use the rotating body-fixed frame expected by the environmental model. 1. **Configure the atmosphere core.** Set its epoch, environmental indices, data paths, and planet-specific options. Verify every required runtime data file and kernel before constructing the atmosphere. 1. **Write a thin provider adapter.** Convert body-fixed Cartesian position to the core's native coordinates, evaluate it once, transform wind into atmosphere-frame Cartesian axes, and return `AtmosphereState`. 1. **Create AerodynamicsDamped.** Pass the propagator, surface nodes, atmosphere frame, and unified provider. 1. **Configure surfaces.** Supply geometry and coefficients in each node's local aerodynamic axes. 1. **Validate conventions.** Test coordinate order, altitude datum, wind sign, frame axes, units, and a zero-wind case before running a full trajectory. 1. **Validate lifecycle and physics.** Run a short simulation and cleanup test, then compare the atmosphere and trajectory with independent reference data. Following this recipe keeps vehicle aerodynamics reusable while placing every planet-specific decision in a small, testable atmosphere adapter. ## NESC Flight Dynamics Check Cases ### NESC 2015 kdFlex Check-Case Overview Written by Codex for jain on 2026-08-06. #### Purpose and scope This workspace ports the NASA Engineering and Safety Center (NESC) 2015 six-degree-of-freedom flight-simulation check cases to kdFlex. The cases are verification problems, not one "correct" reference trajectory. NASA collected independent submissions from several simulation tools so that developers could compare equations of motion, frames, environment models, vehicle models, and integration behavior. The workspace contains two related suites: - `atmospheric`: falling bodies, cannonballs, controlled F-16 maneuvers, and a two-stage launch vehicle. - `orbital`: ISS, sphere, and cylinder cases that progressively add gravity harmonics, third bodies, atmosphere, applied loads, rigid-body rotation, and gravity-gradient torque. Together they contain 45 implemented case directories: 19 atmospheric scenarios or subcases and 26 orbital scenarios or subcases. NASA did not complete Atmospheric Case 14 because of resource constraints. The published orbital sequence begins at Case 02. Those numbering gaps are therefore part of the source assessment, not missing folders in this port. The primary source is the [NASA NESC 2015 check-case website](https://nescacademy.nasa.gov/flightsim/2015), which provides scenario descriptions, DAVE-ML models, initial conditions, raw trajectories, plotting tools, specifications, and errata. The accompanying [Volume I report](https://ntrs.nasa.gov/citations/20150001263) gives the assessment overview, while [Volume II](https://nescacademy.nasa.gov/src/flightsim/Reports/NASA-TM-2015-218675-EOM_checkcase_appendices.pdf) contains detailed definitions and comparisons. ```mermaid flowchart LR N[NASA definitions, models,
errata, and trajectories] C[Typed case configuration] E[Shared kdFlex engine] P[Selected physics branches] O[NASA-shaped CSV output] V[Comparison with independent
NASA submissions] N --> C --> E --> P --> O --> V ``` #### Atmospheric cases ##### Simple-body progression Cases 01 through 10 deliberately add one modeling feature at a time. This makes a discrepancy easier to associate with drag, rotational damping, geodesy, gravity, Earth rotation, or wind. | Case | Scenario | Main verification focus | |---|---|---| | 01 | Dragless dropped sphere | WGS-84/J2 free fall | | 02 | Tumbling brick, no damping | 6-DOF rotation without aero torque | | 03 | Tumbling brick with damping | Angular-rate aero moments | | 04 | Dropped sphere, round Earth | Stationary sphere and point gravity | | 05 | Dropped sphere, rotating Earth | Air-relative velocity from rotation | | 06 | Dropped sphere, WGS-84 | Oblate geodesy and J2 gravity | | 07 | Dropped sphere, steady wind | Constant local-east wind | | 08 | Dropped sphere, wind shear | Altitude-dependent local-east wind | | 09 | Eastward cannonball | Horizontal launch and rotating WGS-84 | | 10 | Northward cannonball | Meridian flight and frame conventions | ##### F-16 and launch-vehicle progression | Case | Scenario | Main verification focus | |---|---|---| | 11 | Subsonic F-16 trim | Open-loop low-speed trim flyout | | 12 | Supersonic F-16 trim | High-Mach aero and propulsion tables | | 13.1 | F-16 altitude step | Closed-loop longitudinal response | | 13.2 | F-16 airspeed step | Throttle and speed control | | 13.3 | F-16 heading step | Lateral-directional control | | 13.4 | F-16 lateral step | Navigation and lateral position control | | 15 | F-16 near the North Pole | Polar navigation singularity handling | | 16 | F-16 at the date line | Longitude wrapping and guidance continuity | | 17 | Unguided two-stage rocket | Staging, variable mass, thrust, and ascent | #### Orbital cases The orbital suite starts with an eight-hour ISS propagation and then isolates or combines environmental and rigid-body effects. | Case | Scenario | Main verification focus | |---|---|---| | 02 | ISS, spherical Earth | Kepler propagation and J2000/ITRF93 output | | 03A | ISS, GEM-T1 4×4 | Harmonic coefficient normalization | | 03B | ISS, GEM-T1 8×8 | Higher-degree rotating gravity field | | 04 | ISS with Sun and Moon | Differential third-body gravity | | 05A | Elliptical ISS, minimum MET | Low solar-activity atmosphere output | | 05B | Elliptical ISS, mean MET | Mean solar-activity atmosphere output | | 05C | Elliptical ISS, maximum MET | High solar-activity atmosphere output | | 06A | Sphere with fixed drag | Constant density and Earth-relative speed | | 06B | Sphere with dynamic drag | Live MET density coupled to drag | | 06C | Cylinder plane change | Finite-duration cross-track thrust | | 06D | Cylinder Earth departure | Finite-duration along-track thrust | | 07A | Sphere, 4×4 plus third bodies | Combined gravity without drag | | 07B | Sphere, 8×8 plus third bodies | Higher-degree combined gravity | | 07C | Case 07A plus drag | 4×4 gravity, third bodies, and MET drag | | 07D | Case 07B plus drag | 8×8 gravity, third bodies, and MET drag | | 08A | ISS free rotation, zero rate | Inertially fixed attitude | | 08B | ISS free rotation, initial rate | Torque-free asymmetric rigid body | | 09A | ISS torque, zero rate | Scheduled body torque | | 09B | ISS torque, initial rate | Torque plus existing rotation | | 09C | ISS force and torque, zero rate | Coupled translation and rotation | | 09D | ISS force and torque, initial rate | Coupled loads with body rotation | | 10A | Circular cylinder, zero offset | Gravity-gradient torque | | 10B | Circular cylinder, rate offset | Gravity gradient plus rate offset | | 10C | Elliptical cylinder, zero offset | Time-varying gravity gradient | | 10D | Elliptical cylinder, rate offset | Elliptical case with rate offset | | FULL | Elliptical ISS, all effects | Coupled environment and attitude | #### General kdFlex modeling approach The port uses a common architecture instead of reproducing each case as a large independent script. A small launcher selects an immutable, typed case configuration. A shared engine creates the kdFlex object graph and activates the required physics. This keeps a correction to frames, units, or force application consistent across every affected case. ```mermaid flowchart TD L[Small case launcher] T[Frozen typed configuration] S[Shared simulation engine] R[Reference frames and transforms] B[Multibody, rigid body, and force nodes] M[Gravity, atmosphere, aero,
propulsion, and applied-load models] I[StatePropagator and RK4 integration] CSV[Case CSV writer] L --> T --> S S --> R S --> B S --> M R --> I B --> I M --> I I --> CSV ``` ##### Object graph and lifecycle Each engine creates or reuses a `Sim`, then constructs a `FrameContainer`, `Multibody`, unconstrained `PhysicalBody` objects, and a `StatePropagator`. Force-producing models use explicit kdFlex force nodes. The engines retain mutually referring objects for the whole run and discard the propagator, multibody, and frame container in kdFlex lifecycle order. The launchers also accept kdpy's pre-created `sim` object. The same main module can therefore run under ordinary Python or kdpy without maintaining a second simulation definition. ##### Quantities and transforms Published inputs mix SI and customary units. The scripts attach kdFlex/Pint units to values and let kdFlex's C++ quantity interfaces perform compatible conversion. Output conversion occurs only when populating a NASA-named CSV column. Frame geometry is composed with `HomTran` and `UnitQuaternion` operations. Atmospheric cases use local navigation, Earth-fixed, inertial, and vehicle frames. Orbital cases propagate in J2000 and use SPICE ITRF93 transformations for Earth-fixed output and atmosphere sampling. The NASA [datum and coordinate-system guide](https://nescacademy.nasa.gov/flightsim/2015/datum_coordinate_system) is the controlling web reference for axis and Earth-model conventions. The [NAIF frame documentation](https://naif.jpl.nasa.gov/pub/naif/toolkit_docs/FORTRAN/req/frames.html) supports the J2000 and ITRF93 implementation. ##### Gravity and celestial bodies The `NBodyGravity` model supplies central gravity, spherical harmonics, and differential third-body acceleration. Atmospheric cases select point gravity or WGS-84/J2 through typed configuration. Orbital Cases 03 and 07 supply the published fully normalized GEM-T1 coefficients through degree 4 or 8. SPICE kernels provide Earth orientation and Sun/Moon ephemerides. NASA Case 04 specified DE405; the local port uses the installed NAIF DE442s ephemeris and documents that version difference rather than hiding it. ##### Atmosphere and aerodynamics The atmospheric engines use the U.S. Standard Atmosphere 1976 implementation in the shared support directory. The `AerodynamicsDamped` C++ KModel receives surface kinematics relative to an atmosphere frame and applies spatial forces at configured nodes. Its unified provider keeps planetary atmosphere details outside the vehicle aerodynamics. Orbital MET cases use the standalone Earth-GRAM core. An Earth-specific adapter converts ITRF93 Cartesian position to WGS-84 longitude, latitude, and altitude, then returns one atmosphere state to `AerodynamicsDamped`. The same live density drives drag in Cases 06B, 07C, 07D, and FULL. The [Earth-GRAM 2024 user guide](https://ntrs.nasa.gov/citations/20240008432) documents the underlying atmosphere product. ```mermaid flowchart LR K[Vehicle surface kinematics] AF[Atmosphere frame] AP[Atmosphere adapter] AC[USSA76 or Earth-GRAM core] AD[AerodynamicsDamped] FN[kdFlex force node] K --> AD AF --> AD AD -->|time and Cartesian position| AP AP --> AC --> AP AP -->|density and wind| AD --> FN ``` ##### DAVE-ML vehicle models The F-16 and rocket use NASA-supplied DAVE-ML data rather than replacement hand-fit equations. The shared reader evaluates MathML expressions and multidimensional tables, interpolating within the published domain and clamping where extrapolation is prohibited. For the F-16, separate files provide inertia, aerodynamic coefficients, propulsion, open-loop control, and closed-loop guidance. Case 17 uses distinct rocket inertia, propulsion, and aerodynamic files. DAVE-ML is an exchange format for dynamic models; background and tools are available at the [DAVE-ML website](https://daveml.org/). ##### Integration, events, and output Most cases use fixed-step fourth-order Runge-Kutta integration. Reporting cadence is independent of the internal integration step. Orbital thrust cases locally refine integration near discontinuous burn boundaries so an RK4 stage does not bias total impulse. FULL uses a smaller step for coupled asymmetric attitude and gravity-gradient response. Every result uses NASA-compatible headings and units. Comparison helpers read all supplied independent trajectories, interpolate them onto the kdFlex output times, and report maximum absolute differences for shared fields. The spread among NASA submissions is meaningful evidence, especially for nonlinear atmospheric cases; it is not reduced to a single presumed truth trajectory. #### Case-specific implementation branches ##### Atmospheric Cases 01–03 `drop_sim.py` handles the shared WGS-84/J2 drop and tumble problem. `DropCaseConfig` selects mass properties, initial position and angular rate, and whether aerodynamic damping is active. Case 03 alone creates the damped aerodynamic surface. ##### Atmospheric Cases 04–10 `sphere_sim.py` handles spherical or WGS-84 Earth, point or J2 gravity, rotating or stationary Earth, and still, steady, or sheared wind. A typed `SphereCaseConfig` selects these features without tests on a case-number string inside the physics functions. ##### Atmospheric Cases 11–16 `f16_sim.py` builds the common six-DOF aircraft. Cases 11 and 12 use fixed trim commands. Cases 13.1–13.4 add command schedules and closed-loop control. Cases 15 and 16 select special guidance modes for polar and date-line navigation. ##### Atmospheric Case 17 Case 17 has its own engine because mass, inertia, propulsion, and aerodynamic reference points change during flight. A derivative callback updates inertia as propellant is consumed, and staging removes the empty first stage before the second-stage burn. ##### Orbital Cases 02–FULL `orbital_case.py` owns the common ISS/sphere/cylinder construction, J2000 and ITRF93 frames, gravity model, atmosphere adapter, scheduled loads, and output. `CaseConfig` activates the requested harmonic degree, third bodies, orbit shape, solar activity, drag, thrust, initial body rate, or gravity gradient. ```mermaid flowchart TD O[Orbital common engine] O --> G{Gravity} G --> G0[Spherical] G --> GH[GEM-T1 4×4 or 8×8] G --> G3[Sun and Moon] O --> E{Environment and loads} E --> D[Fixed or Earth-GRAM drag] E --> F[Scheduled force and torque] O --> A{Attitude} A --> FR[Free rotation] A --> GG[Gravity-gradient torque] G0 --> X[Selected case] GH --> X G3 --> X D --> X F --> X FR --> X GG --> X ``` #### Running and validating a case The normal recipe is: 1. Enter the desired case directory. 1. Run its `nesc_atmos_*.py` or `nesc_orbit_*.py` launcher. 1. Use `--duration` and `--integration-step` for a short smoke or sensitivity run when those options are provided. 1. Run the case comparison helper against its retained `reference_data`. 1. Run the case tests and the applicable shared-engine lifecycle tests. For example: ```sh cd atmospheric/NESC2015checkcases_07 python3 nesc_atmos_07.py python3 compare_reference.py python3 -m pytest -q test_nesc_atmos_07.py ``` The shared orbital suite can be checked with: ```sh cd orbital/NESC2015checkcases_common python3 -m pytest -q test_orbital_cases.py ``` Further local implementation detail is available in the [atmospheric shared guide](atmospheric/NESC2015checkcases_common/README.md), the [orbital shared guide](orbital/NESC2015checkcases_common/README.md), and the [atmosphere/aerodynamics tutorial](../../KModels/aerodynamicsDamped/AerodanamicsUsage.md). #### Web references - [NASA NESC 2015 check-case website](https://nescacademy.nasa.gov/flightsim/2015) - [NASA/TM-2015-218675 Volume I](https://ntrs.nasa.gov/citations/20150001263) - [NASA/TM-2015-218675 Volume II](https://nescacademy.nasa.gov/src/flightsim/Reports/NASA-TM-2015-218675-EOM_checkcase_appendices.pdf) - [NASA datum and coordinate systems](https://nescacademy.nasa.gov/flightsim/2015/datum_coordinate_system) - [NASA Earth-GRAM 2024 user guide](https://ntrs.nasa.gov/citations/20240008432) - [NAIF SPICE reference frames](https://naif.jpl.nasa.gov/pub/naif/toolkit_docs/FORTRAN/req/frames.html) - [DAVE-ML model exchange format](https://daveml.org/) Local case README files document exact inputs, known qualifications, commands, and comparison results for each implementation. ### NESC 2023 Lunar kdFlex Check-Case Overview Written by Codex for jain on 2026-08-06. #### Purpose and scope This workspace ports the NASA Engineering and Safety Center (NESC) lunar six-degree-of-freedom simulation check cases to kdFlex. The assessment extends the 2015 Earth-focused work into cislunar and lunar-orbit regimes relevant to Artemis and Human Landing System analysis. The 17 independently runnable scenarios comprise nine numbered cases and eight lettered subcases. They progress from an analytically checkable two-body orbit through detailed lunar gravity, Earth and Sun perturbations, near-rectilinear halo orbits, offset sensor kinematics, and applied body moments. NASA supplies a family of results from eight independently implemented simulation tools. These are cross-comparison references rather than a claim that one participant file is exact truth. Case 01 additionally has a closed-form Kepler reference. The primary source is the [NASA NESC 2023 Lunar Check-Cases website](https://nescacademy.nasa.gov/flightsim/2023). The assessment overview is [NASA/TM-20240013031 Volume I](https://ntrs.nasa.gov/citations/20240013031), and detailed definitions begin in [Volume II, Part 1](https://ntrs.nasa.gov/citations/20240013556). ```mermaid flowchart LR N[NASA cases, initial states,
kernels, and participant data] C[Typed LunarCaseConfig] E[Shared kdFlex lunar engine] P[Selected gravity, body,
sensor, and load branches] O[NASA-shaped CSV output] V[Cross-comparison with
eight participant solutions] N --> C --> E --> P --> O --> V ``` #### Case inventory | Case | Scenario | Duration | Main verification focus | |---|---|---:|---| | [01][cc01] | Keplerian low orbit | 8 h | Point gravity and analytic solution | | [02][cc02] | Low-fidelity GRAIL | 8 h | Degree/order 8 lunar gravity | | [03][cc03] | High-fidelity GRAIL | 8 h | Degree/order 320 lunar gravity | | [04][cc04] | Circular Apollo | 8 h | Apollo inertia and 8×8 gravity | | [05][cc05] | Perturbed circular | 8 h | Differential Earth/Sun gravity | | [05A][cc05a] | Perturbed tumble | 8 h | Third bodies and initial tumble | | [06][cc06] | Elliptical cylinder | 28 h | Long-period perturbations | | [06A][cc06a] | Zero-rate ellipse | 28 h | Inertially fixed initial attitude | | [07][cc07] | Elliptical Apollo | 28 h | Apollo inertia in same orbit | | [08][cc08] | Baseline NRHO | 7 d | Long-duration cislunar propagation | | [08A][cc08a] | NRHO at anomaly 180° | 7 d | Far-orbit initial geometry | | [08B][cc08b] | NRHO at anomaly 0° | 7 d | Near-Moon initial geometry | | [08C][cc08c] | NRHO radius perturbation | 7 d | Initial-position sensitivity | | [08D][cc08d] | NRHO velocity delta | 7 d | Initial-velocity sensitivity | | [09][cc09] | Polar sensor A | 8 h | Sensor and lunar test points | | [09A][cc09a] | Polar orbit, sensor B | 8 h | Different sensor lever arm | | [09B][cc09b] | Sensor B with moments | 8 h | Piecewise applied body moment | [cc01]: https://nescacademy.nasa.gov/flightsim/2023/cc01 [cc02]: https://nescacademy.nasa.gov/flightsim/2023/cc02 [cc03]: https://nescacademy.nasa.gov/flightsim/2023/cc03 [cc04]: https://nescacademy.nasa.gov/flightsim/2023/cc04 [cc05]: https://nescacademy.nasa.gov/flightsim/2023/cc05 [cc05a]: https://nescacademy.nasa.gov/flightsim/2023/cc05a [cc06]: https://nescacademy.nasa.gov/flightsim/2023/cc06 [cc06a]: https://nescacademy.nasa.gov/flightsim/2023/cc06a [cc07]: https://nescacademy.nasa.gov/flightsim/2023/cc07 [cc08]: https://nescacademy.nasa.gov/flightsim/2023/cc08 [cc08a]: https://nescacademy.nasa.gov/flightsim/2023/cc08a [cc08b]: https://nescacademy.nasa.gov/flightsim/2023/cc08b [cc08c]: https://nescacademy.nasa.gov/flightsim/2023/cc08c [cc08d]: https://nescacademy.nasa.gov/flightsim/2023/cc08d [cc09]: https://nescacademy.nasa.gov/flightsim/2023/cc09 [cc09a]: https://nescacademy.nasa.gov/flightsim/2023/cc09a [cc09b]: https://nescacademy.nasa.gov/flightsim/2023/cc09b #### How the cases increase fidelity The early cases isolate lunar gravity fidelity. The middle cases add vehicle inertia and third-body perturbations. The later cases stress long-duration cislunar propagation and local point kinematics. ```mermaid flowchart LR C1[Case 01
point Moon] C2[Case 02
GRAIL 8×8] C3[Case 03
GRAIL 320×320] C4[Cases 04–07
Apollo/cylinder and
Earth/Sun perturbations] C8[Cases 08–08D
seven-day NRHO variants] C9[Cases 09–09B
sensor points and moments] C1 --> C2 --> C3 --> C4 --> C8 --> C9 ``` The progression should not be read as one vehicle gaining every feature. Individual cases deliberately change body type, initial orbit, attitude, or duration to isolate a particular implementation concern. #### General kdFlex modeling approach Cases 03 through 09B use one shared simulation engine and immutable typed configuration. Cases 01 and 02 retain their extensively documented standalone implementations because they establish the point-gravity and low-order-gravity foundations. All three implementations follow the same kdFlex construction and lifecycle conventions. ##### Shared object graph The simulation creates or reuses a `Sim`, then constructs a `FrameContainer`, `Multibody`, one unconstrained lunar vehicle, and a `StatePropagator`. Vehicle mass and inertia select either NASA's cylinder or Apollo values. An explicit force node is added when Case 09B requires the applied moment. ```mermaid flowchart TD L[Case launcher] T[Typed case configuration] S[Shared lunar_case engine] F[J2000 and lunar PA frames] B[Cylinder or Apollo rigid body] G[NBodyGravity and harmonic field] TB[Optional Earth/Sun bodies] SN[Optional sensor node] M[Optional moment force node] I[StatePropagator and RK4] CSV[Result writer] L --> T --> S S --> F S --> B S --> G S --> TB S --> SN S --> M F --> I B --> I G --> I TB --> I SN --> CSV M --> I I --> CSV ``` Each launcher can run with ordinary Python or reuse kdpy's already-created `sim` object. Cleanup discards the propagator, multibody, and frame container in kdFlex lifecycle order. ##### Units and transforms Inputs remain kdFlex/Pint quantities until they cross a quantity-aware C++ interface. No custom magnitude-conversion helpers are needed. Result values are converted only when populating columns with NASA's specified units. Attitudes use `UnitQuaternion`. Frame poses and position-vector transformations use `HomTran` composition. Participant quaternion signs are compared sign-invariantly because $q$ and $-q$ represent the same orientation. ##### Lunar frames and SPICE Vehicle translation is propagated in the Moon-centered inertial frame aligned with J2000. Harmonic coefficients are fixed in the DE440 lunar principal-axis frame, which rotates relative to J2000. kdFlex frame paths supply that live orientation without manually constructing rotation matrices. ```mermaid flowchart LR J[J2000-aligned
Moon inertial frame] PA[DE440 lunar
principal-axis frame] V[Vehicle body frame] S[Offset sensor frame] H[GRAIL harmonic field] EP[Earth and Sun
DE440 positions] J <-->|SPICE lunar orientation| PA PA --> H H -->|gravity in inertial dynamics| J J --> V --> S EP -->|differential acceleration| J ``` The shared kernel set contains the DE440 planetary ephemeris, DE440 lunar orientation binary PCK, and lunar frame kernel. NASA's [output specification](https://nescacademy.nasa.gov/flightsim/2023/output_specification) defines the reported frames and fields. ##### Lunar gravity Case 01 uses NASA's specified DE400 point-mass gravitational parameter so it can be checked directly against the closed-form Kepler solution. The broader suite uses GRGM660PRIM, a GRAIL-derived fully normalized spherical-harmonic gravity model. The local parser loads the official coefficient file and truncates it to the configured degree and order: - Case 02 and Cases 04–09B use degree and order 8. - Case 03 uses degree and order 320. The NASA report identifies the PDS GRAIL spherical-harmonic archive without giving one exact filename. Testing candidate products against the published initial acceleration identified `GGGRX_0660PM_SHA.TAB` as the matching GRGM660PRIM product. The [PDS GRAIL gravity documentation](https://pds-geosciences.wustl.edu/grail/grail-l-lgrs-2-edr-v1/grail_0001/document/shbdr.htm) describes this solution, and the [PDS SHADR directory](https://pds-geosciences.wustl.edu/grail/grail-l-lgrs-5-rdr-v1/grail_1001/shadr/) hosts the coefficient product. ##### Earth and Sun perturbations Cases 05 through 09B add differential Earth and Sun gravity from DE440. Differential gravity is essential: kdFlex subtracts the acceleration that the third body gives the Moon from the acceleration it gives the spacecraft. The remaining term changes the Moon-relative trajectory without incorrectly accelerating the coordinate origin. ##### Vehicle attitude and inertia The cylinder is used for the low-orbit and selected elliptical cases. Other cases use NASA's Apollo mass and full inertia tensor. NASA lists products of inertia using the engineering convention, so their signs are reversed when placed in the mathematical inertia matrix. Initial attitude is either built from the published vehicle-orbit frame or read from a published quaternion. Body angular velocity is carefully distinguished from rate relative to a moving orbit frame. Case 05A adds tumble, while Case 06A deliberately starts at zero inertial angular velocity. ##### Sensor and applied-moment cases Cases 09, 09A, and 09B attach a sensor at a fixed body-frame offset. kdFlex frame kinematics provide sensor position, velocity, and specific acceleration, including rotational lever-arm effects. Case 09A changes only the sensor location. Case 09B also applies NASA's piecewise body moment. At a switching instant, the callback uses the average of the left and right limiting values so classical fixed-step RK4 integrates the discontinuity without an endpoint bias. #### Case-specific implementation groups ##### Cases 01-03: gravity verification ladder Case 01 checks state propagation and output against an analytical point-gravity solution. Case 02 introduces the Moon-fixed 8x8 GRAIL field and SPICE lunar orientation. Case 03 raises the same gravity machinery to 320x320, stressing coefficient ingestion, normalization, indexing, and computational scale. ##### Cases 04-07: bodies, perturbations, and attitude Case 04 introduces the Apollo vehicle in a high circular orbit. Case 05 adds Earth and Sun perturbations; 05A adds initial tumble. Cases 06 and 06A move a cylinder to a 28-hour elliptical orbit and isolate the initial-rate convention. Case 07 repeats that orbit with Apollo mass properties. ##### Cases 08A-08D: NRHO sensitivity All five NRHO cases run for seven days with 8x8 lunar gravity and Earth/Sun perturbations. Cases 08A and 08B select different true-anomaly starting locations. Cases 08C and 08D perturb the baseline initial radius or velocity, making them sensitive checks of exact initial-state precision and long-term integration. ##### Cases 09A-09B: local kinematics and loads These polar-orbit cases verify that a simulation can report a point away from the vehicle center of mass. Case 09 uses sensor A, 09A uses sensor B, and 09B adds the scheduled moment profile. #### Known Case 09 terrain boundary Vehicle dynamics, attitude, sensor kinematics, and Case 09B moments are fully implemented. NASA Case 09 additionally requests two independent lunar test points and a height lookup in the 80 m/pixel south-pole LOLA raster `LDEM_80S_80M`. That raster and a corresponding `TerrainGrid` ingestion path are not present in this workspace. Five Case 09 test-point/DEM fields are therefore omitted. This does not affect the propagated vehicle or sensor results, and Cases 09A and 09B do not request those fields. #### Integration, output, and validation The models use fixed-step fourth-order Runge-Kutta integration and normally report every 60 seconds. Case duration ranges from eight hours to seven days. NASA-compatible field names and units are retained so a generated file can be compared directly with participant data or uploaded to NASA's plotting page. The shared comparison tool selects and reports agreement with the closest of the eight participant files. Across the generated full-duration results, the largest Moon-inertial position difference from the closest participant is less than 6.3 mm, including the seven-day NRHO cases. This is evidence of cross-tool agreement, not a replacement for sensitivity and model-configuration review. #### Running and validating a case The normal recipe is: 1. Enter the desired case directory. 1. Run its `nesc_lunar_*.py` launcher with ordinary Python or kdpy. 1. Use `--duration`, `--output-step`, and `--integration-step` for a short run when appropriate. 1. Compare the generated CSV with the retained NASA participant data. 1. Run the shared smoke and lifecycle tests. For example: ```sh cd NESC2023checkcases_08 python3 nesc_lunar_08.py --duration 600 --output-step 60 \ --output /tmp/jain/lunar_08_short.csv cd ../NESC2023checkcases_common python3 compare_reference.py 08 python3 -m pytest -q test_lunar_cases.py ``` Further local detail is available in the [suite README](README.md), the [shared implementation guide](NESC2023checkcases_common/README.md), and each case directory's README. #### Web references - [NASA NESC 2023 Lunar Check-Cases](https://nescacademy.nasa.gov/flightsim/2023) - [NASA/TM-20240013031 Volume I](https://ntrs.nasa.gov/citations/20240013031) - [NASA/TM-20240013031 Volume II, Part 1](https://ntrs.nasa.gov/citations/20240013556) - [NASA lunar output specification](https://nescacademy.nasa.gov/flightsim/2023/output_specification) - [PDS GRAIL gravity documentation](https://pds-geosciences.wustl.edu/grail/grail-l-lgrs-2-edr-v1/grail_0001/document/shbdr.htm) - [PDS GRAIL SHADR products](https://pds-geosciences.wustl.edu/grail/grail-l-lgrs-5-rdr-v1/grail_1001/shadr/) - [NAIF SPICE reference frames](https://naif.jpl.nasa.gov/pub/naif/toolkit_docs/FORTRAN/req/frames.html) The NASA umbrella page and linked case pages remain the controlling sources for corrected inputs, definitions, downloadable data, and errata.