Rudimentary Landing Gear: Unsteady Wheel Surface Pressure¶

The Boeing Rudimentary Landing Gear (BRLG) is a standard four-wheel configuration used as a benchmark for validating unsteady CFD methods. This case runs an unsteady k-ω SST DDES simulation of the flow around the gear at a freestream speed of 40 m/s and 180° incidence, and validates the predicted circumferential surface-pressure coefficient on the front and rear wheels against DLR BRLG revA wind-tunnel pressure-port measurements (Miquel et al.).

The study is a fork chain of three runs on the same mesh L1 DDES setup: a base run, a CTU40 warm-up continuation, and a CTU100 continuation (6000 steps, Δt = 0.0001 s). Only the final CTU100 case is used for the comparison, with the circumferential Cp averaged over the 40-100 CTU window.

The notebook runs top-to-bottom in a single kernel: load the geometry, submit the fork chain, wait for completion, then post-process the in-kernel CTU100 case into the two published figures.

Requirements¶

This notebook runs against the Flow360 Python API. Install it and configure your API key by following the installation & setup guide.

This notebook targets the latest Flow360 Python client. It fetches the case's root asset and reference data with flow360.examples.download_benchmark_assets, added in 25.10.4. Install the client (25.10.4 or newer) before running:

pip install "flow360>=25.10.4"

Beyond flow360, the notebook needs a few plotting/data packages:

pip install numpy pandas matplotlib

Run the cells top to bottom in a single kernel.

Imports¶

In [ ]:
import os
import re
import tarfile
import tempfile
from pathlib import Path

import flow360 as fl
from flow360 import u
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

from flow360.examples import download_benchmark_assets

Input data¶

The post-processing compares Flow360 against the DLR BRLG revA experimental Cp measurements. Fetch those reference CSVs from the public benchmark bucket; this recreates a local ./ref_data/ directory that the plotting cells read from.

In [ ]:
download_benchmark_assets("Rudimentary_Landing_Gear", "ref_data")

Load project¶

The root asset for this case is the landing-gear geometry. Download it from the public benchmark bucket, start a fresh project from it, and group the geometry faces by their faceId tag so the boundary conditions can address them by name.

In [ ]:
root_asset_files = download_benchmark_assets("Rudimentary_Landing_Gear", "root_assets")
# The snapshot includes processed copies under results/; select the source
# CAD file(s) to rebuild the geometry project from.
geometry_files = [
    f for f in root_asset_files
    if "/results/" not in f and "/logs/" not in f
    and f.lower().endswith((".csm", ".egads", ".stp", ".step", ".iges", ".igs", ".stl"))
]
project = fl.Project.from_geometry(geometry_files, name="Rudimentary Landing Gear")
geometry = project.geometry
geometry.group_faces_by_tag("faceId")

Case constants¶

Operating point, wheel geometry, and the geometry face groups. The gear surface is split into the viscous wheel walls, the freestream far-field faces, and the tunnel side (slip) walls; the remaining body00001_face* faces are the wheels. Surface and slice output fields are also fixed here.

In [ ]:
VMAG = 40
AOA = 180
CTU = 0.01
CLF = 1
DIAMETER = 0.4064
WIDTH = 0.37 * DIAMETER
WHEEL_CENTER_X = 0.235712
WHEEL_CENTER_Y = 0.5047 * DIAMETER / 2 + WIDTH / 2
WHEEL_CENTER_Z = 0.749808

SIDE_SURFACES = [
    "body00001_face00015",
    "body00001_face00018",
    "body00001_face00019",
    "body00001_face00042",
    "body00001_face00045",
    "body00001_face00046",
]
FREESTREAM_SURFACES = [
    "body00001_face00016",
    "body00001_face00017",
    "body00001_face00043",
    "body00001_face00044",
]
WHEEL_SURFACES = [
    f"body00001_face{i:05d}"
    for i in range(1, 73)
    if f"body00001_face{i:05d}" not in {*SIDE_SURFACES, *FREESTREAM_SURFACES}
]

SURF_FIELDS = ["Cp", "yPlus", "Cf", "CfVec"]
SLICE_FIELDS = ["Cp", "velocity_m_per_s"]

Meshing setup¶

Build the meshing parameters from the geometry. A boundary layer is grown on the wheels (first-layer thickness 3 µm) while it is switched off on the side and freestream faces via a passive unchanged spacing. The wheel surfaces get a finer surface refinement, and two boxed wake refinements (near and far) resolve the downstream flow. A y-plane mesh slice is emitted at each wheel plane.

In [ ]:
def build_meshing(farfield_zone, slice_y1, slice_y2, geometry):
    with fl.SI_unit_system:
        return fl.MeshingParams(
            defaults=fl.MeshingDefaults(
                surface_max_edge_length=0.032 * u.m,
                curvature_resolution_angle=12 * u.deg,
                boundary_layer_first_layer_thickness=0.000003 * u.m,
            ),
            volume_zones=[farfield_zone],
            refinements=[
                fl.PassiveSpacing(
                    name="BL_OFF",
                    faces=[geometry[name] for name in SIDE_SURFACES]
                    + [geometry[name] for name in FREESTREAM_SURFACES],
                    type="unchanged",
                ),
                fl.SurfaceRefinement(
                    name="Wheel_Refinement",
                    faces=[geometry[name] for name in WHEEL_SURFACES],
                    max_edge_length=0.006 * fl.u.m,
                    curvature_resolution_angle=4 * fl.u.deg,
                ),
                fl.UniformRefinement(
                    name="Wake_Near_Refinement",
                    entities=[
                        fl.Box.from_principal_axes(
                            name="Wake_Near",
                            center=(-150, 0, 750) * fl.u.mm,
                            size=(1600, 1000, 800) * fl.u.mm,
                            axes=[(1, 0, 0), (0, 1, 0)],
                        )
                    ],
                    spacing=0.008 * fl.u.m,
                ),
                fl.UniformRefinement(
                    name="Wake_Far_Refinement",
                    entities=[
                        fl.Box.from_principal_axes(
                            name="Wake_Far",
                            center=(800, 0, 750) * fl.u.mm,
                            size=(5000, 1600, 1600) * fl.u.mm,
                            axes=[(1, 0, 0), (0, 1, 0)],
                        )
                    ],
                    spacing=0.032 * fl.u.m,
                ),
            ],
            outputs=[fl.MeshSliceOutput(slices=[slice_y1, slice_y2])],
        )

Physics setup¶

Boundary conditions and the flow solver. The wheels are viscous walls, the far-field faces are freestream, and the tunnel sides are slip walls. The fluid is solved with a small numerical-dissipation factor and the k-ω SST turbulence model running in DDES (Detached Eddy Simulation) hybrid mode for the unsteady, separated wake.

In [ ]:
def build_models(geometry):
    with fl.SI_unit_system:
        return [
            fl.Wall(surfaces=[geometry[name] for name in WHEEL_SURFACES], name="Wall"),
            fl.Freestream(
                surfaces=[geometry[name] for name in FREESTREAM_SURFACES],
                name="Freestream",
            ),
            fl.SlipWall(
                surfaces=[geometry[name] for name in SIDE_SURFACES], name="SlipWall"
            ),
            fl.Fluid(
                navier_stokes_solver=fl.NavierStokesSolver(
                    numerical_dissipation_factor=0.2
                ),
                turbulence_model_solver=fl.KOmegaSST(
                    hybrid_model=fl.DetachedEddySimulation(shielding_function="DDES")
                ),
            ),
        ]

Outputs setup¶

The requested solver outputs: surface fields on the wheels, y-plane slice fields, and four rings of surface probes (front/rear wheel, each side) that sample Cp around the wheel circumference. _probe_locations places 180 probe points on a circle of radius DIAMETER/2 about each wheel center; these feed the SurfaceProbeOutput monitors the post-processing reads.

In [ ]:
def _probe_locations():
    radius = DIAMETER / 2
    theta = np.linspace(0, (360 - 2) * np.pi / 180, 180)
    fw1 = [
        [WHEEL_CENTER_X + radius * np.cos(t), WHEEL_CENTER_Y, WHEEL_CENTER_Z + radius * np.sin(t)]
        for t in theta
    ]
    fw2 = [
        [WHEEL_CENTER_X + radius * np.cos(t), -WHEEL_CENTER_Y, WHEEL_CENTER_Z + radius * np.sin(t)]
        for t in theta
    ]
    rw1 = [
        [-WHEEL_CENTER_X + radius * np.cos(t), WHEEL_CENTER_Y, WHEEL_CENTER_Z + radius * np.sin(t)]
        for t in theta
    ]
    rw2 = [
        [-WHEEL_CENTER_X + radius * np.cos(t), -WHEEL_CENTER_Y, WHEEL_CENTER_Z + radius * np.sin(t)]
        for t in theta
    ]
    return fw1, fw2, rw1, rw2


def build_outputs(
    slice_y1, slice_y2, geometry, fw1_probes, fw2_probes, rw1_probes, rw2_probes
):
    with fl.SI_unit_system:
        return [
            fl.SurfaceOutput(
                output_format="tecplot",
                output_fields=SURF_FIELDS,
                surfaces=[geometry[name] for name in WHEEL_SURFACES],
            ),
            fl.SliceOutput(
                slices=[slice_y1, slice_y2],
                output_format="tecplot",
                output_fields=SLICE_FIELDS,
            ),
            fl.SurfaceProbeOutput(
                name="FW1_Probes",
                entities=[
                    fl.Point(name=f"FW1_Point_{i}", location=tuple(p) * fl.u.m)
                    for i, p in enumerate(fw1_probes)
                ],
                target_surfaces=[
                    geometry["body00001_face00062"],
                    geometry["body00001_face00066"],
                ],
                output_fields=["Cp"],
            ),
            fl.SurfaceProbeOutput(
                name="FW2_Probes",
                entities=[
                    fl.Point(name=f"FW2_Point_{i}", location=tuple(p) * fl.u.m)
                    for i, p in enumerate(fw2_probes)
                ],
                target_surfaces=[
                    geometry["body00001_face00009"],
                    geometry["body00001_face00010"],
                ],
                output_fields=["Cp"],
            ),
            fl.SurfaceProbeOutput(
                name="RW1_Probes",
                entities=[
                    fl.Point(name=f"RW1_Point_{i}", location=tuple(p) * fl.u.m)
                    for i, p in enumerate(rw1_probes)
                ],
                target_surfaces=[
                    geometry["body00001_face00003"],
                    geometry["body00001_face00071"],
                ],
                output_fields=["Cp"],
            ),
            fl.SurfaceProbeOutput(
                name="RW2_Probes",
                entities=[
                    fl.Point(name=f"RW2_Point_{i}", location=tuple(p) * fl.u.m)
                    for i, p in enumerate(rw2_probes)
                ],
                target_surfaces=[
                    geometry["body00001_face00037"],
                    geometry["body00001_face00038"],
                ],
                output_fields=["Cp"],
            ),
        ]

Simulation Params¶

Assemble the full SimulationParams: the meshing, reference geometry (based on the wheel diameter and width), the operating condition (40 m/s at 180° incidence, DLR tunnel thermal state), the models, and the outputs. Time stepping is unsteady with an adaptive CFL, 30 pseudo-steps per physical step, a step size of one hundredth of a CTU, and 2000 physical steps for the base run (the CTU100 continuation later raises this to 6000). We then instantiate the far-field zone, the two wheel-plane slices, and the probe rings, and build the base run_params.

In [ ]:
def make_run_params(
    geometry,
    farfield_zone,
    slice_y1,
    slice_y2,
    fw1_probes,
    fw2_probes,
    rw1_probes,
    rw2_probes,
):
    with fl.SI_unit_system:
        return fl.SimulationParams(
            meshing=build_meshing(farfield_zone, slice_y1, slice_y2, geometry),
            reference_geometry=fl.ReferenceGeometry(
                moment_center=(0.0, 0.0, 0.0) * u.m,
                moment_length=(DIAMETER, DIAMETER, DIAMETER) * u.m,
                area=DIAMETER * WIDTH * u.m * u.m,
            ),
            operating_condition=fl.AerospaceCondition(
                velocity_magnitude=VMAG * u.m / u.s,
                alpha=AOA * u.deg,
                thermal_state=fl.ThermalState(
                    temperature=306.15 * fl.u.K,
                    density=1.153 * fl.u.kg / fl.u.m**3,
                ),
            ),
            time_stepping=fl.Unsteady(
                max_pseudo_steps=30,
                step_size=CTU / 100,
                steps=2000,
                CFL=fl.AdaptiveCFL(convergence_limiting_factor=CLF),
            ),
            models=build_models(geometry),
            outputs=build_outputs(
                slice_y1,
                slice_y2,
                geometry,
                fw1_probes,
                fw2_probes,
                rw1_probes,
                rw2_probes,
            ),
        )


farfield_zone = fl.UserDefinedFarfield()
slice_y1 = fl.Slice(name="SliceY1", normal=(0, 1, 0), origin=(0, WHEEL_CENTER_Y, 0) * u.m)
slice_y2 = fl.Slice(name="SliceY2", normal=(0, 1, 0), origin=(0, -WHEEL_CENTER_Y, 0) * u.m)
fw1_probes, fw2_probes, rw1_probes, rw2_probes = _probe_locations()

run_params = make_run_params(
    geometry,
    farfield_zone,
    slice_y1,
    slice_y2,
    fw1_probes,
    fw2_probes,
    rw1_probes,
    rw2_probes,
)

Submit cases¶

Submit the fork chain to Flow360. The base run advances 2000 physical steps (20 CTU) from the freestream initial condition. CTU40 forks from the base and runs another 2000 steps to warm up the wake. CTU100 forks from CTU40 with the step count raised to 6000 (60 more CTU) and is the case used for the comparison. All three share the same mesh, solver, and outputs and use the beta mesher.

In [ ]:
SOLVER_VERSION = os.environ.get("SOLVER_VERSION_OVERRIDE", "release-25.9")

base_case = project.run_case(
    params=run_params,
    name=f"RLG_MeshL1_kwSST-DDES_{SOLVER_VERSION}",
    use_beta_mesher=True,
    solver_version=SOLVER_VERSION,
)

ctu40_case = project.run_case(
    params=run_params,
    name=f"RLG_MeshL1_kwSST-DDES_CTU40_{SOLVER_VERSION}",
    fork_from=base_case,
    use_beta_mesher=True,
    solver_version=SOLVER_VERSION,
)

run_params.time_stepping.steps = 6000
ctu100_case = project.run_case(
    params=run_params,
    name=f"RLG_MeshL1_kwSST-DDES_CTU100_{SOLVER_VERSION}",
    fork_from=ctu40_case,
    use_beta_mesher=True,
    solver_version=SOLVER_VERSION,
)

submitted_cases = [base_case, ctu40_case, ctu100_case]

Wait for completion¶

Block until every case in the fork chain has finished. This can take a long time (the chain is 100 CTU of an unsteady DDES solve), since each fork must complete before the next continues from it.

In [ ]:
for case in submitted_cases:
    case.wait()

Postprocessing¶

The published results are two figures: the circumferential surface-pressure coefficient on the front wheel and on the rear wheel, each compared against the DLR BRLG revA experiment.

The workflow reads the surface-probe monitors from the finished CTU100 case (ctu100_case, reused directly from the sweep, no project is re-opened). For each wheel the two probe rings (FW1/FW2, RW1/RW2) are averaged, the Cp signal is time-averaged over the 40-100 CTU window, and the result is plotted as Cp versus circumferential angle θ.

Setup¶

Styling colors, the averaging window, the experimental reference file locations, and the helper functions that turn the downloaded probe monitors into a time-averaged Cp-vs-θ curve for each wheel: get_data_at_last_pseudo_step keeps one row per physical step, average_monitor_data averages the two rings of a wheel, average_data_over_indices time-averages over the CTU window, and extract_cp_vs_theta orders the probe points into a Cp(θ) curve.

In [ ]:
FLOW360_COLOR = "#00643c"   # Flow360 results
REFERENCE_COLOR = "black"   # experimental reference

TIME_STEP = 0.0001
CTU_TIME = 0.01
AVG_CTU_START = 40
AVG_CTU_END = 100

REF_FW_PATH = "ref_data/miquel_revA_experiments_cp_cprms_front_wheel.csv"
REF_RW_PATH = "ref_data/miquel_revA_experiments_cp_cprms_rear_wheel_modified.csv"


def get_data_at_last_pseudo_step(filename):
    dataframe = pd.read_csv(filename, skipinitialspace=True)
    dataframe = dataframe.loc[:, ~dataframe.columns.str.contains("^Unnamed")]
    data_raw = dataframe.to_dict("list")
    data = {key: [] for key in data_raw}
    n = len(data_raw["physical_step"])
    for i in range(0, n - 1):
        if data_raw["physical_step"][i] != data_raw["physical_step"][i + 1]:
            for key in data_raw:
                data[key].append(data_raw[key][i])
    for key in data_raw:
        data[key].append(data_raw[key][n - 1])
    return data


def average_monitor_data(data_1, data_2, tag_1, tag_2, output_tag):
    averaged_data = {}
    for key in data_1:
        key_2 = key.replace(tag_1, tag_2)
        output_key = key.replace(tag_1, output_tag)
        values_1 = np.asarray(data_1[key], dtype=float)
        values_2 = np.asarray(data_2[key_2], dtype=float)
        averaged_data[output_key] = ((values_1 + values_2) / 2.0).tolist()
    return averaged_data


def average_data_over_indices(data, indices):
    if len(indices) == 0:
        raise ValueError("No indices found for averaging.")
    averaged_data = {}
    for key, values in data.items():
        values_array = np.asarray(values, dtype=float)
        averaged_data[key] = float(np.mean(values_array[indices]))
    return averaged_data


def extract_cp_vs_theta(avg_data, tag):
    pattern = re.compile(rf"^{tag}_Probes_{tag}_Point_(\d+)_Cp$")
    indexed_cp = []
    for key, value in avg_data.items():
        match = pattern.match(key)
        if match:
            indexed_cp.append((int(match.group(1)), float(value)))
    if not indexed_cp:
        raise ValueError(f"No Cp point data found for {tag}.")
    indexed_cp.sort(key=lambda x: x[0])
    cp = np.asarray([value for _, value in indexed_cp], dtype=float)
    theta = np.linspace(0.0, 360.0, len(cp), endpoint=False)
    return theta, cp


def find_monitor_file(root, stem):
    candidates = sorted(root.rglob(f"monitor_{stem}_v2.csv"))
    if not candidates:
        candidates = sorted(root.rglob(f"monitor_{stem}.csv"))
    if not candidates:
        raise FileNotFoundError(f"monitor file for {stem!r} not found under {root}")
    return candidates[0]


def download_monitors(case, destination):
    destination.mkdir(parents=True, exist_ok=True)
    monitors_tar = destination / "monitors.tar.gz"
    case.results.monitors.download(str(monitors_tar), overwrite=True)
    with tarfile.open(monitors_tar) as archive:
        archive.extractall(destination)


def cp_curves_from_monitors(monitor_root):
    fw1_data = get_data_at_last_pseudo_step(find_monitor_file(monitor_root, "FW1_Probes"))
    fw2_data = get_data_at_last_pseudo_step(find_monitor_file(monitor_root, "FW2_Probes"))
    rw1_data = get_data_at_last_pseudo_step(find_monitor_file(monitor_root, "RW1_Probes"))
    rw2_data = get_data_at_last_pseudo_step(find_monitor_file(monitor_root, "RW2_Probes"))

    fw_data = average_monitor_data(fw1_data, fw2_data, "FW1", "FW2", "FW")
    rw_data = average_monitor_data(rw1_data, rw2_data, "RW1", "RW2", "RW")

    physical_time = np.linspace(
        fw_data["physical_step"][0] * TIME_STEP,
        fw_data["physical_step"][-1] * TIME_STEP,
        len(fw_data["physical_step"]),
    )
    ctu = physical_time / CTU_TIME
    avg_idx = np.where((ctu > AVG_CTU_START) & (ctu <= AVG_CTU_END))[0]
    avg_fw_data = average_data_over_indices(fw_data, avg_idx)
    avg_rw_data = average_data_over_indices(rw_data, avg_idx)

    return extract_cp_vs_theta(avg_fw_data, "FW") + extract_cp_vs_theta(avg_rw_data, "RW")


def plot_wheel_cp(theta_flow, cp_flow, exp_df, title, flip_theta=False):
    theta_plot = 360.0 - theta_flow if flip_theta else theta_flow
    sort_idx = np.argsort(theta_plot)
    theta_sorted = theta_plot[sort_idx]
    cp_sorted = cp_flow[sort_idx]

    exp_theta = exp_df["theta"].to_numpy(dtype=float)
    exp_cp = exp_df["cp"].to_numpy(dtype=float)
    if flip_theta:
        exp_theta = 360.0 - exp_theta
    exp_sort_idx = np.argsort(exp_theta)
    exp_theta_sorted = exp_theta[exp_sort_idx]
    exp_cp_sorted = exp_cp[exp_sort_idx]

    plt.figure(figsize=(5, 4))
    plt.plot(theta_sorted, cp_sorted, color=FLOW360_COLOR, linewidth=2, label="Flow360")
    plt.plot(
        exp_theta_sorted,
        exp_cp_sorted,
        color=REFERENCE_COLOR,
        linewidth=2,
        marker="o",
        markersize=2.5,
        label="Experiment",
    )
    plt.title(title)
    plt.xlabel("Theta (deg)")
    plt.ylabel("Cp")
    plt.xlim(0, 360)
    plt.xticks([0, 90, 180, 270, 360])
    plt.ylim(-1.5, 1.5)
    plt.grid(True, alpha=0.3)
    plt.legend(fontsize=8, loc="upper right")
    plt.tight_layout()


os.makedirs("results", exist_ok=True)

Load the CTU100 monitors¶

Download the surface-probe monitors from the in-kernel CTU100 case and reduce them to a time-averaged Cp(θ) curve per wheel, then read the two experimental reference CSVs.

In [ ]:
with tempfile.TemporaryDirectory() as temp_dir:
    monitor_root = Path(temp_dir)
    download_monitors(ctu100_case, monitor_root)
    fw_theta, fw_cp, rw_theta, rw_cp = cp_curves_from_monitors(monitor_root)

exp_fw = pd.read_csv(REF_FW_PATH)
exp_rw = pd.read_csv(REF_RW_PATH)

Front-wheel Cp vs θ¶

Flow360's front-wheel circumferential Cp against the DLR revA experiment.

In [ ]:
plot_wheel_cp(fw_theta, fw_cp, exp_fw, "Front-wheel Cp vs θ", flip_theta=False)
plt.savefig("results/fw_cp_theta.png", dpi=400)
plt.show()

Rear-wheel Cp vs θ¶

Flow360's rear-wheel circumferential Cp against the DLR revA experiment. The angle convention is mirrored (flip_theta=True) so the rear wheel matches the experimental reference orientation.

In [ ]:
plot_wheel_cp(rw_theta, rw_cp, exp_rw, "Rear-wheel Cp vs θ", flip_theta=True)
plt.savefig("results/rw_cp_theta.png", dpi=400)
plt.show()