The src subpackage

This subpackage contains the main code for the BLDFM model. It includes modules for the PBL model, the solver, and utility functions.

It is organized into three main modules:

  1. src.pbl_model: Computes vertical profiles of horizontal wind and eddy diffusivity in the planetary boundary layer using Monin-Obukhov Similarity Theory (MOST).

  2. src.solver: Solves the steady-state advection-diffusion equation for scalar concentration fields using FFT-based methods and numerical integration schemes.

  3. src.utils: Provides utility functions for preprocessing, synthetic input generation, and diagnostics.

Modules

src.pbl_model module

src.pbl_model.phi(x)[source]

Stability correction function for eddy diffusivity in Monin-Obukhov Similarity Theory (MOST).

Parameters:

x (float or numpy.ndarray) – Dimensionless stability parameter (z / mol).

Returns:

Stability correction factor for eddy diffusivity.

Return type:

numpy.ndarray

Notes

  • For stable conditions (x > 0), a linear relationship is used.

  • For unstable conditions (x < 0), the Businger-Dyer formulation is applied.

src.pbl_model.psi(x)[source]

Stability correction function for Monin-Obukhov Similarity Theory (MOST).

Parameters:

x (float or numpy.ndarray) – Dimensionless stability parameter (z / mol).

Returns:

Stability correction factor for momentum.

Return type:

numpy.ndarray

Notes

  • For stable conditions (x > 0), a linear relationship is used.

  • For unstable conditions (x < 0), the Businger-Dyer formulation is applied.

src.pbl_model.vertical_profiles(n, meas_height, wind, ustar=0.2, mol=1000000000.0, prsc=1.0, closure='MOST', z0=-1000000000.0, z0_min=0.001, z0_max=2.0, tke=None)[source]

Computes vertical profiles of horizontal wind components and eddy diffusivity in the planetary boundary layer (PBL) based on Monin-Obukhov Similarity Theory (MOST) or other closure models.

Parameters:

  • n (int) – Number of vertical grid points.

  • meas_height (float) – Measurement height above the ground.

  • wind (tuple of floats) – Zonal (u) and meridional (v) wind components at the measurement height.

  • ustar (float or numpy.ndarray) – Friction velocity [m/s].

  • mol (float or numpy.ndarray, optional) – Monin-Obukhov length. Default is 1e9 (neutral conditions).

  • prsc (float, optional) – Prandtl or Schmidt number. Default is 1.0.

  • closure (str, optional) – Closure model to use. Options are “MOST”, “CONSTANT”, or “OAAHOC”. Default is “MOST”.

  • z0 (float or numpy.ndarray, optional) – Roughness length. Default is -1e9 (auto-calculated).

  • z0_min (float, optional) – Minimum allowable roughness length. Default is 0.001.

  • z0_max (float, optional) – Maximum allowable roughness length. Default is 2.0.

  • tke (float or numpy.ndarray, optional) – Turbulent kinetic energy [m²/s²] for the “OAAHOC” closure. If not provided, it will default to 1.0.

Returns:

  • z (numpy.ndarray): 1D array of vertical grid points.

  • profiles (tuple of numpy.ndarray): 1D arrays of horizontal wind components (u, v) and eddy diffusivity (K) at each vertical grid point.

Return type:

tuple

Raises:

ValueError – If invalid closure type is provided.

Notes

  • The “OAAHOC” closure uses a one-and-a-half order closure model based on Schumann-Lilly.

  • The “MOST” closure uses Monin-Obukhov Similarity Theory.

References

  • Kormann, R., & Meixner, F. X. (2001). An analytical footprint model for non-neutral stratification. Boundary-Layer Meteorology, 99(2), 207–224. https://doi.org/10.1023/A:1018991015119

  • Schumann, U. (1991). Subgrid length-scales for large-eddy simulation of stratified turbulence. Theoretical and Computational Fluid Dynamics, 2(5), 279–290. https://doi.org/10.1007/BF00271468

src.solver module

src.solver.ivp_solver(fftpq, profiles, z, Lx, Ly, method='SIE')[source]

Solves the initial value problem resulting from the discretization of the steady-state advection-diffusion equation using the Fast Fourier Transform.

Parameters:

  • fftpq (tuple of ndarray) – Tuple containing the initial Fourier-transformed pressure and flux fields (fftp0, fftq0).

  • profiles (tuple of ndarray) – Tuple containing 1D arrays of vertical profiles of zonal wind, meridional wind [m/s], and eddy diffusivity [m²/s].

  • z (ndarray of float) – 1D array of vertical grid points from z0 to zm [m].

  • Lx (ndarray of float) – 2D array of zonal wavenumbers.

  • Ly (ndarray of float) – 2D array of meridional wavenumbers.

  • method (str, optional) –

    Method for solving the initial value problem. Options are:

    • SIE (Semi-Implicit Euler, default)

    • EI (Exponential Integrator)

    • TSEI3 (Taylor Series Exponential Integrator, 3rd order)

    • EE (Explicit Euler)

Returns:

  • fftp (ndarray of complex) – Fourier-transformed pressure field after solving the IVP.

  • fftq (ndarray of complex) – Fourier-transformed flux field after solving the IVP.

src.solver.steady_state_transport_solver(srf_flx, z, profiles, domain, modes=(512, 512), meas_pt=(0.0, 0.0), srf_bg_conc=0.0, footprint=False, analytic=False, fetch=-1000000000.0, ivp_method='SIE')[source]

Solves the steady-state advection-diffusion equation for a concentration with flux boundary condition given vertical profiles of wind and eddy diffusivity using the Fourier, linear shooting, and semi-implicit Euler methods.

Parameters:

  • srf_flx (ndarray of float) – 2D field of surface kinematic flux at z=z0 [m/s].

  • z (ndarray of float) – 1D array of vertical grid points from z0 to zm [m].

  • profiles (tuple of ndarray) – Tuple containing 1D arrays of vertical profiles of zonal wind, meridional wind [m/s], and eddy diffusivity [m²/s].

  • domain (tuple of float) – Tuple containing domain sizes (xmax, ymax) [m].

  • modes (tuple of int, optional) – Tuple containing the number of zonal and meridional Fourier modes (nlx, nly). Default is (512, 512).

  • meas_pt (tuple of float, optional) – Coordinates of the measurement point (xm, ym) [m] relative to srf_flx, where the origin is at srf_flx[0, 0]. Default is (0.0, 0.0).

  • srf_bg_conc (float, optional) – Surface background concentration at z=z0 [scalar_unit]. Default is 0.0.

  • footprint (bool, optional) – If True, activates footprints (Green’s function) for output. Default is False.

  • analytic (bool, optional) – If True, uses the analytic solution for constant wind and eddy diffusivity. Default is False.

  • fetch (float, optional) – Width of the zero-flux halo around the domain [m]. Default is -1e9, which sets the fetch to max(xmax, ymax).

  • ivp_method (str, optional) – Method for solving the initial value problem. Options are “SIE” (default), “EI”, “TSEI3”, or “EE”.

Returns:

  • srf_conc (ndarray of float) – 2D field of surface concentrations at z=z0.

  • bg_conc (float) – Background concentration at z=zm.

  • conc (ndarray of float) – 2D field of concentration at z=zm or Green’s function.

  • flx (ndarray of float) – 2D field of kinematic flux at z=zm or footprint.

src.utils module

src.utils.compute_wind_fields(u_rot, wind_dir)[source]

Computes the zonal (u) and meridional (v) wind components from a rotated wind speed and direction.

Parameters:

  • u_rot (float) – Rotated wind speed.

  • wind_dir (float) – Wind direction in degrees (clockwise from north).

Returns:

A tuple (u, v) where:

  • u (float): Zonal wind component (east-west).

  • v (float): Meridional wind component (north-south).

Return type:

tuple

src.utils.ideal_source(nxy, domain, shape='diamond')[source]

Creates a synthetic source field in the shape of a circle or diamond. Useful for testing purposes.

Parameters:

  • nxy (tuple) – Number of grid points in the x and y directions (nx, ny).

  • domain (tuple) – Physical dimensions of the domain (xmax, ymax).

  • shape (str) – Shape of the source field. Options are “circle” or “diamond”. Default is “diamond”.

Returns:

A 2D array representing the source field.

Return type:

numpy.ndarray

src.utils.point_measurement(f, g)[source]

Computes the convolution of two 2D arrays evaluated at a specific point.

Parameters:
Returns:

The result of the convolution at the specified point.

Return type:

float

src.utils.point_source(nxy, domain, src_pt)[source]

Generates a point source field in Fourier space and transforms it back to the spatial domain.

Parameters:

  • nxy (tuple) – Number of grid points in the x and y directions (nx, ny).

  • domain (tuple) – Physical dimensions of the domain (xmax, ymax).

  • src_pt (tuple) – Coordinates of the source point (xs, ys).

Returns:

A 2D array representing the point source field in the spatial domain.

Return type:

numpy.ndarray

src.utils.setup_logging()[source]

Set up logging configuration.