PALM
The model PALM is based on the non-hydrostatic, filtered, incompressible Navier-Stokes equations in Boussinesq-approximated form (an anelastic approximation is available as an option for simulating deep convection). By default, PALM has at least six prognostic quantities: the velocity components u, v, w on a Cartesian grid, the potential temperature θ, water vapor mixing ratio qv and possibly a passive scalar s. Furthermore, an additional equation is solved for either the subgrid-scale turbulent kinetic energy (SGS-TKE) e (LES mode, default) or the total turbulent kinetic energy (RANS mode, see PALM-4U components).
In the LES mode, the filtering process yields four subgrid-scale (SGS) covariance terms. These SGS terms are parametrized using a 1.5-order closure after Deardorff (1980). PALM uses the modified version of Moeng and Wyngaard (1988) and Saiki et al. (2000). The closure is based on the assumption that the energy transport by SGS eddies is proportional to the local gradients of the mean quantities.
The following sections outline the most important specifics of PALM. For a more detailed description, see PALM documentation and related publications (Raasch and Schröter, 2001; Maronga et al., 2015).
Discretization in time and space
The model domain in PALM is discretized in space using finite differences and equidistant horizontal grid spacings. The Arakawa staggered C-grid (Harlow and Welch, 1965; Arakawa and Lamb, 1977) is used, where scalar quantities are defined at the center of each grid volume, whereas velocity components are shifted by half a grid width in their respective direction so that they are defined at the edges of the grid volumes. By default, the advection terms in the prognostic equations are discretized using an upwind-biased 5th-order differencing scheme (Wicker and Skamarock, 2002) in combination with a 3rd-order Runge–Kutta time-stepping scheme after Williamson (1980).
Pressure solver
The Boussinesq approximation requires incompressibility of the flow, but the integration of the governing equations does not provide this feature. Divergence of the flow field is thus inherently produced. Hence, a predictor corrector method is used where an equation is solved for the modified perturbation pressure after every time step (e.g., Patrinos and Kistler, 1977). In case of cyclic lateral boundary conditions, the solution of the Poisson equation is achieved by using a direct fast Fourier transform (FFT). PALM provides the inefficient but less restrictive Singleton-FFT (Singleton, 1969) and the well optimized Temperton-FFT (Temperton, 1992). External FFT libraries can be used as well, with the FFTW (Frigo and Johnson, 1998) being the most efficient one. Alternatively, the iterative multigrid scheme can be used (e.g., Hackbusch, 1985). This scheme uses the Gauss–Seidel method for the inner iterations on each grid level.
Boundary conditions
PALM offers a variety of boundary conditions. Dirichlet or Neumann boundary conditions can be chosen for u, v, θ, qv, and p∗ at the bottom and top of the model. For the horizontal velocity components the choice of Neumann (Dirichlet) boundary conditions yields free-slip (no-slip) conditions. Neumann boundary conditions are also used for the SGS-TKE. Kinematic fluxes of heat and moisture can be prescribed at the surface instead (Neumann conditions) of temperature and humidity (Dirichlet conditions). At the top of the model, Dirichlet boundary conditions can be used with given values of the geostrophic wind. Vertical velocity is assumed to be zero at the surface and top boundaries, which implies using Neumann conditions for pressure.
Following Monin-Obukhov similarity theory (MOST) a constant flux layer can be assumed as boundary condition between the surface and the first grid level where scalars and horizontal velocities are defined. In PALM we assume that MOST can be also applied locally and we therefore calculate local fluxes, velocities, and scaling parameters. This scheme involves calculation of the Obukhov length L, which can be either done based on variables of the previous time step ("circular"), via a Newton iteration method, or via a look-up table for the stability parameter z/L.
Topography parameterization
The Cartesian topography (complex terrain and buildings) in PALM is generally based on the mask method (Briscolini and Santangelo, 1989) and allows for explicitly resolving solid obstacles such as buildings and orography. The implementation makes use of the following simplifications:
- the obstacle shape is approximated by (an appropriate number of) full grid cells to fit the grid, i.e., a grid cell is either 100% fluid or 100% obstacle,
- the obstacles are fixed (not moving).
Topography is realized in 3-D, e.g., overhanging structures as for example bridges, ceilings, or tunnels, are allowed, i.e. topography does not necessarily be surface-mounted. If no overhanging structures are present, the 3-D obstacle dimension reduces to a 2.5-D topography format, which is conform to the Digital Elevation Model (DEM) format (DEMs of city morphologies have become increasingly available worldwide due to advances in remote sensing technologies). In case of overhanging structures, however, 3-D topography information is required to mask obstacles and their faces in PALM.
The model domain is then separated into three subdomains:
- grid points in free fluid without adjacent surfaces, where the standard PALM code is executed,
- grid points next to surface that require extra code (e.g., surface parametrization), and
- grid points within obstacles, where the standard PALM code is executed but multiplied by zero.
Additional topography code is executed in grid volumes of subdomain B. The faces of the obstacles are always located where the respective surface-normal velocity components u, v, and w are defined so that the impermeability boundary condition can be implemented by setting the respective surface-normal velocity component to zero.
In case of 5th-order advection scheme, the numerical stencil at grid points adjacent to obstacles would require data which is located within the obstacle. In order to avoid this, the order of the advection scheme is successively degraded at respective grid volumes adjacent to obstacles, i.e., from the 5th-order to 3rd-order at the second grid point above/beside an obstacle and from 3rd-order to 1st-order at grid points directly adjacent to an obstacle.
Simulations with topography require the application of MOST between each surface and the first computational grid point outside of the topography. For vertical and horizontal downward-facing surfaces, neutral stratification is assumed for MOST.
Buildings are primarily realized as obstacles that react to the flow dynamics via form drag and friction forces by assuming a constant flux layer between the building surface and the adjacent air volume. A simple thermodynamic coupling is also possible by prescribing surface fluxes of sensible (and latent heat) at any of the building surface grid elements. Thermodynamic interactive buildings are realized within the PALM-4U components.
The technical realization of the topography and treatment of surface-bounded grid cells is be outlined in Section topography implementation.
Parallelization and scaling
The parallelization of the code is achieved by a 2-D domain decomposition method along the x and y direction with equally sized subdomains. Ghost layers are added at the side boundaries of the subdomains in order to account for the local data dependencies, which are caused by the need to compute finite differences at these positions. The number of ghost layers that are used in PALM depend on the order of the advection scheme, with three layers for the 5th-order Wicker-Skamarock scheme. Ghost layer data are exchanged after every time step. Data exchange between proecessor cores is realized using the Message Passing Interface (MPI). Additional loop vectorization via OpenMP is realized which also allows a so-called hybrid parallelization.
PALM shows excellent scaling which was tested for up to 50,000 processor cores (details).
External forcing and nesting
Usually, PALM is used to simulate the flow in the boundary layer which is a certain part of the atmosphere. Processes occurring on larger scales than those in the boundary layer including large scale advection of scalars, large scale pressure gradients or large scale subsidence have also to be considered in the model, especially when focusing on realistic situations observed during measurement campaigns. In limited domain models with non-cyclic boundary conditions the large scale state enters through the boundary conditions at the lateral walls, and is usually taken from larger-scale models. An additional possibility to account for tendencies in the LES model resulting from larger scales than those in the boundary layer is the usage of nudging. Nudging is a (Newtonian) relaxation technique which adjusts the large-eddy simulation to a given, larger scale flow situation (Anthes, 1974; Stauffer and Bao, 1993). Alternatively, a nesting system is available which allows for self-nesting of PALM, i.e. running a parent domain run at coarse grid resolution and large model domain with a nested child model domain run at finer resolution and smaller domain.
Ocean option and coupling to atmosphere
PALM allows for studying the ocean mixed layer (OML) by using an ocean option where the sea surface is defined at the top of the model, so that negative values of z indicate the depth.
A coupled mode for the atmospheric and oceanic versions of PALM has been developed in order to allow for studying the interaction between turbulent processes in the ABL and OML. The coupling is realized by the online exchange of information at the sea surface (boundary conditions) between two PALM runs (one atmosphere and one ocean). The atmospheric model uses a constant flux layer and transfers the kinematic surface fluxes of heat and moisture as well as the momentum fluxes to the oceanic model. Flux conservation between the ocean and the atmosphere requires an adjustment of the fluxes for the density of water.
Embedded models
PALM comes with an increasing number of embedded models:
- Cloud microphysics
- Lagrangian particle model (LPM)
- Lagrangian cloud model (LCM)
- Canopy model
- 1D model for precursor runs
- Land surface model
- Radiation models
- Wind turbine model
For more details, see PALM documentation.
Data output and handling
Due to the enormous amount of data that comes along with computationally expensive LES, the data handling plays a key role for the performance of LES models and for data analysis during post-processing. PALM is optimized to pursue the strategy of performing data operations to great extent online during the simulation instead of postpone these operations to the post-processing. In this way, the data output (e.g., of huge 4-D data, or temporal averages) can be significantly reduced. In order to allow the user to perform own calculations during runtime, the user interface offers a wide range of possibilities, e.g., for defining user-defined output quantities (see Sect. user interface).
PALM allows data output for different quantities as time series, (horizontally-averaged) vertical profiles, 2-D cross sections, 3-D volume data, and masked data (see Sect. user interface ). All data output files are in netCDF format, which can be processed by different public domain and commercial software. NetCDF data can also be easily read from Fortran programs, provided that a netCDF library is available. The netCDF libraries currently support three different binary formats for netCDF files: classic, 64-bit offset, and netCDF-4. The latter was introduced in netCDF version 4.0 and is based on the HDF5 (http://www.hdfgroup.org/HDF5) data format. PALM is able to handle all three netCDF formats and also supports parallel I/O for netCDF-4.