8.2
general documentation
Atmospheric flows module

Table of Contents

Data files

When using the atmospheric module, a file called meteo may be added to a case's DATA directorys in order to provide vertical profiles of the main variables.

Atmospheric domain mesh requirements

An atmospheric mesh has the following specific features:

  • The boundary located at the top of the domain should be a plane. So, horizontal wind speed at a given altitude can be prescribed at the top face as an inlet boundary.
  • Cells may have very different sizes, from very small (near ground or buildings) to very large (near the top of domain or far from zone of interest).
  • Vertical resolution: from tiny cells (e.g. $\Delta$\upshape z = 1 m) near the ground to a few hundreds of meters at the top.
  • Horizontal resolution: from a few meters to hundreds of meters.
  • The length ratio between two adjacent cells (in each direction) should preferably be between $0.7$ and $1.3$.
  • The z axis represents the vertical axis.

A topography map can be used to generate a mesh. In this case, the preprocessor mode is particularly useful to check the quality of the mesh (run type Mesh quality criteria).

Atmospheric flow model and steady/unsteady algorithm

The GUI may be used to enable the atmospheric flow module and set up the following calculation parameters in the Thermophysical models - Calculation features page see fig:steady.

The atmospheric flow model

The user can choose one of the following atmospheric flow models:

  • Constant density: To simulate neutral atmosphere.
  • Dry atmosphere: To simulate dry, thermally-stratified atmospheric flows (enables Potential temperature as thermal model).
  • Humid atmosphere: To simulate thermally stratified atmospheric flows (air-water mixture) with phase changes (enables Liquid potential temperature as thermal model). The model is described in [6] and [5].

Allowed time-stepping options

  • The pseudo-steady time-stepping is usually chosen. It sets a time step variable in space and time. It can be selected if constant boundary conditions are used, and usually provides fastest and smoothest convergence.
  • The unsteady time-stepping algorithm must be used in presence of time varying boundary conditions or source terms (the time step can then be variable in time or constant).

Selection of atmospheric model

Selection of steady/unsteady flow algorithm
Warning
The following points have to be considered when setting the parameters described above:
  • The potential temperature thermal model and the liquid potential temperature one (see the paragraph Atmospheric main variables for the definition) requires that the vertical component of the gravity is set to $g_z=-9.81 m.s^{-2}$ ( $g_x=g_y=0 m.s^{-2}$), otherwise pressure and density won't be correctly computed.
  • As well, the use of scalar with drift for atmospheric dispersion requires the gravity to be set to $g_z=-9.81$ ( $g_x=g_y=0 m.s^{-2}$), even if the density is constant.

Physical properties

The specific heat value has to be set to the atmospheric value $C_{p}=1005 J/kg/K$.

Parameters Constant
density 		Dry atmosphere Humid atmosphere Explanation
pressure boundary
condition 		Neumann first
order 		Extrapolation Extrapolation In case of Extrapolaion,
the pressure gradient is assumed (and set) constant, whereas in case of
Neumann first order, the pressure gradient is assumed (and set) to zero.
Improved pressure no yes yes If yes, exact balance between the hydrostatic part of the pressure gradient and the gravity term $\rho$g is numerically ensured.
Gravity (gravity is assumed aligned with the z-axis) $g_z=0$ or $g_z=-9.81 m.s^{-2}$ (the latter is useful for scalar with drift) $g_z=-9.81 m.s^{-2}$ $g_z=-9.81 m.s^{-2}$
Thermal variable no potential temperature liquid potential temperature
Others variables no no total water content, droplets number

Boundary and initial conditions

The meteo file can be used to define initial conditions for the different fields and to set up the inlet boundary conditions. For the velocity field, code_saturne can automatically detect if the boundary is an inlet boundary or an outflow boundary, according to the wind speed components given in the meteo file with respect to the boundary face orientation. This is often used for the lateral boundaries of the atmospheric domain, especially if the profile is evolving in time. In the case of inlet flow, the data given in the meteo file will be used as the input data (Dirichlet boundary condition) for velocity, temperature, humidity and turbulent variables. In the case of outflow, a Neumann boundary condition is automatically imposed (except for the pressure). The unit of temperature in the meteo file is the degree Celsius whereas the unit in the GUI is the kelvin.

To be taken into account, the meteo file has to be selected in the GUI (Atmospheric flows page, see fig:meteo) and the check box on the side ticked. This file gives the profiles of prognostic atmospheric variables containing one or a list of time stamps. The file has to be put in the DATA directory.

Selection of the *meteo* file

An example of file meteo is given in the directory data/user/meteo. The file format has to be strictly respected. The horizontal coordinates are not used at the present time (except when boundary conditions are based on several meteorological vertical profiles) and the vertical profiles are defined with the altitude above sea level. The highest altitude of the profile should be above the top of the simulation domain and the lowest altitude of the profile should be below or equal to the lowest level of the simulation domain. The line at the end of the meteo file should not be empty.

If the boundary conditions are variable in time, the vertical profiles for the different time stamps have to be written sequentially in the meteo file.

You can also set the profiles of atmospheric variables directly in the GUI. The following boundary conditions can be selected in the GUI:

  • Inlet/Outlet is automatically calculated for lateral boundaries (e.g. North, West\textellipsis ) of the computational domain (see fig:inlet).
  • Inlet for the top of the domain (see fig:top).
  • Rough wall for building walls (see fig:walls) or for the ground (see fig:ground). The user has to enter the roughness length. In case of variable roughness length, the user has to provide the land use data and the association between the roughness length values and land use categories.

Selection of automatic inlet/ outlet for boundary conditions

Selection of the boundary condition for the top of the domain

Selection of the boundary condition for building walls

Selection of the boundary condition for the ground
Remarks
If a meteorological file is given, it is used by default to initialize the variables. If a meteorological file is not given, the user can use the standard code_saturne initial and boundary conditions set up but has to be aware that even small inconsistencies can create very large buoyancy forces and spurious circulations.

Boundary conditions based on several meteorological vertical profiles

In some cases, especially when outputs of a mesoscale model are used, you need to build input boundary conditions from several meteorological vertical wind profiles. Cressman interpolation is then used to create the boundary conditions. The following files need to be put in the \texttt{DATA} directory: \item All meteo files giving the different vertical profiles of prognostic variables (wind, temperature, turbulent kinetic energy and dissipation).

  • A file called imbrication_files_list.txt which is a list of the meteo files used.
  • A separate meteo file which is used for the initial conditions and to impose inlet boundary conditions for the variables for which Cressman interpolation is not used (for example: temperature, turbulent kinetic energy). This file must follow the rules indicated previously.

The following files should be put in the SRC directory:

  • The user source file cs_user_parameters.f90. In this file, set the cressman_flag of each variable, for which the Cressman interpolation should be enabled, to true.

User-defined functions

User-defined functions may be used when the graphical user interface is not sufficient to set up the calculation. We provide some examples of user file for atmospheric application:

Physical models

Atmospheric dispersion of pollutants

To simulate the atmospheric dispersion of pollutant, one first need to define the source(s) term(s). That is to say the location i.e. the list of cells or boundary faces, the total air flow, the emitted mass fraction of pollutant, the emission temperature and the speed with the associated turbulent parameters. The mass fraction of pollutant is simulated through a user added scalar that could be a scalar with drift if wanted (aerosols for example).

The simulations can be done using 2 different methods:

  • Prescribing a boundary condition code total imposed mass flux for some boundary faces using the cs_user_boundary_conditions.c user function.
  • Using a scalar source term. In this case, the air inflow is not taken into account. The user has to add an explicit part to the equations for the scalar through the cs_user_source_terms.c file. This is done by selecting the cells and adding the source term st_exp which equals to the air flux multiplied by the mass fraction, while the implicit part st_imp is set to zero.

With the first method, the same problem of sources interactions appears, and moreover standard Dirichlet conditions should not be used (use itypfb=i_convective_inlet and icodcl=13 instead) as the exact emission rate cannot be prescribed because the diffusive part (usually negligible) cannot be quantified. Additionally, it requires that the boundary faces of the emission are explicitly represented in the mesh.
Finally the second method does not take into account the jet effect of the emission and so must be used only if it is sure that the emission does not modify the flow.
Whatever solution is chosen, the mass conservation should be verified by using for example the cs_user_extra_operations-scalar_balance.c file.

Soil/atmosphere interaction model

This model is based on the force restore model ([9]). It takes into account heat and humidity exchanges between the ground and the atmosphere at daily scale and the time evolution of ground surface temperature and humidity. Surface temperature is calculated with a prognostic equation whereas a 2-layers model is used to compute surface humidity.

The parameter iatsoil in the file atini0.f90 needs to be equal to one to activate the model. Then, the source file solvar.f90 is used.

Three variables need to be initialized in the file atini0.f90: deep soil temperature, surface temperature and humidity.

The user needs to give the values of the model constants in the file solcat.f90: roughness length, albedo, emissivity...

In case of a 3D simulation domain, land use data has to be provided for the domain. Values of model constants for the land use categories have also to be provided.

Radiative model (1D)

The 1D-radiative model calculates the radiative exchange between different atmospheric layers and the surface radiative fluxes.

The radiative exchange is computed separately for two wave lengths intervals

  • Calculation in the infrared spectral domain (file rayir.f90)
  • Calculation in the spectral range of solar radiation (file rayso.f90)

This 1D-radiative model is needed if the soil/atmosphere interaction model is activated.

This model is activated if the parameter iatra1 is equal to one in the file cs_users_parameters.f90.

Atmospheric main variables

For more details on the topic of atmospheric boundary layers, see [18].

  • Definition of the potential temperature:

    \[ \theta =T\left(\frac{P}{P_{r}}\right)^{-\frac{R_{d}}{C_{p}}} \]

  • Definition of liquid potential temperature:

    \[ \theta_{l} = \theta \left( 1-\frac{L}{C_{p}T} q_{l} \right) \]

  • Definition of virtual temperature:

    \[ T_{v} = \left(1+0.61q\right)T \]

  • Gas law:

    \[ P = \rho \frac{R}{M_{d}}\left(1+0,61q\right)T \]

    with $R=R_{d} M_{d}$.
  • Hydrostatic state:

    \[ \frac{\partial P}{\partial z} = -\rho g \]

Constant name Symbol Values Unit
Gravity acceleration at sea level $g$ $9.81$ $m.s^{-2}$
Effective Molecular Mass for dry air $M_{d}$ $28.97$ $kg.kmol^{-1}$
Standard reference pressure $P_{r}$ $10^{5}$ $Pa$
Universal gas constant $R$ $8.3143$ $J.K^{-1}.mol$
Gas constant for dry air $R_{d}$ $287$ $J.kg^{-1}.K^{-1}$
Variable name Symbol
Specific heat capacity of dry air $C_{p}$
Atmospheric pressure $P$
Specific humidity $q$
Specific content for liquid water $q_{l}$
Temperature $T$
Virtual temperature $T_{v}$
Potential temperature $\theta$
Liquid potential temperature $\theta_{l}$
Latent heat of vaporization $L$
Density $\rho $
Altitude $z$

Recommendations

This part is a list of recommendations for atmospheric numerical simulations.

  • Enough probes at different vertical levels in the domain should be used to check the convergence of the calculation.
  • An inflow boundary condition at the top level of the domain should be set (symmetry and automatic inlet/outlet are not appropriate).
  • A Courant number too small or too big has to be avoided (see code_saturne Best Practice Guidelines). That is the reason why the variable time step in space and in time option is recommended for steady simulations when there are large differences of cell size inside the domain (which is generally the case for atmospheric simulations). With this option, it can be necessary to change the reference time step and the time step maximal increase (by default, the time step increase rate is $10$).

In some cases, results can be improved with the following modifications:

  • In some case, the turbulent eddy viscosity can drop to unrealistically low values (especially with $k-\varepsilon$ model in stable atmospheric condition). In those cases, it is suggested to put an artificial molecular viscosity around $0.1 m^{2}.s^{-1}$.
  • If the main direction of wind is parallel to the boundary of your computing domain, try to set symmetry boundary conditions for the lateral boundaries to avoid inflow and outflow on the same boundary zone (side of your domain). Another possibility is to use a cylindrical mesh.
  • To avoid inflow and outflow on the same boundary zone (side of your domain), avoid the case of vertical profile in the input data \texttt{meteo} file with changes of the sign of velocity of wind ( $V_x$ or/and $V_y$).