There have been several fixes / modifications in the model, to run properly in single precision:
In COSMO-Model Version 5.01, some work has been invested by MeteoSwiss to implement different variants of epsilons into TERRA:
When implementing the ICON version of TERRA, this work was lost, because the ICON version was based on an earlier COSMO version, which only had the eps_soil. This work has been resurrected now.
Erwan Brisson noticed large differences appearing between a double precision and a single precision simulation for the year 1999 over Europe, especially over mountainous areas and the beginning of spring. After replacing some appearances of eps_soil by eps_temp, where temperatures where checked against the freezing point, these differences vanished.
Now we also checked all other occurences of eps_soil and replaced them, where necessary, according to the work in src_soil_multlay.f90.
For single precision runs, this will change the results.
For double precision results, all three epsilons have the same value as eps_soil, so there will be no changes of results.
Introduced clipping to 0.0 for computation of cloud variables to avoid negative
values when running in single precision.
(by M. Baldauf; also implemented in 5.05b1)
The RTTOV-library (and therefore also the interface library radiance) works only in double precision. The RTTOV library is used for the following tasks in the COSMO-Model:
Now the synthetic satellite pictures (lsynsat) can also be computed when COSMO is running in single precision and RTTOV10 or RTTOV12 is used. But the other tasks, which are necessary during data assimilation, cannot be performed then.
Computing the LPI requires solving a Newton method. The implemented convergence criteria (ABS (f/f') < 1.0E-2) is not practicable for single precision. There are many cases where the Newton method does not converge then.
Nevertheless, the value of ABS(f) is going below 1.0E-6 in all these cases. Therefore the criterion ABS(f) < eps_for_null = 1.0E-6 has been implemented now.
There are some power functions with q_rain and q_graupel, where the argument can have very small values (< 1.0E-15), which lead to floating point exceptions with Cray cce (altough everything runs fine for GNU).
Such small values can occur, because all q-values are only clipped at the end of the time loop. The computations have been modified in a way, that the result is set to 0.0, if the q-values are below 1.0E-15, otherwise, the power function a**b is computed with EXP (b * LOG(a)).
Approach from Alla Yurova et al., 2014, https://doi.org/10.1002/2013WR014624
New namelist variable in /PHYCTL/:
| Group | Name | Meaning | Default | |
|---|---|---|---|---|
| /PHYCTL/ | itype_mire | NEW |
|
0 |
NOTE: The mire parameterization only works with soil moisture dependent heat conductivity, i.e. itype_heatcond ≥ 2!
Changes of results due to a modification in sfc_terra:
An expression (eai(i)/sai(i)) was replaced by a local variable.
Using the local variable instead of this expression in another computation
changes the order of execution, which gives a different rounding.
The turbulence modules have been aligned between COSMO and ICON and some cleanups have been done. The current versions are now implemented into ICON code.
A bug has been fixed in the COSMO turbulence interface for itype_vdif = 1:
Variables, which are computed in turbdiff (zvari, rhon) have not been saved to
pass them to the vertical diffusion in subroutine vertdiff, which is called after
the physics for itype_vdif = 1.
This has been fixed now by using full 3D variables to save all values.
The temperature tendencies from the microphysics are now also gathered in case of lgsp_first and given to the dynamics in src_runge_kutta.f90 for further treatment and vertical diffusion for itype_vdif = -1. A treatment of these tendencies for itype_vdif = 0 or = 1 has not been considered yet.
Add CLAW directives to enhance performance of some GPU ported parts: (in turb_utilities.f90)
Boundary fields are now only saved once to GPU.
The transfer to GPU is controlled by bd_list (TYPE bd_field) in src_input.f90 and
triggered in acc_utilities.f90 with:
GNSS delays are non-local observations which are integrated along the
whole signal path. However, the LETKF requires a position for each
observation. The namelist variable Href is used to estimate the ZTD/STD
position. In case of ZTDs the latitude and longitude of the observation
is copied from the GNSS station coordinates, the height above ground
in m is specified by Href. In case of STDs the position on the signal
path with h=Href is determined and used as the location of the
observation. This leads to a small horizontal displacement from the GNSS
station into the direction the the GNSS satellite.
For Href = 0.0 m both, ZTDs and STDs, are located at the GNSS station.
New namelist entry in /STD_OBS/:
| Group | Name | Meaning | Default | |
|---|---|---|---|---|
| /STD_OBS/ | Href | NEW | GNSS delay reference height above ground in [m] | 0.0 |
The STD operator checks the computed ZTDs/STDs and rejects observations which lead to invalid quantities. Furthermore, GNSS stations located far below or above the model orography can be rejected. In case of STDs it is also checked if parts of the signal path proceed below the model surface. The valid ranges can be specified in the /STD_OBS/ namelist. A detailed description of the namelist parameters can be found in the the User's Guide (from Version 5.06 on).
New namelist entries in /STD_OBS/:
| Group | Name | Meaning | Default | |
|---|---|---|---|---|
| /STD_OBS/ | ZTDminUse | NEW | minimum ZTD used for assimilation, delay in [m] | 0.0 |
| ZTDmaxUse | NEW | maximum ZTD used for assimilation, delay in [m] | 3.0 | |
| StatBelowSurface | NEW | GNSS station height below model surface in [m] | 500.0 | |
| StatAboveSurface | NEW | GNSS station height above model surface in [m] | 1000.0 | |
| StatBelowColumn | NEW | GNSS station height below any model grid column used for interpolation in [m] | 500.0 | |
The STD operator works with different weather models and cannot access model fields like the HHL field which defines the structure of the vertical grid nodes. Nevertheless it is desirable to distribute the supporting points on the signal path in a way which resembles the vertical model grid structure.
In the previous version of the STD operator one single scale height HScaleP was used to approximate the vertical grid structure. In the current version a much more flexible option is introduced. The vertical model grid levels (full levels in COSMO) are used to define several intervals which will be approximated in separate steps. For each interval the number of model levels and the desired number of supporting points needs to be specified. Based on these data individual scale heights are computed and used to distribute the supporting points. Up to 10 intervals can be defined.
Example, COSMO-DE
New namelist entries in /STD_OBS/:
| Group | Name | Meaning | Default | |
|---|---|---|---|---|
| /STD_OBS/ | NStepOpt | NEW |
Option for vertical distribution of supporting points along the signal path
(integer); valid options are:
|
2 |
| Heights | NEW | Vertical model grid levels defining the upper bounds of the height intervals (REAL array of size 10); only the desired number of intervals needs to be specified. | 2829.47, 6773.21, 12877.68, 21500.0 | |
| Hlevel | NEW | Number of model levels within the intervals (INTEGER array of size 10) | 20, 10, 10, 10 | |
| Hpoints | NEW | Desired number of supporting points within each interval, should be greater than Hlevel (INTEGER array of size 10) | 21, 11, 11, 11 | |
If this option (lredn_repro) is set, redundancy is constrained to occur only between reports assigned to the same grid point. This limits the validity of the redundancy limits set in data_obs_cdfin.f90 and increases the risk that redundancy of the two reports, e.g. one from an original ASCII TAC report and the other one from an original BUFR report, is not detected and that both reports may be used actively.
If this option is not set (which coincides with the implementation in the previous code versions) then the results of the redundancy check may depend on the domain decomposition (which was a known bug). As a consequence of this, the fof feedback files from different ensemble members running on different processor configurations could then contain different sets of reports, and using these fof files in the LETKF as input could then lead to a fatal crash.
New namelist entry in /NUDGING/:
| Group | Name | Meaning | Default | |
|---|---|---|---|---|
| /NUDGING/ | lredn_repro | NEW | if .true., then reproducibility of redundancy check ensured irrespective of domain decomposition by allowing for redundancy only between reports assigned to the same grid point. | .FALSE. |
None for the use of conventional data, unless lredn_repro is switched on.
Changes in src_lheat_nudge.f90, data_lheat_nudge.f90 and organize_assimilation.f90.
Add new features to LHN approach, which have been proofen beneficial when testing LHN within ICON model
The new profile can be defined by the new namelist parameters: tt_clim_max, zlev_clim_max, std_clim_max.
To use the new definition, set the new namelist switch lhn_tclim_old=.FALSE.
Default values of
tt_clim_max = 0.00015_wp and
zlev_clim_max = 3500.0_wp,
std_clim_max = 4.0_wp
will get similar results as the old approach.
New namelist entry in /NUDGING/:
| Group | Name | Meaning | Default | |
|---|---|---|---|---|
| /NUDGING/ | lhn_tclim_only | NEW | If .TRUE., the old model climatological
heating profile is used, if no surrounding grid point with suffient
precipitation is found. If .FALSE., an artificial heating profile can be determined with 3 additional new namelist variables zlev_clim_max, tt_clim_max and std_clim_max. |
.TRUE. |
| lhn_artif_only | NEW | If .TRUE., apply only increments calculated with the climatological / artificial heating profile. | .FALSE. | |
| lhn_no_ttend | NEW | If .TRUE., only take into account increments in humidity, calculated when lhn_hum_adj=.TRUE.. | .FALSE. | |
| zlev_clim_max | NEW | Height of maximum value of atificial temperature profile. | 3500.0 | |
| tt_clim_max | NEW | Maximum value of atificial temperature profile. | 0.00015 | |
| std_clim_max | NEW | Parameter to define length of artifical temperature profile. | 4.0 | |
This new tuning can be activated by setting the namelist variable itype_diag_gusts = 5 (new option) in the group /DIACTL/.
This parameterization is an empirically trained statistical linear model with four predictors: the model mean wind at 10m and three different gust estimates from the Brasseur wind gust parameterization (the main gust estimate and the lower and upper bounds of the bounding interval). The coefficients in the linear model were trained using COSMO at a horizontal grid-spacing of 1.1 km over Central Europe. itype_diag_gust=5 can thus be expected to work best under these conditions but to some degree the parameterization should generalize to other grid-spacing and to other areas. This parameterization was created with the goal to overcome the systematic underestimation of strong wind gusts as it occurs with the simple model itype_diag_gusts=1.
There are two cases of extrapolation that the statistical model does not handle well because it was not trained under such conditions:
Using the lockfile mechanism is an alternative approach (compared to READY files) to wait for input files, which is used by MeteoSwiss. It has been implemented into INT2LM some time ago and now also in the COSMO-Model.
Files are not written directly to a file called filename, but first to a file with name filename.lck. Only if the file is fully written and is closed, the name is changed to filename. For this the (non-standard) Fortran intrinsic RENAME is used. But this intrinsic nowadays is offered by (nearly) all compilers.
Using lockfiles is activated by setting the namelist variable (in /IOCTL/) llockfiles = .TRUE..
| Group | Name | Meaning | Default | |
|---|---|---|---|---|
| /PHYCTL/ | itype_mire | NEW |
|
0 |
| /DIACTL/ | itype_diag_gusts | new option | 5: empirically trained statistical linear model. | 1 |
| /IOCTL/ | llockfiles | NEW | indicates whether to use lockfiles or not. | .FALSE. |
| /NUDGING/ | lredn_repro | NEW | if .true., then reproducibility of redundancy check ensured irrespective of domain decomposition by allowing for redundancy only between reports assigned to the same grid point. | .FALSE. |
| lhn_tclim_only | NEW | If .TRUE., the old model climatological
heating profile is used, if no surrounding grid point with suffient
precipitation is found. If .FALSE., an artificial heating profile can be determined with 3 additional new namelist variables zlev_clim_max, tt_clim_max and std_clim_max. |
.TRUE. | |
| lhn_artif_only | NEW | If .TRUE., apply only increments calculated with the climatological / artificial heating profile. | .FALSE. | |
| lhn_no_ttend | NEW | If .TRUE., only take into account increments in humidity, calculated when lhn_hum_adj=.TRUE.. | .FALSE. | |
| zlev_clim_max | NEW | Height of maximum value of atificial temperature profile. | 3500.0 | |
| tt_clim_max | NEW | Maximum value of atificial temperature profile. | 0.00015 | |
| std_clim_max | NEW | Parameter to define length of artifical temperature profile. | 4.0 | |
| /STD_OBS/ | Href | NEW | GNSS delay reference height above ground in [m] | 0.0 |
| ZTDminUse | NEW | minimum ZTD used for assimilation, delay in [m] | 0.0 | |
| ZTDmaxUse | NEW | maximum ZTD used for assimilation, delay in [m] | 3.0 | |
| StatBelowSurface | NEW | GNSS station height below model surface in [m] | 500.0 | |
| StatAboveSurface | NEW | GNSS station height above model surface in [m] | 1000.0 | |
| StatBelowColumn | NEW | GNSS station height below any model grid column used for interpolation in [m] | 500.0 | |
| NStepOpt | NEW |
Option for vertical distribution of supporting points along the signal path
(integer); valid options are:
|
2 | |
| Heights | NEW | Vertical model grid levels defining the upper bounds of the height intervals (REAL array of size 10); only the desired number of intervals needs to be specified. | 2829.47, 6773.21, 12877.68, 21500.0 | |
| Hlevel | NEW | Number of model levels within the intervals (INTEGER array of size 10) | 20, 10, 10, 10 | |
| Hpoints | NEW | Desired number of supporting points within each interval, should be greater than Hlevel (INTEGER array of size 10) | 21, 11, 11, 11 | |
If not activating the new namelist switches, there are only numerical changes: