Documentation of the Changes in the COSMO-Model
Version 4.24

22.06.2012

4.24 is another version towards the new unified version 5.0. It contains again some fixes for 4.23, but again a few new features. It contains some modifications that do influence the results of the simulations. Hopefully, these are the last changes which do change the results, so this version is a candidate for performing meteorological tests. If everything works out well, the next versions only contain technical changes with no effects on the results.

1. New Fast Waves Solver

(by Michael Baldauf)

Most of the operational setups of the COSMO model now use the so-called Runge-Kutta time integration scheme (Wicker, Skamarock, 2002). The basic idea of this time-splitting procedure is to treat the slow parts like advection or Coriolis force with a large time step, whereas the fast waves modes sound and gravity wave expansion are treated with a small time step.

As in the original proposal of Wicker, Skamarock (2002, 1998) the fast waves are treated horizontally with a backward-forward scheme and vertically implicit to allow larger values for the small time step. An additional filter process must stabilize this whole time-splitting procedure; usually a divergence damping is used (e.g. Skamarock, Klemp, 1992) in the fast waves solver.

Properties of the new fast waves solver compared to the current one are:

Improvement of the accuracy of all vertical derivatives and averages:
During the evaluation of the Runge-Kutta dynamical core for the COSMO-EU (7km)-setup it became obvious that the proper treatment of the vertical discretizations in strongly vertically stretched grids in particular improves the behaviour of the pressure bias. Now all explicit and implicit terms are discretized in a consistent way with appropriate weights to take into account the grid stretching.
Use of the divergence operator in strong conservation form:
This form of the divergence has a benefit in the calculation of the metric terms in vertically stretched grids, too, and in the formulation of the lower boundary condition for the vertical velocity w.
Isotropic treatment of the artificial divergence damping:
The current divergence damping considers only the horizontal gradients of the 3D-divergence in the momentum equations. This seems appropriate in the vicinity of the lower boundary, where the grid is highly anisotropic, too. But for smaller scale setups of COSMO (e.g. COSMO-DE with 2.8 km grid mesh size) the grid becomes more isotropic in the troposphere. Therefore an isotropic treatment of the artificial divergence damping seems to be reasonable. Due to high divergence damping Courant numbers near the ground, this process must be also treated vertically implicit like the sound and buoyancy terms. The stability of this treatment is shown in Baldauf (2010).
This isotropic treatment can be used optionally. By default it is switched off (internal switch l_3D_div_damping).
Also optionally a preliminary version of the Mahrer-discretization of the horizontal pressure gradient terms is available. By default it is switched off (internal switch l_hor_p_grad_Mahrer).

A new namelist switch itype_fast_waves has been introduced in /DYNCTL/, to activate the old (Default) or new fast-waves solver:

Group	Name	Meaning	Default
`/DYNCTL/`	`itype_fast_waves`	Type of fast waves solver for Runge-Kutta dynamics Old scheme (from module `fast_waves_rk.f90)` (Default) New scheme (from module `fast_waves_sc.f90)` with proper weightings for all vertical discretizations, divergence operator in strong conservation form, optionally: isotropic (fully 3D) divergence damping, optionally: Mahrer (1984) discretization	`1`

New options for the bottom boundary condition itype_bbc_w can be set which are only valid for itype_fast_waves=2. These will be elaborated further and will be documented in a later release.

Recommended Namelist settings for the new fast-waves solver in /DYNCTL/:

itype_fast_waves = 2
itype_bbc_w = 14
ldyn_bbc = .FALSE.

Back to Contents

2. Changes to Subroutines computing the Reference Atmospheres

(by Michael Baldauf)

There have been some changes to the subroutines reference_atmosphere, reference_atmosphere_2 and reference_atmosphere_BV in module meteo_utilities.f90:

They have been extended to also compute the temperature of the reference state t0hl on half levels.
The computation of the k-indices for special pressure levels has been removed from these subroutines and has been put to a new subroutine k_index_of_pressure_levels (also in meteo_utilities).
For irefatm=1 and using the new fast_waves_sc-solver, the pressure and temperature of the reference state on full levels are now computed analytically, instead of just averaging the values from the half levels.
For irefatm=2 (subroutine reference_atmosphere_2) the reference temperature on full levels, t0, is now calculated directly as a function of the height (instead of using p0), which is a bit more intuitive.
NOTE:
This slightly different calculation can change results of applications with irefatm=2 (e.g. COSMO_EU) due to round off errors.

The interfaces in src_input.f90 and src_artifdata.f90 have been adapted accordingly.

Back to Contents

3. Changes for sub-hourly Analysis Update with LETKF

(by Hendrik Reich)

When using the LETKF method for data assimilation in the future, there will be the need to have an assimilation cycle of only few minutes, and not of hours. Therefore the COSMO-Model has been adapted to read and write files, where also minutes and seconds are included, e.g.

laf20130504121500: for 04th of May 2013, 12 UTC + 15 Minutes and 00 seconds

instead of the hitherto existing laf2013050412.

The main changes are:

Subroutine make_fn from io_utilities:
The file names, that are created by this subroutine, can have 4 additional digits. This is controlled by the Namelist switch ydate_ini:
- If ydate_ini is specified with 10 digits (as it was up to now), also the file names are as they were before.
- If ydate_ini is specified with 14 digits (including minutes and seconds), the file names contain 4 additional digits. But even then old files will be recognized for input.
Subroutine get_utc_date: to also compute the minutes and seconds in a precise manner, the computation of the return parameter acthour has been modified, and therefore the results are slightly different.

NOTE:
Because the Radiation module uses this return parameter to calculate the zenith angle, also the results of the simulations are changing. But the differences are numerically negligible.

the date variables that are returned by the routine get_utc_date from module utilities have been extended:

Variable	Before	After
`yactdate1`	`2013050412`	`20130504121532`
`yactdate2`	`SAT 04.05.2013 12 UTC`	`SAT 04.05.2013 12:15:32 UTC`

This also changes some formats in the ASCII output files

All routines that use the date variables have been adapted.

Back to Contents

4. Problems in FLake Model with GRIB packing

(by Dmitrii Mironov)

There were problems with the FLake Model related to Grib packing. At the end of each time step the FLake Model checks all computed values for consistency. Then the variables are written to Grib and read again in the next assimilation cycle or at the start of a forecast. Due to grib packing the consistency of the data could be lost.

Therefore, a consistency check has now been added also at the beginning of a simulation in the subroutine flake_init.

Moreover, another security procedure has been implemented:
To limit the values of EXP-arguments, a maximum value for EXP arguments has been introduced in data_flake.f90: c_maxearg_flk = 1.0E+02_ireals.

NOTE: These modifications do change the results of the simulations!

Back to Contents

5. Technical Changes and Bug Fixes

io_utilities.f90:
- Reading restart files: Adaptations are necessary in check_input_grid to deal with binary restart files.
- Adaptations for reading NetCDF files: Unification of IDs with INT2LM: Changed dim_id for topo corrections from 11 to 15 (by Burkhardt Rockel).
data_io.f90:
Unification of NetCDF dimension ID's with INT2LM (by Burkhardt Rockel, CLM).
organize_data.f90:
- Corrected check for itype_albedo in line 1074: was erroneously checked two times for value = 2.
- SR input_ioctl: replaced some WARNINGs for async I/O by ERROR, because error status is set now.
src_input.f90:
- Reading restart files: Also check the possibility of irefatm=3 for restarts.
- Corrected an error message in SR fill_realarray.
- Adapted definition of dim_ids to INT2LM: changed ID for topo corrections from 11 to 15 (by Burkhardt Rockel).
src_output.f90: Conditional compilation for GRIBDWD in SR makegds:
The vertical coordinate parameters can only be written, if GRIBDWD is set and the Grib library is available, because of packing of REALs to INTEGERs. In case of restart these parameters are not written, but also not needed.
src_sat_rttov.f90:
Added a NOMOVE Compiler directive for the SX-compiler.
src_slow_tendencies_rk.f90:
Order of summations for humidity tracer qvtens in advection: was different than for all the other humidity tracers and has been adapted now
NOTE: This changes the results because of numerical reasons.

Back to Contents

6. Changes to the Namelists

There were the following changes for the Namelist variables:

Group	Name	Meaning	Default
`/RUNCTL/`	`ydate_ini`	The character string for this namelist switch can now be specified in two ways: `ydate_ini = 'yyyymmddhh'`: All file names are specified in the way it was before. `ydate_ini = 'yyyymmddhhmmss'`: (with additional digits for the minutes and seconds): All file names for the initial files or for `ytunit='d'` (in day format) are also specified with the additional digits.	`' '`
`/DYNCTL/`	`itype_fast_waves`	Type of fast waves solver for Runge-Kutta dynamics Old scheme (from module `fast_waves_rk.f90)` (Default) New scheme (from module `fast_waves_sc.f90)` with proper weightings for all vertical discretizations, divergence operator in strong conservation form, optionally: isotropic (fully 3D) divergence damping, optionally: Mahrer (1984) discretization	`1`

Back to Contents

7. Changes of Results

The following changes influence the results:

Computation of t0 in subroutine reference_atmosphere_2 in a slightly different way: Affects only applications that use irefatm=2
Use of precomputed t0 in module src_runge_kutta instead of re-computing (because of numerical reasons).
Changes in the FLake Model.
Changes to the routine get_utc_date and the return parameter acthour: Slight changes of the values because of numerical reasons. Because this parameter is used in the radiation scheme for determining the zenith angle of the sun, the values of the simulations can change.

Back to Contents