GABLS1 LES Intercomparison Study for Stable Boundary Layers: Description
Last updated: May 2026
Warning
All GABLS1 runs violate the near-surface resolution guideline \(z_1 > 50\,z_0 = 5\) m (\(z_0 = 0.1\) m; Basu & Lacser, 2017). The first half-model level \(z_1\) for each resolution is:
\(64^3\): \(z_1 \approx 3.13\) m
\(128^3\): \(z_1 \approx 1.57\) m
\(256^3\): \(z_1 \approx 0.78\) m
\(384^3\): \(z_1 \approx 0.52\) m
A Roughness Sublayer (RSL) parameterization is planned for a future release of JAX-ALFA; all these runs should be re-simulated once RSL support is available. Results should be interpreted accordingly.
Basu, S., and Lacser, A. (2017). Boundary-Layer Meteorology, 163, 197–207. https://doi.org/10.1007/s10546-016-0225-y
References
Beare, R. J., and Coauthors. (2006). An intercomparison of large-eddy simulations of the stable boundary layer. Boundary-Layer Meteorology, 118, 247–272.
Basu, S. and Porté-Agel, F. (2006). Large-eddy simulation of stably stratified atmospheric boundary layer turbulence: A scale-dependent dynamic modeling approach. Journal of the Atmospheric Sciences, 63, 2074–2091.
Dai, Y., Basu, S., Maronga, B., and de Roode, S. R. (2020). Addressing the grid-size sensitivity issue in large-eddy simulations of stable boundary layers. Boundary-Layer Meteorology, 178, 63–89.
Case Description
The GABLS1 case simulates a stably stratified boundary layer driven by prescribed surface cooling. Starting from an initially (conventionally) neutral profile, the boundary layer approaches a quasi-steady stable state.
Parameter |
Value |
|---|---|
Reference run |
|
Domain |
400 m x 400 m x 400 m |
Reference grid |
\(128^3\) |
Run matrix |
\(64^3\), \(128^3\), \(256^3\), \(384^3\) x LASDD-SM, LASDD-WL, LAD-SM, LAD-WL x SP, DP |
Simulation time |
32400 s (9 h) |
Reference time step |
0.1 s |
Geostrophic wind |
\(U_g = 8\) m s\(^{-1}\), \(V_g = 0\) m s\(^{-1}\) |
Coriolis parameter |
\(f = 1.39 \times 10^{-4}\) s\(^{-1}\) |
Roughness lengths |
\(z_{0m} = 0.1\) m, \(z_{0T} = 0.1\) m |
Thermal surface option |
|
Sensible heat flux |
0 K m s\(^{-1}\) |
Reference temperature |
\(T_0 = 265\) K |
Inversion strength |
0.01 K m\(^{-1}\) |
SGS model |
selected by the user |
Filter-to-grid ratio |
|
Damping layer |
above \(z = 300\) m |
Output statistics interval |
60 s |
3D output interval |
SimTime s |
Geostrophic forcing option |
|
Advection forcing option |
|
Moisture option |
|
Surface boundary condition |
Prescribed surface cooling |
Reference averaging window |
hours 8-9 |
The representative configuration below is copied from examples/SBL_GABLS1/runs/128x128x128_LASDD_SM_DP/Config.py.
Reference Configuration
[1]:
# Copyright (C) 2025 Sukanta Basu
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program. If not, see <https://www.gnu.org/licenses/>.
"""
File: Config.py
===============
:Author: Sukanta Basu
:AI Assistance: Claude.AI (Anthropic) is used for documentation,
code restructuring, and performance optimization
:Date: 2026-05-20
:Description: Stable BL benchmark after Beare et al. (2006), GABLS1.
Domain 400x400x400 m; surface cooling 0.25 K/hr.
Grid: 128x128x128, SGS: LASDD-SM,
Precision: double.
"""
# ============================================================
# Imports
# ============================================================
import numpy as np
# ============================================================
# User Input
# ============================================================
# ------------------------------------------------------------
# Platform options
# ------------------------------------------------------------
use_double_precision = True
# 0: use CPU, 1: use GPU
optGPU = 1
GPU_ID = 0
# ------------------------------------------------------------
# Domain configuration
# ------------------------------------------------------------
# Domain size (m)
l_x = 400
l_y = 400
l_z = 400
# Number of grid points
nx = 128
ny = 128
nz = 128
# ------------------------------------------------------------
# Time integration configuration
# ------------------------------------------------------------
# Change this if it is a restart run
istep = 1
# Time stepping and simulation time
dt = 0.1 # unit: sec
SimTime = 9 * 3600 # unit: sec
# Galilean transformation (m/s)
Ugal = 5
# ------------------------------------------------------------
# Surface configuration
# ------------------------------------------------------------
# optSurfFlux: 0 = homogeneous, 1 = heterogeneous
optSurfFlux = 0
# optSurfBC: 0 = constant flux
# 1 = time-varying flux (from SurfaceBCFile)
# 2 = time-varying surface temperature (from SurfaceBCFile)
optSurfBC = 2
# Roughness lengths (m)
z0m = 0.1
z0T = 0.1
# Screen-level temperature reference height (m); 0 = use z0T
zTemperature = 0.0
# SensibleHeatFlux: not used when optSurfBC >= 1; set to 0 as placeholder
SensibleHeatFlux = 0.0 # K m/s
# Path to surface BC file (relative to run directory)
SurfaceBCFile = 'input/SurfaceBC.npz'
# ------------------------------------------------------------
# Forcing configuration
# ------------------------------------------------------------
# Geostrophic wind option:
# 0 = constant Ug2, Vg2 from Config
# 1 = time + height varying, loaded from GeoWindFile
optGeoWind = 0
# Constant geostrophic wind (m/s)
Ug2 = 8
Vg2 = 0
# Path to geostrophic wind file (relative to run directory)
GeoWindFile = 'input/GeoWind.npz'
# Coriolis parameter (1/s)
f_coriolis = 0.000139
# Potential temperature lapse rate above domain top (K/m)
inversion = 0.01
# Buoyancy calculation: 0 = use reference T_0, 1 = use local THv
optBuoyancy = 1
# Reference temperature (K)
T_0 = 265
# ------------------------------------------------------------
# Subgrid-scale configuration
# ------------------------------------------------------------
# SGS model: 0 = Static SM, 1 = LASDD-SM, 2 = LASDD-WL, 3 = LAD-SM, 4 = LAD-WL
optSgs = 1
# Dynamic SGS update frequency (every N steps)
dynamicSGS_call_time = 1
# Filter to grid ratio (FGR=1: implicit + dealiasing; FGR>=2: explicit)
FGR = 2
# Initial SGS coefficients (used before first dynamic update)
Cs2 = 0.1 ** 2 # SM models: initial Cs^2
Cwl = 0.1 ** 2 # WL models: initial C_WL
Cs2PrRatio = Cs2 / 1.0
CwlPrRatio = Cwl / 1.0
# ------------------------------------------------------------
# Damping layer configuration
# ------------------------------------------------------------
optDamping = 1 # 1: activate Rayleigh damping
z_damping = 300 # unit: m
RelaxTime = 60 # unit: s
# ------------------------------------------------------------
# Statistics computation
# ------------------------------------------------------------
SampleInterval_sec = 10.0 # collect a sample every N s
OutputInterval_sec = 60.0 # output averaged stats every N s
Output3DInterval_sec = SimTime # output 3D fields only at the end of the run
# ------------------------------------------------------------
# Large-scale advection forcing
# ------------------------------------------------------------
# 0: none, 1: time/height-varying (from AdvectionFile)
optAdvection = 0
AdvectionFile = 'input/AdvForcing.npz'
# ------------------------------------------------------------
# Moisture configuration
# ------------------------------------------------------------
# 0: dry run, 1: prognostic specific humidity Q
optMoisture = 0
# Screen-level moisture reference height (m); 0 = use z0T
zMoisture = 0.0
# Surface moisture flux (kg/kg m/s); used when optMoistureSurfBC = 0
MoistureFlux = 0.0
# 0: constant flux, 1: time-varying flux, 2: time-varying surface Q
optMoistureSurfBC = 0
MoistureSurfaceBCFile = 'input/MoistureSurfaceBC.npz'
# Specific humidity lapse rate above domain top (kg/kg/m); 0 = zero gradient
q_inversion = 0.0
# Pressure solver: 0 = LU (original), 1 = Thomas (tridiagonal, faster)
optPressureSolver = 1