---
# System prepended metadata

title: '[Greenline] Cleanup stencils in ICON4py'
tags: [cycle 18]

---

# Cleanup stencils in ICON4py



- Shaped by: Magdalena, Nikki, Abishek
- Appetite (FTEs, weeks):1-2 devs, 2-3 weeks
- Developers: <!-- Filled in at the betting table unless someone is specifically required here -->
- Note: Anurag should be involved at least for the stencil renaming, if not other tasks. Please communicate with him during the project.

## Problem
<!-- The raw idea, a use case, or something we’ve seen that motivates us to work on this -->
The stencils in the dycore and 
were ported from ICON quite some time ago, since then gt4py has evolved and we also learned a lot about the dycore. This is not yet reflected in the stencil code. At the time we went through the Fortran code replacing one loop-nest after the other, and did not pay much attention as to what the code actually does. This resulted eventually in the same code being ported several times. 

## Solution
<!-- The core elements we came up with, presented in a form that’s easy for people to immediately understand -->
- remove duplicated code by reusing calling other field_operators
- use fused stencils in green line
- rename functions

### update to recent gt4py features
- Slicing of fields instead of different bounds in fused stencils (TODO havogt: if I understand this correctly it's a misunderstanding. Currently, you cannot slice fields in GT4Py.)

### use fused stencils where available
Depends on previouse project [from cycle 17](https://hackmd.io/J0jTgdbKROeb0P-VaUOVDg)
We want to use the fused stencil in the greenline granules as well:

- use fused stencils from project  where possible in greenline granules (diffusion, solve_nonhydro)


### Naming of stencils
The stencils were originally named after the Fortran module they appear in and numbered in order of appearance there. 

The diffusion stencils have be renamed from their original numbering to descriptive names, for the stencils in solve_nonhydro we collected a list of descriptive names at DWD offices in Offenbach.

Renaming needs to be done in `icon4py` and `icon-exclaim`:.

#### descriptive names:
##### velocity advection names:
- [x] mo_velocity_advection_stencil_01 - compute_tangential_wind
- [x] mo_velocity_advection_stencil_02 - interpolate_vn_to_ie_and_compute_ekin_on_edges
- [x] mo_velocity_advection_stencil_03 - interpolate_vt_to_interface_edges
- [x] mo_velocity_advection_stencil_04 - compute_contravariant_correction
- [x] mo_velocity_advection_stencil_05 - compute_horizontal_kinetic_energy
- [x] mo_velocity_advection_stencil_06 - extrapolate_at_top
- [x] mo_velocity_advection_stencil_07 - compute_horizontal_advection_term_for_vertical_velocity
- [x] mo_velocity_advection_stencil_08 - interpolate_to_cell_center
- [x] mo_velocity_advection_stencil_09 - interpolate_to_cell_center
- [x] mo_velocity_advection_stencil_10 - interpolate_to_half_levels_vp
- [x] mo_velocity_advection_stencil_11 - copy_cell_field_to_vp
- [x] mo_velocity_advection_stencil_12 - set_cell_kdim_field_to_zero_vp
- [x] mo_velocity_advection_stencil_13 - correct_contravariant_vertical_velocity
- [x] mo_velocity_advection_stencil_14 - compute_maximum_cfl_and_clip_contravariant_vertical_velocity
- [x] mo_velocity_advection_stencil_15 - interpolate_contravariant_vertical_verlocity_to_full_levels
- [x] mo_velocity_advection_stencil_16 - compute_advective_vertical_wind_tendency
- [x] mo_velocity_advection_stencil_17 - add_interpolated_horizontal_advection_of_w
- [x] mo_velocity_advection_stencil_18 - add_extra_diffusion_for_w_con_approaching_cfl 
- [x] mo_velocity_advection_stencil_19 - compute_advective_normal_wind_tendency
- [x] mo_velocity_advection_stencil_20 - add_extra_diffusion_for_wn_approaching_cfl

##### solve non hydro names:
- [x] mo_solve_nonhydro_stencil_01 - set_two_cell_kdim_fields_to_zero_vp
- [x] mo_solve_nonhydro_stencil_02 - extrapolate_temporally_exner_pressure (mentioned in dycore paper, numerical aspect, increases stability)
- [x] mo_solve_nonhydro_stencil_03 - set_cell_kdim_field_to_zero_vp
- [x] mo_solve_nonhydro_stencil_04 - interpolate_to_surface
- [x] mo_solve_nonhydro_stencil_05 - interpolate_to_half_levels_vp
- [x] mo_solve_nonhydro_stencil_06 - compute_first_vertical_derivative
- [x] mo_solve_nonhydro_stencil_07 - compute_pertubation_of_rho_and_theta
- [x] mo_solve_nonhydro_stencil_08 - compute_pertubation_of_rho_and_theta_and_rho_at_ic
- [x] mo_solve_nonhydro_stencil_09 - def compute_virtual_potential_temperatures_and_pressure_gradient
- [x] mo_solve_nonhydro_stencil_10 - compute_rho_virtual_potential_temperatures_and_pressure_gradient
- [x] mo_solve_nonhydro_stencil_11_lower - set_cell_kdim_field_to_zero_vp
- [x] mo_solve_nonhydro_stencil_11_upper - set_theta_v_prime_ic_at_lower_boundary
- [x] mo_solve_nonhydro_stencil_12 - compute_approx_of_2nd_vertical_derivative_of_exner
- [x] mo_solve_nonhydro_stencil_13 - compute_pertubation_of_rho_and_theta
- [x] mo_solve_nonhydro_stencil_14 - set_two_edge_kdim_fields_to_zero_wp
- [x] mo_solve_nonhydro_stencil_15 - set_two_edge_kdim_fields_to_zero_wp
- [x] mo_solve_nonhydro_stencil_16 - compute_horizontal_advection_of_rho_and_theta
- [x] mo_solve_nonhydro_stencil_17 - add_vertical_wind_derivative_to_divergence_damping
- [x] mo_solve_nonhydro_stencil_18 - compute_horizontal_gradient_of_exner_pressure_for_flat_coordinates
- [x] mo_solve_nonhydro_stencil_19 - compute_horizontal_gradient_of_exner_pressure_for_nonflat_coordinates
- [x] mo_solve_nonhydro_stencil_20 - compute_horizontal_gradient_of_exner_pressure_for_multiple_levels
- [x] mo_solve_nonhydro_stencil_21 - compute_hydrostatic_correction_term
- [x] mo_solve_nonhydro_stencil_22 - apply_hydrostatic_correction_to_horizontal_gradient_of_exner_pressure
- [x] mo_solve_nonhydro_stencil_23 - add_temporal_tendencies_to_vn_by_interpolating_between_time_levels
- [x] mo_solve_nonhydro_stencil_24 - add_temporal_tendencies_to_vn
- [x] mo_solve_nonhydro_stencil_25 - compute_graddiv2_of_vn
- [x] mo_solve_nonhydro_stencil_26 - apply_2nd_order_divergence_damping
- [x] mo_solve_nonhydro_stencil_27 - apply_weighted_2nd_and_4th_order_divergence_damping
- [x] mo_solve_nonhydro_stencil_28 - add_analysis_increments_to_vn
- [x] mo_solve_nonhydro_stencil_29 - compute_vn_on_lateral_boundary
- [x] mo_solve_nonhydro_stencil_30 - compute_avg_vn_and_graddiv_vn_and_vt
- [x] mo_solve_nonhydro_stencil_31 - compute_avg_vn
- [x] mo_solve_nonhydro_stencil_32 - compute_mass_flux
- [x] mo_solve_nonhydro_stencil_33 - set_two_edge_kdim_fields_to_zero_wp
- [x] mo_solve_nonhydro_stencil_34 - accumulate_prep_adv_fields
- [x] mo_solve_nonhydro_stencil_35 - compute_contravariant_correction
- [x] mo_solve_nonhydro_stencil_36 - interpolate_vn_and_vt_to_ie_and_compute_ekin_on_edges
- [x] mo_solve_nonhydro_stencil_37 - compute_horizontal_kinetic_energy
- [x] mo_solve_nonhydro_stencil_38 - extrapolate_at_top
- [x] mo_solve_nonhydro_stencil_39 - compute_contravariant_correction_of_w
- [x] mo_solve_nonhydro_stencil_40 - compute_contravariant_correction_of_w_for_lower_boundary
- [x] mo_solve_nonhydro_stencil_41 - compute_divergence_of_fluxes_of_rho_and_theta
- [x] mo_solve_nonhydro_stencil_42 - compute_explicit_vertical_wind_from_advection_and_vertical_wind_density
- [x] mo_solve_nonhydro_stencil_43 - compute_explicit_vertical_wind_speed_and_vertical_wind_times_density
- [x] mo_solve_nonhydro_stencil_44 - compute_solver_coefficients_matrix
- [x] mo_solve_nonhydro_stencil_45 - set_two_cell_kdim_fields_index_to_zero_vp
- [x] mo_solve_nonhydro_stencil_46 - set_two_cell_kdim_fields_to_zero_wp
- [x] mo_solve_nonhydro_stencil_47 - set_lower_boundary_condition_for_w_and_contravariant_correction
- [x] mo_solve_nonhydro_stencil_48 - compute_explicit_part_for_rho_and_exner
- [x] mo_solve_nonhydro_stencil_49 - compute_explicit_part_for_rho_and_exner
- [x] mo_solve_nonhydro_stencil_50 - add_analysis_increments_from_data_assimilation
- [ ] solve nonhydro stencil 51 - mo_solve_nonhydro_stencil_51 (isn't called in icon-exclaim)
- [x] mo_solve_nonhydro_stencil_52 - solve_tridiagonal_matrix_for_w_forward_sweep
- [x] mo_solve_nonhydro_stencil_53 - solve_tridiagonal_matrix_for_w_back_substitution
- [x] mo_solve_nonhydro_stencil_54 -  apply_rayleigh_damping_mechanism
- [x] mo_solve_nonhydro_stencil_55 - compute_results_for_thermodynamic_variables
- [x] mo_solve_nonhydro_stencil_56_63 - compute_dwdz_for_divergence_damping
- [x] mo_solve_nonhydro_stencil_57 - set_cell_kdim_field_to_zero_wp
- [x] mo_solve_nonhydro_stencil_58 - update_mass_flux
- [x] mo_solve_nonhydro_stencil_59 - copy_cell_kdim_field_to_vp
- [x] mo_solve_nonhydro_stencil_60 - update_dynamical_exner_time_increment
- [x] mo_solve_nonhydro_stencil_61 - update_density_exner_wind
- [x] mo_solve_nonhydro_stencil_62 - update_wind
- [x] mo_solve_nonhydro_stencil_64 - set_cell_kdim_field_to_zero_wp
- [x] mo_solve_nonhydro_stencil_65 - update_mass_flux_weighted
- [x] mo_solve_nonhydro_stencil_66 - compute_theta_and_exner
- [x] mo_solve_nonhydro_stencil_67 - compute_exner_from_rhotheta
- [x] mo_solve_nonhydro_stencil_68 - update_theta_v
- [x] mo_solve_nonhydro_4th_order_divdamp - apply_4th_order_divergence_damping


On the other hand people got used to the numbered namings and it reflects where the stencils is used in the original code.

Naming convention proposal: rename the stencils, keep old name with numbering in doc string.



### Duplicated stencils
There is quite some duplicated code in the stencils. Duplicates should be removed and code reuse enforced by composition of field operators. A short scan of 
the velocity advection stencils shows that: 


- [x] `mo_velocity_advection_stencil_01` is the same as the 3d line of `mo_solve_nonhydro_stencil_30`
```
vt = _mo_velocity_advection_stencil_01(...)
```

- [x] `mo_solve_nonhydro_stencil_36` reuses `mo_velocity_advection_stencil_02` and `mo_velocity_advection_stencil_03`
- [x] `mo_solve_nonhydro_stencil_35` is the same as `mo_velocity_advection_stencil_04`
- [x] `mo_velocity_advection_stencil_05` line 3 reuses line 1 from `mo_velocity_advection_stencil_02`
- [x] `mo_velocity_advection_stencil_06`is the same as `mo_solve_nonhydro_stencil_38`
- [x] `mo_velocity_advection_stencil_08`and `mo_velocity_advection_stencil_09` do the same thing on different fields
- [x] `mo_velocity_advection_stencil_11` just copies a Cell Field, why not call it that, copying fields is also a part of `mo_velocity_advection_stencil_05` (lines 1 and 2: *dh: both are EdgeFields, plus second line is purely `vp`?) `mo_solve_nonhydro_stencil_46`(line 1: it is a init, am I looking at the wrong thing?) But it can be the same as mo_solve_nonhydro_stencil_59*
- [x] `mo_velocity_advection_stencil_12` broadcasts 0 to a Cell, K field, why not call it that? ,  `mo_solve_nonhydro_stencil_01`, `mo_solve_nonhydro_stencil_03`, `mo_solve_nonhydro_stencil_11_lower`,`mo_solve_nonhydro_stencil_33`, `mo_solve_nonhydro_stencil_45`,`mo_solve_nonhydro_stencil_45_b`, `mo_solve_nonhydro_stencil_46`, `mo_solve_nonhydro_stencil_47` (line 2) , `mo_solve_nonhydro_stencil_57`, `mo_solve_nonhydro_stencil_64` (the last two 57 and 64 being completely identical)
- [x] `mo_velocity_advection_stencil_13` computes the difference of Cell, K fields., why not call it such? *dh: see list from above*
- [x] `mo_velocity_advection_stencil_15`computes the average between k layers k and k+1 *dh: see list from above*
- [ ] `mo_velocty_advection_stencil_17` adds the result of `mo_velocity_advection_stencil_08` to a field *dh: the 08 stencil returns `vp`, the addition in 17 is done in `wp`?*

- [x] `mo_solve_nonhydro_stencil_08` computes `z_rth_pr_1` and `z_rth_pr_2` in the same way as `mo_solve_nonhydro_stencil_07`. These stencils are fused 
- [x] `mo_solve_nonhydro_stencil_09` computes `z_theta_v_pr_ic` in the same way as `mo_solve_nonhydro_stencil_05` only with `z_rth_pr_2` instead of `z_exner_ex_pr` 
- [x] `mo_solve_nonhydro_stencil_13` is the same as `mo_solve_nonhydro_stencil_07`
- [ ] `mo_solve_nonhydro_stencil_19` has part of the expression the same as `mo_solve_nonhydro_stencil_18` *dh: 18 returns a vp but in the computation of 19 is a wp needed.*
- [x] `mo_solve_nonhydro_stencil_31` computes `z_vn_avg` exactly as `mo_solve_nonhydro_stencil_30`
- [x] `mo_solve_nonhydro_stencil_36`: z_vt_ie is computed in the same way as in `mo_velocity_advection_stencil_03`
- [x] `mo_solve_nonhydro_stencil_38` is the same as `mo_velocity_advection_stencil_06`
- [ ] `mo_solve_nonhydro_stencil_39` and `mo_solve_nonhydro_stencil_40` have some analogous parts *dh: some are analogous but majority is different, not sure that added complexity is worth it? Am I overseeing something?*
- [x] `mo_solve_nonhydro_stencil_59` only assigns an input field to another variable. I think that this can be fused with either another stencil
- [x] `mo_solve_nonhydro_stencil_61` and `mo_solve_nonhydro_stencil_62` compute `w_new` exactly the same way
- [ ] stencils 46 to 55 and `set_zero_c_k` are the same for both predictor and corrector. Perhaps we can create a function with this stencils

### Cleanup stencils
- [x] The compute order in fortran and the icon4py stencil is sometimes different (`a+c+b` instead of `a+b+c`). This should be checked and fixed. For example: mo_solve_nonhydro_stencil_24, mo_velocity_advection_stencil_09, mo_velocity_advection_stencil_19
- [x] Sometimes there is still `**2.0` used instead of `**2`. 

### consolidate switches in solve_nonhydro
There are a couple of switches in `SolveNonhydro` that can be removed. (list composed by Christoph and Anurag):

The following dycore switches can be kept/ignored:
* [x] only support idiv_method 1, delete control flow (was removed in ICON)
* [ ] itime_scheme we keep flexible: requires porting of 1 additional stencil for itime_scheme 6
* [ ] divdamp_type: we support 32 and 3
* [ ] divdamp_order: we support 2, 4, and 24
* [ ] igradp_method we use 3, throw everything else out
* [x] l_open_ubc: has been removed from ICON-NWP and from ICON 2.6.7-rc
* ~~[ ] nesting should be ported sometime in the future, for now remove all switches~~

## Appetite
~ 3 weeks, needs coordination with other projects (fusion and renaming) Can easily be split and done by several people.



## Rabbit holes
<!-- Details about the solution worth calling out to avoid problems -->
Don't get lost in too many details trying to squeeze out the last bit. A first round removing obvious duplication will already help a lot.
Needs coordination with other projects.
## No-gos
<!-- Anything specifically excluded from the concept: functionality or use cases we intentionally aren’t covering to fit the ## appetite or make the problem tractable -->
