diff --git a/model/atmosphere/dycore/src/icon4py/model/atmosphere/dycore/solve_nonhydro.py b/model/atmosphere/dycore/src/icon4py/model/atmosphere/dycore/solve_nonhydro.py index 8b939ed6d6..12780f8554 100644 --- a/model/atmosphere/dycore/src/icon4py/model/atmosphere/dycore/solve_nonhydro.py +++ b/model/atmosphere/dycore/src/icon4py/model/atmosphere/dycore/solve_nonhydro.py @@ -506,6 +506,10 @@ def __init__( "iau_wgt_dyn": self._config.iau_wgt_dyn, "is_iau_active": self._config.is_iau_active, "limited_area": self._grid.limited_area, + "divdamp_order": gtx.int32(self._config.divdamp_order), + "mean_cell_area": self._grid.global_properties.mean_cell_area, + "max_nudging_coefficient": self._config.max_nudging_coefficient, + "dbl_eps": constants.DBL_EPS, }, variants={ "apply_2nd_order_divergence_damping": [False, True], @@ -684,16 +688,6 @@ def __init__( "vertical_end": gtx.int32(self._grid.num_levels), }, ) - self._calculate_divdamp_fields = setup_program( - backend=backend, - program=dycore_utils.calculate_divdamp_fields, - constant_args={ - "divdamp_order": gtx.int32(self._config.divdamp_order), - "mean_cell_area": self._grid.global_properties.mean_cell_area, - "max_nudging_coefficient": self._config.max_nudging_coefficient, - "dbl_eps": constants.DBL_EPS, - }, - ) self._compute_rayleigh_damping_factor = setup_program( backend=backend, program=dycore_utils.compute_rayleigh_damping_factor, @@ -947,18 +941,6 @@ def _allocate_local_fields(self, allocator: gtx_allocators.FieldBufferAllocation self._grid, dims.KDim, dtype=ta.wpfloat, allocator=allocator ) """ - Declared as enh_divdamp_fac in ICON. - """ - self.reduced_fourth_order_divdamp_coeff_at_nest_boundary = data_alloc.zero_field( - self._grid, dims.KDim, dtype=ta.wpfloat, allocator=allocator - ) - """ - Declared as bdy_divdamp in ICON. - """ - self.fourth_order_divdamp_scaling_coeff = data_alloc.zero_field( - self._grid, dims.KDim, dtype=ta.wpfloat, allocator=allocator - ) - """ Declared as scal_divdamp in ICON. """ self.intermediate_fields = IntermediateFields.allocate(grid=self._grid, allocator=allocator) @@ -1296,13 +1278,6 @@ def run_corrector_step( second_order_divdamp_factor * self._grid.global_properties.mean_cell_area ) - self._calculate_divdamp_fields( - interpolated_fourth_order_divdamp_factor=self.interpolated_fourth_order_divdamp_factor, - fourth_order_divdamp_scaling_coeff=self.fourth_order_divdamp_scaling_coeff, - reduced_fourth_order_divdamp_coeff_at_nest_boundary=self.reduced_fourth_order_divdamp_coeff_at_nest_boundary, - second_order_divdamp_factor=second_order_divdamp_factor, - ) - log.debug("corrector run velocity advection") self.velocity_advection.run_corrector_step( diagnostic_state=diagnostic_state_nh, @@ -1352,8 +1327,8 @@ def run_corrector_step( normal_wind_iau_increment=diagnostic_state_nh.normal_wind_iau_increment, theta_v_at_edges_on_model_levels=z_fields.theta_v_at_edges_on_model_levels, horizontal_pressure_gradient=z_fields.horizontal_pressure_gradient, - reduced_fourth_order_divdamp_coeff_at_nest_boundary=self.reduced_fourth_order_divdamp_coeff_at_nest_boundary, - fourth_order_divdamp_scaling_coeff=self.fourth_order_divdamp_scaling_coeff, + interpolated_fourth_order_divdamp_factor=self.interpolated_fourth_order_divdamp_factor, + second_order_divdamp_factor=second_order_divdamp_factor, second_order_divdamp_scaling_coeff=second_order_divdamp_scaling_coeff, dtime=dtime, apply_2nd_order_divergence_damping=apply_2nd_order_divergence_damping, diff --git a/model/atmosphere/dycore/src/icon4py/model/atmosphere/dycore/solve_nonhydro_stencils.py b/model/atmosphere/dycore/src/icon4py/model/atmosphere/dycore/solve_nonhydro_stencils.py index bdb18e294c..2f9bf6efcd 100644 --- a/model/atmosphere/dycore/src/icon4py/model/atmosphere/dycore/solve_nonhydro_stencils.py +++ b/model/atmosphere/dycore/src/icon4py/model/atmosphere/dycore/solve_nonhydro_stencils.py @@ -6,33 +6,13 @@ # Please, refer to the LICENSE file in the root directory. # SPDX-License-Identifier: BSD-3-Clause import gt4py.next as gtx -from gt4py.next.experimental import concat_where from icon4py.model.atmosphere.dycore.dycore_utils import ( _broadcast_zero_to_three_edge_kdim_fields_wp, ) -from icon4py.model.atmosphere.dycore.stencils.compute_contravariant_correction import ( - _compute_contravariant_correction, -) -from icon4py.model.atmosphere.dycore.stencils.compute_horizontal_kinetic_energy import ( - _compute_horizontal_kinetic_energy, -) -from icon4py.model.atmosphere.dycore.stencils.compute_perturbation_of_rho_and_theta import ( - _compute_perturbation_of_rho_and_theta, -) -from icon4py.model.atmosphere.dycore.stencils.compute_perturbation_of_rho_and_theta_and_rho_interface_cell_centers import ( - _compute_perturbation_of_rho_and_theta_and_rho_interface_cell_centers, -) -from icon4py.model.atmosphere.dycore.stencils.compute_virtual_potential_temperatures_and_pressure_gradient import ( - _compute_virtual_potential_temperatures_and_pressure_gradient, -) -from icon4py.model.atmosphere.dycore.stencils.extrapolate_at_top import _extrapolate_at_top from icon4py.model.atmosphere.dycore.stencils.init_cell_kdim_field_with_zero_wp import ( _init_cell_kdim_field_with_zero_wp, ) -from icon4py.model.atmosphere.dycore.stencils.interpolate_vn_and_vt_to_ie_and_compute_ekin_on_edges import ( - _interpolate_vn_and_vt_to_ie_and_compute_ekin_on_edges, -) from icon4py.model.atmosphere.dycore.stencils.update_density_exner_wind import ( _update_density_exner_wind, ) @@ -63,165 +43,6 @@ def init_test_fields( ) -@gtx.field_operator -def _compute_pressure_gradient_and_perturbed_rho_and_potential_temperatures( - rho: fa.CellKField[float], - z_rth_pr_1: fa.CellKField[float], - z_rth_pr_2: fa.CellKField[float], - rho_ref_mc: fa.CellKField[float], - theta_v: fa.CellKField[float], - theta_ref_mc: fa.CellKField[float], - rho_ic: fa.CellKField[float], - wgtfac_c: fa.CellKField[float], - vwind_expl_wgt: fa.CellField[float], - exner_pr: fa.CellKField[float], - d_exner_dz_ref_ic: fa.CellKField[float], - ddqz_z_half: fa.CellKField[float], - z_theta_v_pr_ic: fa.CellKField[float], - theta_v_ic: fa.CellKField[float], - z_th_ddz_exner_c: fa.CellKField[float], -) -> tuple[ - fa.CellKField[float], - fa.CellKField[float], - fa.CellKField[float], - fa.CellKField[float], - fa.CellKField[float], - fa.CellKField[float], -]: - (z_rth_pr_1, z_rth_pr_2) = concat_where( - dims.KDim == 0, - _compute_perturbation_of_rho_and_theta(rho, rho_ref_mc, theta_v, theta_ref_mc), - (z_rth_pr_1, z_rth_pr_2), - ) - - (rho_ic, z_rth_pr_1, z_rth_pr_2) = concat_where( - dims.KDim >= 1, - _compute_perturbation_of_rho_and_theta_and_rho_interface_cell_centers( - wgtfac_c, rho, rho_ref_mc, theta_v, theta_ref_mc - ), - (rho_ic, z_rth_pr_1, z_rth_pr_2), - ) - - (z_theta_v_pr_ic, theta_v_ic, z_th_ddz_exner_c) = concat_where( - dims.KDim >= 1, - _compute_virtual_potential_temperatures_and_pressure_gradient( - wgtfac_c, - z_rth_pr_2, - theta_v, - vwind_expl_wgt, - exner_pr, - d_exner_dz_ref_ic, - ddqz_z_half, - ), - (z_theta_v_pr_ic, theta_v_ic, z_th_ddz_exner_c), - ) - - return z_rth_pr_1, z_rth_pr_2, rho_ic, z_theta_v_pr_ic, theta_v_ic, z_th_ddz_exner_c - - -@gtx.field_operator -def _predictor_stencils_35_36( - vn: fa.EdgeKField[float], - ddxn_z_full: fa.EdgeKField[float], - ddxt_z_full: fa.EdgeKField[float], - vt: fa.EdgeKField[float], - z_w_concorr_me: fa.EdgeKField[float], - wgtfac_e: fa.EdgeKField[float], - vn_ie: fa.EdgeKField[float], - z_vt_ie: fa.EdgeKField[float], - z_kin_hor_e: fa.EdgeKField[float], - k_field: fa.KField[gtx.int32], - nflatlev_startindex: gtx.int32, -) -> tuple[ - fa.EdgeKField[float], - fa.EdgeKField[float], - fa.EdgeKField[float], - fa.EdgeKField[float], -]: - z_w_concorr_me = concat_where( - dims.KDim >= nflatlev_startindex, - _compute_contravariant_correction(vn, ddxn_z_full, ddxt_z_full, vt), - z_w_concorr_me, - ) - (vn_ie, z_vt_ie, z_kin_hor_e) = concat_where( - dims.KDim >= 1, - _interpolate_vn_and_vt_to_ie_and_compute_ekin_on_edges(wgtfac_e, vn, vt), - (vn_ie, z_vt_ie, z_kin_hor_e), - ) - return z_w_concorr_me, vn_ie, z_vt_ie, z_kin_hor_e - - -@gtx.program(grid_type=gtx.GridType.UNSTRUCTURED) -def predictor_stencils_35_36( - vn: fa.EdgeKField[float], - ddxn_z_full: fa.EdgeKField[float], - ddxt_z_full: fa.EdgeKField[float], - vt: fa.EdgeKField[float], - z_w_concorr_me: fa.EdgeKField[float], - wgtfac_e: fa.EdgeKField[float], - vn_ie: fa.EdgeKField[float], - z_vt_ie: fa.EdgeKField[float], - z_kin_hor_e: fa.EdgeKField[float], - k_field: fa.KField[gtx.int32], - nflatlev_startindex: gtx.int32, - horizontal_start: gtx.int32, - horizontal_end: gtx.int32, - vertical_start: gtx.int32, - vertical_end: gtx.int32, -): - _predictor_stencils_35_36( - vn, - ddxn_z_full, - ddxt_z_full, - vt, - z_w_concorr_me, - wgtfac_e, - vn_ie, - z_vt_ie, - z_kin_hor_e, - k_field, - nflatlev_startindex, - out=(z_w_concorr_me, vn_ie, z_vt_ie, z_kin_hor_e), - domain={ - dims.EdgeDim: (horizontal_start, horizontal_end), - dims.KDim: (vertical_start, vertical_end), - }, - ) - - -@gtx.program(grid_type=gtx.GridType.UNSTRUCTURED) -def predictor_stencils_37_38( - vn: fa.EdgeKField[float], - vt: fa.EdgeKField[float], - vn_ie: fa.EdgeKField[float], - z_vt_ie: fa.EdgeKField[float], - z_kin_hor_e: fa.EdgeKField[float], - wgtfacq_e_dsl: fa.EdgeKField[float], - horizontal_start: gtx.int32, - horizontal_end: gtx.int32, - vertical_start: gtx.int32, - vertical_end: gtx.int32, -): - _compute_horizontal_kinetic_energy( - vn, - vt, - out=(vn_ie, z_vt_ie, z_kin_hor_e), - domain={ - dims.EdgeDim: (horizontal_start, horizontal_end), - dims.KDim: (vertical_start, vertical_start + 1), - }, - ) - _extrapolate_at_top( - vn, - wgtfacq_e_dsl, - out=vn_ie, - domain={ - dims.EdgeDim: (horizontal_start, horizontal_end), - dims.KDim: (vertical_end - 1, vertical_end), - }, - ) - - @gtx.program(grid_type=gtx.GridType.UNSTRUCTURED) def stencils_61_62( rho_now: fa.CellKField[float], diff --git a/model/atmosphere/dycore/src/icon4py/model/atmosphere/dycore/stencils/apply_4th_order_divergence_damping.py b/model/atmosphere/dycore/src/icon4py/model/atmosphere/dycore/stencils/apply_4th_order_divergence_damping.py index 79c86fa3a0..d468c9f4bc 100644 --- a/model/atmosphere/dycore/src/icon4py/model/atmosphere/dycore/stencils/apply_4th_order_divergence_damping.py +++ b/model/atmosphere/dycore/src/icon4py/model/atmosphere/dycore/stencils/apply_4th_order_divergence_damping.py @@ -6,8 +6,11 @@ # Please, refer to the LICENSE file in the root directory. # SPDX-License-Identifier: BSD-3-Clause import gt4py.next as gtx -from gt4py.next import astype, broadcast +from gt4py.next import astype +from icon4py.model.atmosphere.dycore.dycore_utils import ( + _calculate_fourth_order_divdamp_scaling_coeff, +) from icon4py.model.common import field_type_aliases as fa from icon4py.model.common.dimension import EdgeDim, KDim from icon4py.model.common.type_alias import vpfloat, wpfloat @@ -15,31 +18,45 @@ @gtx.field_operator def _apply_4th_order_divergence_damping( - scal_divdamp: fa.KField[wpfloat], + interpolated_fourth_order_divdamp_factor: fa.KField[wpfloat], z_graddiv2_vn: fa.EdgeKField[vpfloat], vn: fa.EdgeKField[wpfloat], + divdamp_order: gtx.int32, + mean_cell_area: float, + second_order_divdamp_factor: float, ) -> fa.EdgeKField[wpfloat]: - """Formelry known as _mo_solve_nonhydro_4th_order_divdamp.""" + """Formerly known as _mo_solve_nonhydro_4th_order_divdamp.""" + scal_divdamp = _calculate_fourth_order_divdamp_scaling_coeff( + interpolated_fourth_order_divdamp_factor, + divdamp_order, + mean_cell_area, + second_order_divdamp_factor, + ) z_graddiv2_vn_wp = astype(z_graddiv2_vn, wpfloat) - scal_divdamp = broadcast(scal_divdamp, (EdgeDim, KDim)) vn_wp = vn + (scal_divdamp * z_graddiv2_vn_wp) return vn_wp @gtx.program(grid_type=gtx.GridType.UNSTRUCTURED) def apply_4th_order_divergence_damping( - scal_divdamp: fa.KField[wpfloat], + interpolated_fourth_order_divdamp_factor: fa.KField[wpfloat], z_graddiv2_vn: fa.EdgeKField[vpfloat], vn: fa.EdgeKField[wpfloat], + divdamp_order: gtx.int32, + mean_cell_area: float, + second_order_divdamp_factor: float, horizontal_start: gtx.int32, horizontal_end: gtx.int32, vertical_start: gtx.int32, vertical_end: gtx.int32, ): _apply_4th_order_divergence_damping( - scal_divdamp, + interpolated_fourth_order_divdamp_factor, z_graddiv2_vn, vn, + divdamp_order, + mean_cell_area, + second_order_divdamp_factor, out=vn, domain={ EdgeDim: (horizontal_start, horizontal_end), diff --git a/model/atmosphere/dycore/src/icon4py/model/atmosphere/dycore/stencils/apply_weighted_2nd_and_4th_order_divergence_damping.py b/model/atmosphere/dycore/src/icon4py/model/atmosphere/dycore/stencils/apply_weighted_2nd_and_4th_order_divergence_damping.py index b8d69dadeb..f99353e001 100644 --- a/model/atmosphere/dycore/src/icon4py/model/atmosphere/dycore/stencils/apply_weighted_2nd_and_4th_order_divergence_damping.py +++ b/model/atmosphere/dycore/src/icon4py/model/atmosphere/dycore/stencils/apply_weighted_2nd_and_4th_order_divergence_damping.py @@ -6,47 +6,69 @@ # Please, refer to the LICENSE file in the root directory. # SPDX-License-Identifier: BSD-3-Clause import gt4py.next as gtx -from gt4py.next import astype, broadcast +from gt4py.next import astype +from icon4py.model.atmosphere.dycore.dycore_utils import ( + _calculate_fourth_order_divdamp_scaling_coeff, + _calculate_reduced_fourth_order_divdamp_coeff_at_nest_boundary, +) from icon4py.model.common import dimension as dims, field_type_aliases as fa from icon4py.model.common.type_alias import vpfloat, wpfloat @gtx.field_operator def _apply_weighted_2nd_and_4th_order_divergence_damping( - scal_divdamp: fa.KField[wpfloat], - bdy_divdamp: fa.KField[wpfloat], + interpolated_fourth_order_divdamp_factor: fa.KField[wpfloat], nudgecoeff_e: fa.EdgeField[wpfloat], z_graddiv2_vn: fa.EdgeKField[vpfloat], vn: fa.EdgeKField[wpfloat], + divdamp_order: gtx.int32, + mean_cell_area: float, + second_order_divdamp_factor: float, + max_nudging_coefficient: float, + dbl_eps: float, ) -> fa.EdgeKField[wpfloat]: - """Formelry known as _mo_solve_nonhydro_stencil_27.""" + """Formerly known as _mo_solve_nonhydro_stencil_27.""" + scal_divdamp = _calculate_fourth_order_divdamp_scaling_coeff( + interpolated_fourth_order_divdamp_factor, + divdamp_order, + mean_cell_area, + second_order_divdamp_factor, + ) + bdy_divdamp = _calculate_reduced_fourth_order_divdamp_coeff_at_nest_boundary( + scal_divdamp, max_nudging_coefficient, dbl_eps + ) z_graddiv2_vn_wp = astype(z_graddiv2_vn, wpfloat) - - scal_divdamp = broadcast(scal_divdamp, (dims.EdgeDim, dims.KDim)) - bdy_divdamp = broadcast(bdy_divdamp, (dims.EdgeDim, dims.KDim)) vn_wp = vn + (scal_divdamp + bdy_divdamp * nudgecoeff_e) * z_graddiv2_vn_wp return vn_wp @gtx.program(grid_type=gtx.GridType.UNSTRUCTURED) def apply_weighted_2nd_and_4th_order_divergence_damping( - scal_divdamp: fa.KField[wpfloat], - bdy_divdamp: fa.KField[wpfloat], + interpolated_fourth_order_divdamp_factor: fa.KField[wpfloat], nudgecoeff_e: fa.EdgeField[wpfloat], z_graddiv2_vn: fa.EdgeKField[vpfloat], vn: fa.EdgeKField[wpfloat], + divdamp_order: gtx.int32, + mean_cell_area: float, + second_order_divdamp_factor: float, + max_nudging_coefficient: float, + dbl_eps: float, horizontal_start: gtx.int32, horizontal_end: gtx.int32, vertical_start: gtx.int32, vertical_end: gtx.int32, ): _apply_weighted_2nd_and_4th_order_divergence_damping( - scal_divdamp, - bdy_divdamp, + interpolated_fourth_order_divdamp_factor, nudgecoeff_e, z_graddiv2_vn, vn, + divdamp_order, + mean_cell_area, + second_order_divdamp_factor, + max_nudging_coefficient, + dbl_eps, out=vn, domain={ dims.EdgeDim: (horizontal_start, horizontal_end), diff --git a/model/atmosphere/dycore/src/icon4py/model/atmosphere/dycore/stencils/compute_edge_diagnostics_for_dycore_and_update_vn.py b/model/atmosphere/dycore/src/icon4py/model/atmosphere/dycore/stencils/compute_edge_diagnostics_for_dycore_and_update_vn.py index 341de7eb15..a8e2e79385 100644 --- a/model/atmosphere/dycore/src/icon4py/model/atmosphere/dycore/stencils/compute_edge_diagnostics_for_dycore_and_update_vn.py +++ b/model/atmosphere/dycore/src/icon4py/model/atmosphere/dycore/stencils/compute_edge_diagnostics_for_dycore_and_update_vn.py @@ -280,8 +280,6 @@ def _apply_divergence_damping_and_update_vn( corrector_normal_wind_advective_tendency: fa.EdgeKField[ta.vpfloat], normal_wind_tendency_due_to_slow_physics_process: fa.EdgeKField[ta.vpfloat], normal_wind_iau_increment: fa.EdgeKField[ta.vpfloat], - reduced_fourth_order_divdamp_coeff_at_nest_boundary: fa.KField[ta.wpfloat], - fourth_order_divdamp_scaling_coeff: fa.KField[ta.wpfloat], second_order_divdamp_scaling_coeff: ta.wpfloat, theta_v_at_edges_on_model_levels: fa.EdgeKField[ta.wpfloat], horizontal_pressure_gradient: fa.EdgeKField[ta.vpfloat], @@ -290,6 +288,7 @@ def _apply_divergence_damping_and_update_vn( inv_dual_edge_length: fa.EdgeField[ta.wpfloat], nudgecoeff_e: fa.EdgeField[ta.wpfloat], geofac_grdiv: gtx.Field[[dims.EdgeDim, dims.E2C2EODim], ta.wpfloat], + interpolated_fourth_order_divdamp_factor: fa.KField[ta.wpfloat], advection_explicit_weight_parameter: ta.wpfloat, advection_implicit_weight_parameter: ta.wpfloat, dtime: ta.wpfloat, @@ -298,6 +297,11 @@ def _apply_divergence_damping_and_update_vn( limited_area: bool, apply_2nd_order_divergence_damping: bool, apply_4th_order_divergence_damping: bool, + divdamp_order: gtx.int32, + mean_cell_area: float, + second_order_divdamp_factor: float, + max_nudging_coefficient: float, + dbl_eps: float, ) -> fa.EdgeKField[ta.wpfloat]: # add dw/dz for divergence damping term. In ICON, this stencil starts from k = kstart_dd3d until k = nlev - 1. # Since scaling_factor_for_3d_divdamp is zero when k < kstart_dd3d, it is meaningless to execute computation @@ -337,17 +341,24 @@ def _apply_divergence_damping_and_update_vn( ) if limited_area: next_vn = _apply_weighted_2nd_and_4th_order_divergence_damping( - scal_divdamp=fourth_order_divdamp_scaling_coeff, - bdy_divdamp=reduced_fourth_order_divdamp_coeff_at_nest_boundary, + interpolated_fourth_order_divdamp_factor=interpolated_fourth_order_divdamp_factor, nudgecoeff_e=nudgecoeff_e, z_graddiv2_vn=squared_horizontal_gradient_of_total_divergence, vn=next_vn, + divdamp_order=divdamp_order, + mean_cell_area=mean_cell_area, + second_order_divdamp_factor=second_order_divdamp_factor, + max_nudging_coefficient=max_nudging_coefficient, + dbl_eps=dbl_eps, ) else: next_vn = _apply_4th_order_divergence_damping( - scal_divdamp=fourth_order_divdamp_scaling_coeff, + interpolated_fourth_order_divdamp_factor=interpolated_fourth_order_divdamp_factor, z_graddiv2_vn=squared_horizontal_gradient_of_total_divergence, vn=next_vn, + divdamp_order=divdamp_order, + mean_cell_area=mean_cell_area, + second_order_divdamp_factor=second_order_divdamp_factor, ) if is_iau_active: @@ -534,8 +545,6 @@ def apply_divergence_damping_and_update_vn( corrector_normal_wind_advective_tendency: fa.EdgeKField[ta.vpfloat], normal_wind_tendency_due_to_slow_physics_process: fa.EdgeKField[ta.vpfloat], normal_wind_iau_increment: fa.EdgeKField[ta.vpfloat], - reduced_fourth_order_divdamp_coeff_at_nest_boundary: fa.KField[ta.wpfloat], - fourth_order_divdamp_scaling_coeff: fa.KField[ta.wpfloat], second_order_divdamp_scaling_coeff: ta.wpfloat, theta_v_at_edges_on_model_levels: fa.EdgeKField[ta.wpfloat], horizontal_pressure_gradient: fa.EdgeKField[ta.vpfloat], @@ -544,6 +553,7 @@ def apply_divergence_damping_and_update_vn( inv_dual_edge_length: fa.EdgeField[ta.wpfloat], nudgecoeff_e: fa.EdgeField[ta.wpfloat], geofac_grdiv: gtx.Field[[dims.EdgeDim, dims.E2C2EODim], ta.wpfloat], + interpolated_fourth_order_divdamp_factor: fa.KField[ta.wpfloat], advection_explicit_weight_parameter: ta.wpfloat, advection_implicit_weight_parameter: ta.wpfloat, dtime: ta.wpfloat, @@ -552,6 +562,11 @@ def apply_divergence_damping_and_update_vn( limited_area: bool, apply_2nd_order_divergence_damping: bool, apply_4th_order_divergence_damping: bool, + divdamp_order: gtx.int32, + mean_cell_area: float, + second_order_divdamp_factor: float, + max_nudging_coefficient: float, + dbl_eps: float, horizontal_start: gtx.int32, horizontal_end: gtx.int32, vertical_start: gtx.int32, @@ -574,8 +589,6 @@ def apply_divergence_damping_and_update_vn( - corrector_normal_wind_advective_tendency: horizontal advection tendency of the normal wind at corrector step [m s-2] - normal_wind_tendency_due_to_slow_physics_process: normal wind tendeny due to slow physics [m s-2] - normal_wind_iau_increment: iau increment to normal wind (data assimilation) [m s-1] - - reduced_fourth_order_divdamp_coeff_at_nest_boundary: fourth order divergence damping coefficient at nest boundary [m2 s2] - - fourth_order_divdamp_scaling_coeff: fourth order divergence damping coefficient [m2 s2] - second_order_divdamp_scaling_coeff: second order divergence damping coefficient [m s] - theta_v_at_edges_on_model_levels: virtual potential temperature at edges on model levels [K] - horizontal_pressure_gradient: horizontal pressure gradient at edges on model levels [Pa m-1] @@ -584,18 +597,20 @@ def apply_divergence_damping_and_update_vn( - inv_dual_edge_length: inverse dual edge length - nudgecoeff_e: nudging coefficient for fourth order divergence damping at nest boundary - geofac_grdiv: metric coefficient for computation of horizontal gradient of divergence - - fourth_order_divdamp_factor: scaling factor for fourth order divergence damping - - second_order_divdamp_factor: scaling factor for second order divergence damping + - interpolated_fourth_order_divdamp_factor - advection_explicit_weight_parameter: explicitness weight of normal_wind_advective_tendency - advection_implicit_weight_parameter: implicitness weight of normal_wind_advective_tendency - dtime: time step [s] - iau_wgt_dyn: a scaling factor for iau increment - is_iau_active: option for iau increment analysis - - itime_scheme: ICON itime scheme (see ICON tutorial) - limited_area: option indicating the grid is limited area or not + - apply_2nd_order_divergence_damping: scaling factor for second order divergence damping + - apply_4th_order_divergence_damping: scaling factor for fourth order divergence damping - divdamp_order: divergence damping order (see the class DivergenceDampingOrder) - - start_edge_nudging_level_2: start index of second nudging level zone for edges - - end_edge_local: end index of local zone for edges + - mean_cell_area + - second_order_divdamp_factor + - max_nudging_coefficient + - dbl_eps Returns: - next_vn: normal wind to be updated [m s-1] @@ -611,14 +626,13 @@ def apply_divergence_damping_and_update_vn( normal_wind_iau_increment=normal_wind_iau_increment, theta_v_at_edges_on_model_levels=theta_v_at_edges_on_model_levels, horizontal_pressure_gradient=horizontal_pressure_gradient, - reduced_fourth_order_divdamp_coeff_at_nest_boundary=reduced_fourth_order_divdamp_coeff_at_nest_boundary, - fourth_order_divdamp_scaling_coeff=fourth_order_divdamp_scaling_coeff, second_order_divdamp_scaling_coeff=second_order_divdamp_scaling_coeff, horizontal_mask_for_3d_divdamp=horizontal_mask_for_3d_divdamp, scaling_factor_for_3d_divdamp=scaling_factor_for_3d_divdamp, inv_dual_edge_length=inv_dual_edge_length, nudgecoeff_e=nudgecoeff_e, geofac_grdiv=geofac_grdiv, + interpolated_fourth_order_divdamp_factor=interpolated_fourth_order_divdamp_factor, advection_explicit_weight_parameter=advection_explicit_weight_parameter, advection_implicit_weight_parameter=advection_implicit_weight_parameter, dtime=dtime, @@ -627,6 +641,11 @@ def apply_divergence_damping_and_update_vn( limited_area=limited_area, apply_2nd_order_divergence_damping=apply_2nd_order_divergence_damping, apply_4th_order_divergence_damping=apply_4th_order_divergence_damping, + divdamp_order=divdamp_order, + mean_cell_area=mean_cell_area, + second_order_divdamp_factor=second_order_divdamp_factor, + max_nudging_coefficient=max_nudging_coefficient, + dbl_eps=dbl_eps, out=next_vn, domain={ dims.EdgeDim: (horizontal_start, horizontal_end), diff --git a/model/atmosphere/dycore/tests/dycore/integration_tests/test_solve_nonhydro.py b/model/atmosphere/dycore/tests/dycore/integration_tests/test_solve_nonhydro.py index e2711f3c82..d0118054a1 100644 --- a/model/atmosphere/dycore/tests/dycore/integration_tests/test_solve_nonhydro.py +++ b/model/atmosphere/dycore/tests/dycore/integration_tests/test_solve_nonhydro.py @@ -585,16 +585,6 @@ def test_nonhydro_corrector_step( at_last_substep=at_last_substep, ) - if icon_grid.limited_area: - assert test_utils.dallclose( - solve_nonhydro.reduced_fourth_order_divdamp_coeff_at_nest_boundary.asnumpy(), - init_savepoint.bdy_divdamp().asnumpy(), - ) - - assert test_utils.dallclose( - solve_nonhydro.fourth_order_divdamp_scaling_coeff.asnumpy(), - init_savepoint.scal_divdamp().asnumpy(), - ) # stencil 10 assert test_utils.dallclose( diagnostic_state_nh.rho_at_cells_on_half_levels.asnumpy(), @@ -1625,20 +1615,18 @@ def test_apply_divergence_damping_and_update_vn( corrector_normal_wind_advective_tendency = sp_stencil_init.ddt_vn_apc_ntl(1) normal_wind_tendency_due_to_slow_physics_process = sp_stencil_init.ddt_vn_phy() normal_wind_iau_increment = sp_stencil_init.vn_incr() - reduced_fourth_order_divdamp_coeff_at_nest_boundary = sp_nh_init.bdy_divdamp() - fourth_order_divdamp_scaling_coeff = sp_nh_init.scal_divdamp() + interpolated_fourth_order_divdamp_factor = sp_nh_init.enh_divdamp_fac() theta_v_at_edges_on_model_levels = sp_stencil_init.z_theta_v_e() horizontal_pressure_gradient = sp_stencil_init.z_gradh_exner() current_vn = sp_stencil_init.vn() next_vn = savepoint_nonhydro_init.vn_new() horizontal_gradient_of_normal_wind_divergence = sp_nh_init.z_graddiv_vn() config = definitions.construct_nonhydrostatic_config(experiment) + mean_cell_area = grid_savepoint.mean_cell_area() iau_wgt_dyn = config.iau_wgt_dyn divdamp_order = config.divdamp_order - second_order_divdamp_scaling_coeff = ( - sp_nh_init.divdamp_fac_o2() * grid_savepoint.mean_cell_area() - ) + second_order_divdamp_scaling_coeff = sp_nh_init.divdamp_fac_o2() * mean_cell_area second_order_divdamp_factor = savepoint_nonhydro_init.divdamp_fac_o2() apply_2nd_order_divergence_damping = ( divdamp_order == dycore_states.DivergenceDampingOrder.COMBINED @@ -1668,8 +1656,6 @@ def test_apply_divergence_damping_and_update_vn( normal_wind_iau_increment=normal_wind_iau_increment, theta_v_at_edges_on_model_levels=theta_v_at_edges_on_model_levels, horizontal_pressure_gradient=horizontal_pressure_gradient, - reduced_fourth_order_divdamp_coeff_at_nest_boundary=reduced_fourth_order_divdamp_coeff_at_nest_boundary, - fourth_order_divdamp_scaling_coeff=fourth_order_divdamp_scaling_coeff, second_order_divdamp_scaling_coeff=second_order_divdamp_scaling_coeff, horizontal_mask_for_3d_divdamp=metrics_savepoint.hmask_dd3d(), scaling_factor_for_3d_divdamp=metrics_savepoint.scalfac_dd3d(), @@ -1684,6 +1670,12 @@ def test_apply_divergence_damping_and_update_vn( limited_area=grid_savepoint.get_metadata("limited_area").get("limited_area"), apply_2nd_order_divergence_damping=apply_2nd_order_divergence_damping, apply_4th_order_divergence_damping=apply_4th_order_divergence_damping, + interpolated_fourth_order_divdamp_factor=interpolated_fourth_order_divdamp_factor, + divdamp_order=divdamp_order, + mean_cell_area=mean_cell_area, + second_order_divdamp_factor=second_order_divdamp_factor, + max_nudging_coefficient=config.max_nudging_coefficient, + dbl_eps=constants.DBL_EPS, horizontal_start=start_edge_nudging_level_2, horizontal_end=end_edge_local, vertical_start=gtx.int32(0), diff --git a/model/atmosphere/dycore/tests/dycore/stencil_tests/test_apply_4th_order_divergence_damping.py b/model/atmosphere/dycore/tests/dycore/stencil_tests/test_apply_4th_order_divergence_damping.py index 82148b6461..b4ca1f1d56 100644 --- a/model/atmosphere/dycore/tests/dycore/stencil_tests/test_apply_4th_order_divergence_damping.py +++ b/model/atmosphere/dycore/tests/dycore/stencil_tests/test_apply_4th_order_divergence_damping.py @@ -14,13 +14,15 @@ from icon4py.model.atmosphere.dycore.stencils.apply_4th_order_divergence_damping import ( apply_4th_order_divergence_damping, ) -from icon4py.model.common import dimension as dims +from icon4py.model.common import constants, dimension as dims from icon4py.model.common.grid import base from icon4py.model.common.states import utils as state_utils from icon4py.model.common.type_alias import vpfloat, wpfloat -from icon4py.model.common.utils.data_allocation import random_field +from icon4py.model.common.utils import data_allocation as data_alloc from icon4py.model.testing.stencil_tests import StencilTest +from . import test_dycore_utils + def apply_4th_order_divergence_damping_numpy( scal_divdamp: np.ndarray, @@ -39,24 +41,40 @@ class TestApply4thOrderDivergenceDamping(StencilTest): @staticmethod def reference( connectivities: dict[gtx.Dimension, np.ndarray], - scal_divdamp: np.ndarray, + interpolated_fourth_order_divdamp_factor: np.ndarray, z_graddiv2_vn: np.ndarray, vn: np.ndarray, + divdamp_order: gtx.int32, + mean_cell_area: float, + second_order_divdamp_factor: float, **kwargs: Any, ) -> dict: + scal_divdamp = test_dycore_utils.fourth_order_divdamp_scaling_coeff_numpy( + interpolated_fourth_order_divdamp_factor, + divdamp_order, + second_order_divdamp_factor, + mean_cell_area, + ) vn = apply_4th_order_divergence_damping_numpy(scal_divdamp, z_graddiv2_vn, vn) return dict(vn=vn) @pytest.fixture def input_data(self, grid: base.Grid) -> dict[str, gtx.Field | state_utils.ScalarType]: - scal_divdamp = random_field(grid, dims.KDim, dtype=wpfloat) - z_graddiv2_vn = random_field(grid, dims.EdgeDim, dims.KDim, dtype=vpfloat) - vn = random_field(grid, dims.EdgeDim, dims.KDim, dtype=wpfloat) + interpolated_fourth_order_divdamp_factor = data_alloc.random_field(grid, dims.KDim) + z_graddiv2_vn = data_alloc.random_field(grid, dims.EdgeDim, dims.KDim, dtype=vpfloat) + vn = data_alloc.random_field(grid, dims.EdgeDim, dims.KDim, dtype=wpfloat) + + divdamp_order = 24 + mean_cell_area = 1000.0 + second_order_divdamp_factor = 3.0 return dict( - scal_divdamp=scal_divdamp, + interpolated_fourth_order_divdamp_factor=interpolated_fourth_order_divdamp_factor, z_graddiv2_vn=z_graddiv2_vn, vn=vn, + divdamp_order=divdamp_order, + mean_cell_area=mean_cell_area, + second_order_divdamp_factor=second_order_divdamp_factor, horizontal_start=0, horizontal_end=gtx.int32(grid.num_edges), vertical_start=0, diff --git a/model/atmosphere/dycore/tests/dycore/stencil_tests/test_apply_divergence_damping_and_update_vn.py b/model/atmosphere/dycore/tests/dycore/stencil_tests/test_apply_divergence_damping_and_update_vn.py index ba5d783c5d..94849c4b3b 100644 --- a/model/atmosphere/dycore/tests/dycore/stencil_tests/test_apply_divergence_damping_and_update_vn.py +++ b/model/atmosphere/dycore/tests/dycore/stencil_tests/test_apply_divergence_damping_and_update_vn.py @@ -20,6 +20,8 @@ from icon4py.model.common.grid import base, horizontal as h_grid from icon4py.model.common.utils import data_allocation as data_alloc +from . import test_dycore_utils + divergence_damp_order = DivergenceDampingOrder() @@ -60,22 +62,26 @@ def reference( normal_wind_iau_increment: np.ndarray, theta_v_at_edges_on_model_levels: np.ndarray, horizontal_pressure_gradient: np.ndarray, - reduced_fourth_order_divdamp_coeff_at_nest_boundary: np.ndarray, - fourth_order_divdamp_scaling_coeff: np.ndarray, second_order_divdamp_scaling_coeff: ta.wpfloat, horizontal_mask_for_3d_divdamp: np.ndarray, scaling_factor_for_3d_divdamp: np.ndarray, inv_dual_edge_length: np.ndarray, nudgecoeff_e: np.ndarray, geofac_grdiv: np.ndarray, + interpolated_fourth_order_divdamp_factor: np.ndarray, advection_explicit_weight_parameter: ta.wpfloat, advection_implicit_weight_parameter: ta.wpfloat, dtime: ta.wpfloat, iau_wgt_dyn: ta.wpfloat, - is_iau_active: gtx.int32, - limited_area: gtx.int32, + is_iau_active: bool, + limited_area: bool, apply_2nd_order_divergence_damping: bool, apply_4th_order_divergence_damping: bool, + divdamp_order: gtx.int32, + mean_cell_area: float, + second_order_divdamp_factor: float, + max_nudging_coefficient: float, + dbl_eps: float, horizontal_start: gtx.int32, horizontal_end: gtx.int32, vertical_start: gtx.int32, @@ -113,6 +119,14 @@ def reference( next_vn, ) + if apply_2nd_order_divergence_damping: + next_vn = np.where( + (horizontal_start <= horz_idx) & (horz_idx < horizontal_end), + next_vn + + (second_order_divdamp_scaling_coeff * horizontal_gradient_of_total_divergence), + next_vn, + ) + if apply_4th_order_divergence_damping: e2c2eO = connectivities[dims.E2C2EODim] # verified for e-10 @@ -129,16 +143,17 @@ def reference( ), np.zeros_like(horizontal_gradient_of_total_divergence), ) - - if apply_2nd_order_divergence_damping: - next_vn = np.where( - (horizontal_start <= horz_idx) & (horz_idx < horizontal_end), - next_vn - + (second_order_divdamp_scaling_coeff * horizontal_gradient_of_total_divergence), - next_vn, + fourth_order_divdamp_scaling_coeff = ( + test_dycore_utils.fourth_order_divdamp_scaling_coeff_numpy( + interpolated_fourth_order_divdamp_factor, + divdamp_order, + second_order_divdamp_factor, + mean_cell_area, + ) + ) + reduced_fourth_order_divdamp_coeff_at_nest_boundary = test_dycore_utils.calculate_reduced_fourth_order_divdamp_coeff_at_nest_boundary_numpy( + fourth_order_divdamp_scaling_coeff, max_nudging_coefficient ) - - if apply_4th_order_divergence_damping: if limited_area: next_vn = np.where( (horizontal_start <= horz_idx) & (horz_idx < horizontal_end), @@ -212,12 +227,13 @@ def input_data(self, request: pytest.FixtureRequest, grid: base.Grid) -> dict: theta_v_at_edges_on_model_levels = data_alloc.random_field(grid, dims.EdgeDim, dims.KDim) horizontal_pressure_gradient = data_alloc.random_field(grid, dims.EdgeDim, dims.KDim) geofac_grdiv = data_alloc.random_field(grid, dims.EdgeDim, dims.E2C2EODim) - fourth_order_divdamp_scaling_coeff = data_alloc.random_field(grid, dims.KDim) - reduced_fourth_order_divdamp_coeff_at_nest_boundary = data_alloc.random_field( - grid, dims.KDim - ) + interpolated_fourth_order_divdamp_factor = data_alloc.random_field(grid, dims.KDim) nudgecoeff_e = data_alloc.random_field(grid, dims.EdgeDim) + mean_cell_area = 1000.0 + max_nudging_coefficient = 0.3 + dbl_eps = constants.DBL_EPS + dtime = 0.9 advection_implicit_weight_parameter = 0.75 advection_explicit_weight_parameter = 0.25 @@ -254,14 +270,13 @@ def input_data(self, request: pytest.FixtureRequest, grid: base.Grid) -> dict: normal_wind_iau_increment=normal_wind_iau_increment, theta_v_at_edges_on_model_levels=theta_v_at_edges_on_model_levels, horizontal_pressure_gradient=horizontal_pressure_gradient, - reduced_fourth_order_divdamp_coeff_at_nest_boundary=reduced_fourth_order_divdamp_coeff_at_nest_boundary, - fourth_order_divdamp_scaling_coeff=fourth_order_divdamp_scaling_coeff, second_order_divdamp_scaling_coeff=second_order_divdamp_scaling_coeff, horizontal_mask_for_3d_divdamp=horizontal_mask_for_3d_divdamp, scaling_factor_for_3d_divdamp=scaling_factor_for_3d_divdamp, inv_dual_edge_length=inv_dual_edge_length, nudgecoeff_e=nudgecoeff_e, geofac_grdiv=geofac_grdiv, + interpolated_fourth_order_divdamp_factor=interpolated_fourth_order_divdamp_factor, advection_explicit_weight_parameter=advection_explicit_weight_parameter, advection_implicit_weight_parameter=advection_implicit_weight_parameter, dtime=dtime, @@ -270,6 +285,11 @@ def input_data(self, request: pytest.FixtureRequest, grid: base.Grid) -> dict: limited_area=limited_area, apply_2nd_order_divergence_damping=apply_2nd_order_divergence_damping, apply_4th_order_divergence_damping=apply_4th_order_divergence_damping, + divdamp_order=divdamp_order, + mean_cell_area=mean_cell_area, + second_order_divdamp_factor=second_order_divdamp_factor, + max_nudging_coefficient=max_nudging_coefficient, + dbl_eps=dbl_eps, horizontal_start=start_edge_nudging_level_2, horizontal_end=end_edge_local, vertical_start=0, diff --git a/model/atmosphere/dycore/tests/dycore/stencil_tests/test_apply_weighted_2nd_and_4th_order_divergence_damping.py b/model/atmosphere/dycore/tests/dycore/stencil_tests/test_apply_weighted_2nd_and_4th_order_divergence_damping.py index 31de7f151f..33505b461b 100644 --- a/model/atmosphere/dycore/tests/dycore/stencil_tests/test_apply_weighted_2nd_and_4th_order_divergence_damping.py +++ b/model/atmosphere/dycore/tests/dycore/stencil_tests/test_apply_weighted_2nd_and_4th_order_divergence_damping.py @@ -14,13 +14,15 @@ from icon4py.model.atmosphere.dycore.stencils.apply_weighted_2nd_and_4th_order_divergence_damping import ( apply_weighted_2nd_and_4th_order_divergence_damping, ) -from icon4py.model.common import dimension as dims +from icon4py.model.common import constants, dimension as dims from icon4py.model.common.grid import base from icon4py.model.common.states import utils as state_utils from icon4py.model.common.type_alias import vpfloat, wpfloat -from icon4py.model.common.utils.data_allocation import random_field +from icon4py.model.common.utils import data_allocation as data_alloc from icon4py.model.testing.stencil_tests import StencilTest +from . import test_dycore_utils + def apply_weighted_2nd_and_4th_order_divergence_damping_numpy( scal_divdamp: np.ndarray, @@ -41,13 +43,28 @@ class TestApplyWeighted2ndAnd4thOrderDivergenceDamping(StencilTest): @staticmethod def reference( connectivities: dict[gtx.Dimension, np.ndarray], - scal_divdamp: np.ndarray, - bdy_divdamp: np.ndarray, + interpolated_fourth_order_divdamp_factor: np.ndarray, nudgecoeff_e: np.ndarray, z_graddiv2_vn: np.ndarray, vn: np.ndarray, + divdamp_order: gtx.int32, + mean_cell_area: float, + second_order_divdamp_factor: float, + max_nudging_coefficient: float, + dbl_eps: float, **kwargs: Any, ) -> dict: + scal_divdamp = test_dycore_utils.fourth_order_divdamp_scaling_coeff_numpy( + interpolated_fourth_order_divdamp_factor, + divdamp_order, + second_order_divdamp_factor, + mean_cell_area, + ) + bdy_divdamp = ( + test_dycore_utils.calculate_reduced_fourth_order_divdamp_coeff_at_nest_boundary_numpy( + scal_divdamp, max_nudging_coefficient + ) + ) vn = apply_weighted_2nd_and_4th_order_divergence_damping_numpy( scal_divdamp, bdy_divdamp, @@ -59,18 +76,27 @@ def reference( @pytest.fixture def input_data(self, grid: base.Grid) -> dict[str, gtx.Field | state_utils.ScalarType]: - scal_divdamp = random_field(grid, dims.KDim, dtype=wpfloat) - bdy_divdamp = random_field(grid, dims.KDim, dtype=wpfloat) - nudgecoeff_e = random_field(grid, dims.EdgeDim, dtype=wpfloat) - z_graddiv2_vn = random_field(grid, dims.EdgeDim, dims.KDim, dtype=vpfloat) - vn = random_field(grid, dims.EdgeDim, dims.KDim, dtype=wpfloat) + interpolated_fourth_order_divdamp_factor = data_alloc.random_field(grid, dims.KDim) + nudgecoeff_e = data_alloc.random_field(grid, dims.EdgeDim, dtype=wpfloat) + z_graddiv2_vn = data_alloc.random_field(grid, dims.EdgeDim, dims.KDim, dtype=vpfloat) + vn = data_alloc.random_field(grid, dims.EdgeDim, dims.KDim, dtype=wpfloat) + + divdamp_order = 24 + mean_cell_area = 1000.0 + second_order_divdamp_factor = 3.0 + max_nudging_coefficient = 0.3 + dbl_eps = constants.DBL_EPS return dict( - scal_divdamp=scal_divdamp, - bdy_divdamp=bdy_divdamp, + interpolated_fourth_order_divdamp_factor=interpolated_fourth_order_divdamp_factor, nudgecoeff_e=nudgecoeff_e, z_graddiv2_vn=z_graddiv2_vn, vn=vn, + divdamp_order=divdamp_order, + mean_cell_area=mean_cell_area, + second_order_divdamp_factor=second_order_divdamp_factor, + max_nudging_coefficient=max_nudging_coefficient, + dbl_eps=dbl_eps, horizontal_start=0, horizontal_end=gtx.int32(grid.num_edges), vertical_start=0, diff --git a/model/atmosphere/dycore/tests/dycore/stencil_tests/test_dycore_utils.py b/model/atmosphere/dycore/tests/dycore/stencil_tests/test_dycore_utils.py index bde56aba25..cae87978ba 100644 --- a/model/atmosphere/dycore/tests/dycore/stencil_tests/test_dycore_utils.py +++ b/model/atmosphere/dycore/tests/dycore/stencil_tests/test_dycore_utils.py @@ -27,15 +27,15 @@ # TODO(): apply StencilTest structure to this test -def fourth_order_divdamp_scaling_coeff_for_order_24_numpy( - a: np.ndarray, factor: float, mean_cell_area: float +def fourth_order_divdamp_scaling_coeff_numpy( + a: np.ndarray, divdamp_order: int, factor: float, mean_cell_area: float ) -> np.ndarray: - a = np.maximum(0.0, a - 0.25 * factor) - return -a * mean_cell_area**2 + b = np.maximum(0.0, a - 0.25 * factor) if divdamp_order == 24 else np.full_like(a, factor) + return -b * mean_cell_area**2 def calculate_reduced_fourth_order_divdamp_coeff_at_nest_boundary_numpy( - coeff: float, field: np.ndarray + field: np.ndarray, coeff: float ) -> np.ndarray: return 0.75 / (coeff + constants.DBL_EPS) * np.abs(field) @@ -61,8 +61,9 @@ def test_calculate_fourth_order_divdamp_scaling_coeff_order_24( offset_provider={}, ) - ref = fourth_order_divdamp_scaling_coeff_for_order_24_numpy( + ref = fourth_order_divdamp_scaling_coeff_numpy( interpolated_fourth_order_divdamp_factor.asnumpy(), + divdamp_order, second_order_divdamp_factor, mean_cell_area, ) @@ -106,7 +107,7 @@ def test_calculate_reduced_fourth_order_divdamp_coeff_at_nest_boundary( assert test_utils.dallclose( out.asnumpy(), calculate_reduced_fourth_order_divdamp_coeff_at_nest_boundary_numpy( - coeff, fourth_order_divdamp_scaling_coeff.asnumpy() + fourth_order_divdamp_scaling_coeff.asnumpy(), coeff ), ) @@ -123,13 +124,16 @@ def test_calculate_divdamp_fields(backend: gtx_typing.Backend) -> None: second_order_divdamp_factor = 0.7 max_nudging_coefficient = 0.3 - scaled_ref = fourth_order_divdamp_scaling_coeff_for_order_24_numpy( - divdamp_field.asnumpy(), second_order_divdamp_factor, mean_cell_area + scaled_ref = fourth_order_divdamp_scaling_coeff_numpy( + divdamp_field.asnumpy(), + divdamp_order, + second_order_divdamp_factor, + mean_cell_area, ) reduced_fourth_order_divdamp_coeff_at_nest_boundary_ref = ( calculate_reduced_fourth_order_divdamp_coeff_at_nest_boundary_numpy( - max_nudging_coefficient, scaled_ref + scaled_ref, max_nudging_coefficient ) )