—
Flow effects on axial temperature profiles.
Reproduces Fig. 2 from Mei et al. (2019) Combustion and Flame (doi: 10.1016/j.combustflame.2018.11.010): influence of N2 flow rate on the measured axial gas temperature in a tube furnace.
Usage:
cd examples
python run_tube_furnace_flow_sweep.py
import importlib.util
import io
from pathlib import Path
import matplotlib.pyplot as plt
from ruamel.yaml import YAML
from bloc.batch import import_simulation
_HERE = Path(__file__).parent.absolute()
_RESULTS = _HERE / "results"
_RESULTS.mkdir(exist_ok=True)
MODULE_NAME = "simulation_tube_furnace.py"
## Configuration
Geometry is calibrated from the 0 L/min static measurement in Fig. 2.
The exit cold leg (30 cm) is longer than the inlet (22 cm), consistent
with the published profile. N2 is non-reactive; kappa_grey = 0
disables the gas-phase radiation term (N2 is IR-transparent).
_YAML_CONFIG = """\
flow_sweep:
mechanism: gri30.yaml # proper N2 transport properties
composition: "N2:1" # non-reactive inert gas
kappa_grey: 0.0 # N2 is IR-transparent
diameter: 0.016 # [m] tube inner diameter
T_wall_K: 1673.0 # [K] furnace plateau temperature
T_ambient_K: 300.0 # [K] ambient / cold-leg temperature
entry_leg: 0.22 # [m] inlet cold leg
entry_zone: 0.22 # [m] heating ramp
plateau_zone: 0.40 # [m] heated plateau
exit_leg: 0.30 # [m] outlet cold leg (longer than inlet)
total_length: 1.40 # [m] entry_leg+entry_zone+plateau+exit_zone+exit_leg
# exit_zone is computed automatically: 1.40 - 0.30 - (0.22+0.22+0.40) = 0.26 m
flow_rates_lpm: [0, 7, 13] # [L/min] flow rates to simulate
x_axis_max_cm: 140 # [cm] x-axis upper limit on the plot
y_axis_max_K: 2000 # [K] y-axis upper limit on the plot
"""
## Run
if __name__ == "__main__":
# --- Parse embedded config ------------------------------------------------
yaml_loader = YAML(typ="safe")
cfg = yaml_loader.load(io.StringIO(_YAML_CONFIG))["flow_sweep"]
mechanism = cfg.get("mechanism", "CRECK_2003_TOT_HT_SOOT_kinetics.yaml")
composition = cfg.get("composition", "CH4:1")
kappa_grey = float(cfg.get("kappa_grey", 0.5))
diameter = float(cfg.get("diameter", 0.016))
T_wall_K = float(cfg.get("T_wall_K", 1673.0))
T_ambient_K = float(cfg.get("T_ambient_K", 300.0))
entry_leg = float(cfg.get("entry_leg", cfg.get("leg", 0.2)))
exit_leg = float(cfg.get("exit_leg", cfg.get("leg", 0.2)))
entry_zone = float(cfg.get("entry_zone", cfg.get("ramp", 0.2)))
plateau_zone = float(cfg.get("plateau_zone", cfg.get("plateau", 0.6)))
total_length = float(
cfg.get("total_length", entry_leg + exit_leg + entry_zone + plateau_zone)
)
flow_rates_lpm = cfg.get("flow_rates_lpm", [0, 7, 13])
x_axis_max_cm = float(cfg.get("x_axis_max_cm", 140))
y_axis_max_K = float(cfg.get("y_axis_max_K", 2000))
exit_zone = total_length - exit_leg - entry_leg - entry_zone - plateau_zone
print("=" * 65)
print(" Tube Furnace – Flow effects on axial temperature profiles")
print("=" * 65)
print(f" Gas : {composition} ({mechanism})")
print(f" D = {diameter * 1e3:.0f} mm, L = {total_length * 1e2:.0f} cm")
print(
f" entry_leg={entry_leg * 1e2:.0f} cm entry_zone={entry_zone * 1e2:.0f} cm "
f"plateau={plateau_zone * 1e2:.0f} cm "
f"exit_zone={exit_zone * 1e2:.0f} cm exit_leg={exit_leg * 1e2:.0f} cm"
)
print(f" T_wall = {T_wall_K:.0f} K, flow rates: {flow_rates_lpm} L/min")
print("-" * 65)
# --- Load simulation module -----------------------------------------------
module = import_simulation(MODULE_NAME)
spec = importlib.util.spec_from_file_location("sim_tube", module)
sim_module = importlib.util.module_from_spec(spec)
spec.loader.exec_module(sim_module)
base_params = {
**sim_module.defaults()["parameters"],
"mechanism": mechanism,
"composition": composition,
"kappa_grey": kappa_grey,
"diameter": diameter,
"total_length": total_length,
"entry_leg": entry_leg,
"exit_leg": exit_leg,
"entry_zone": entry_zone,
"plateau_zone": plateau_zone,
"T_wall_K": T_wall_K,
"T_ambient_K": T_ambient_K,
}
# --- Simulate each flow rate ----------------------------------------------
profiles = []
for flow_lpm in flow_rates_lpm:
params = {**base_params, "flow_slm": flow_lpm}
if flow_lpm == 0:
print(" Computing 0 L/min (static wall)...")
x_cm, T_wall_profile = sim_module.get_wall_profile(params)
profiles.append((flow_lpm, x_cm, T_wall_profile, None))
else:
print(f" Running {flow_lpm} L/min...")
_, x_cm, T_wall_C, T_gas_C = sim_module.run_with_profiles(params)
T_wall_profile = T_wall_C + 273.15
T_gas_K = T_gas_C + 273.15
T_out = T_gas_K[-1] if len(T_gas_K) else 0
print(f" -> T_outlet = {T_out:.0f} K")
profiles.append((flow_lpm, x_cm, T_wall_profile, T_gas_K))
# --- Plot – style matching Mei et al. 2019 Fig. 2 ------------------------
# 0 L/min: black solid + squares
# 7 L/min: red dashed + circles
# 13 L/min: blue dash-dot + triangles
fig, ax = plt.subplots(figsize=(8, 5))
ax.set_facecolor("#e8e8e8")
styles = {
0: ("0 L/min (static)", "k-s", "black"),
7: ("7 L/min", "r--o", "red"),
13: ("13 L/min", "b-.^", "blue"),
}
for flow_lpm, x_cm, T_wall_profile, T_gas_K in profiles:
T_plot = T_wall_profile if T_gas_K is None else T_gas_K
label, fmt, color = styles[flow_lpm]
ax.plot(
x_cm,
T_plot,
fmt,
markersize=6,
markevery=max(1, len(x_cm) // 14),
color=color,
fillstyle="none",
label=label,
)
ax.set_xlabel("Distance [cm]")
ax.set_ylabel("Temperature [K]")
ax.set_title("Flow effects on axial temperature profiles in the flow reactor")
ax.set_xlim(0, x_axis_max_cm)
ax.set_ylim(0, y_axis_max_K)
ax.set_yticks([0, 400, 800, 1200, 1600, 2000])
ax.grid(True, linestyle="--", alpha=0.5)
ax.legend(loc="lower right", ncol=2)
ax.annotate(
f"T = {T_wall_K:.0f} K",
xy=(0.04, 0.95),
xycoords="axes fraction",
fontsize=10,
verticalalignment="top",
style="italic",
fontweight="bold",
)
fig.tight_layout()
out_path = _RESULTS / "flow_sweep_temperature_profiles.png"
fig.savefig(out_path, dpi=150)
print(f"\n Saved: {out_path}")
plt.show()
print("\nDone.")