From 8e6733b28191993459659eb07864d4afca37e5b6 Mon Sep 17 00:00:00 2001 From: Sebastien Masson <sebastien.masson@locean.ipsl.fr> Date: Tue, 28 Jun 2022 09:06:17 +0000 Subject: [PATCH] dev_2021_rk3_halo into 68-summer-body-2022 --- cfgs/SHARED/field_def_nemo-oce.xml | 112 +++++++------- cfgs/SHARED/grid_def_nemo.xml | 33 ++-- src/OCE/DIA/diadetide.F90 | 21 ++- src/OCE/DIA/diahsb.F90 | 136 ++++++++++++----- src/OCE/DIA/diahth.F90 | 78 ++++++---- src/OCE/DIA/diamlr.F90 | 9 +- src/OCE/DIA/diaptr.F90 | 237 +++++++++++++++-------------- src/OCE/DIA/diawri.F90 | 4 +- src/OCE/IOM/iom.F90 | 36 +++-- src/OCE/IOM/iom_nf90.F90 | 37 ++--- src/OCE/IOM/prtctl.F90 | 44 ++++-- src/OCE/TRA/eosbn2.F90 | 14 +- src/OCE/TRA/traadv.F90 | 2 +- src/OCE/TRD/trddyn.F90 | 17 ++- src/OCE/TRD/trdglo.F90 | 111 +++++++------- src/OCE/TRD/trdken.F90 | 101 ++++++------ src/OCE/TRD/trdmxl.F90 | 115 +++++--------- src/OCE/TRD/trdmxl_oce.F90 | 48 +++--- src/OCE/TRD/trdpen.F90 | 36 +++-- src/OCE/TRD/trdtra.F90 | 122 +++++++-------- src/OCE/TRD/trdvor.F90 | 10 +- src/OCE/lib_fortran_generic.h90 | 4 +- 22 files changed, 705 insertions(+), 622 deletions(-) diff --git a/cfgs/SHARED/field_def_nemo-oce.xml b/cfgs/SHARED/field_def_nemo-oce.xml index e5bff5cc7..838cd2834 100644 --- a/cfgs/SHARED/field_def_nemo-oce.xml +++ b/cfgs/SHARED/field_def_nemo-oce.xml @@ -25,7 +25,7 @@ that are available in the tidal-forcing implementation (see --> <!-- Time --> - <field id="diamlr_time" grid_ref="diamlr_grid_T_2D" prec="8" /> + <field id="diamlr_time" grid_ref="diamlr_grid_T_2D_inner" prec="8" /> <!-- Regressors for tidal harmonic analysis --> <field id="diamlr_r001" field_ref="diamlr_time" expr="sin( __TDE_M2_omega__ * diamlr_time )" enabled=".TRUE." comment="harmonic:sin:M2" /> @@ -199,22 +199,22 @@ that are available in the tidal-forcing implementation (see <field id="wi_cff" long_name="Fraction of implicit vertical velocity" unit="#" grid_ref="grid_T_3D" /> <!-- next variables available with key_diahth --> - <field id="mlddzt" long_name="Thermocline Depth (depth of max dT/dz)" standard_name="depth_at_maximum_upward_derivative_of_sea_water_potential_temperature" unit="m" /> - <field id="mldr10_3" long_name="Mixed Layer Depth (dsigma = 0.03 wrt 10m)" standard_name="ocean_mixed_layer_thickness_defined_by_sigma_theta" unit="m" /> - <field id="mldr0_1" long_name="Mixed Layer Depth (dsigma = 0.01 wrt sfc)" standard_name="ocean_mixed_layer_thickness_defined_by_sigma_theta" unit="m" /> - <field id="mldr0_3" long_name="Mixed Layer Depth (dsigma = 0.03 wrt sfc)" standard_name="ocean_mixed_layer_thickness_defined_by_sigma_theta" unit="m" /> - <field id="mld_dt02" long_name="Mixed Layer Depth (|dT| = 0.2 wrt 10m)" standard_name="ocean_mixed_layer_thickness_defined_by_temperature" unit="m" /> - <field id="topthdep" long_name="Top of Thermocline Depth (dT = -0.2 wrt 10m)" standard_name="ocean_mixed_layer_thickness_defined_by_temperature" unit="m" /> - <field id="pycndep" long_name="Pycnocline Depth (dsigma[dT=-0.2] wrt 10m)" standard_name="ocean_mixed_layer_thickness_defined_by_sigma_theta" unit="m" /> - <field id="BLT" long_name="Barrier Layer Thickness" unit="m" > topthdep - pycndep </field> - <field id="tinv" long_name="Max of vertical invertion of temperature" unit="degC" /> - <field id="depti" long_name="Depth of max. vert. inv. of temperature" unit="m" /> - <field id="20d" long_name="Depth of 20C isotherm" standard_name="depth_of_isosurface_of_sea_water_potential_temperature" unit="m" axis_ref="iax_20C" /> - <field id="26d" long_name="Depth of 26C isotherm" standard_name="depth_of_isosurface_of_sea_water_potential_temperature" unit="m" axis_ref="iax_26C" /> - <field id="28d" long_name="Depth of 28C isotherm" standard_name="depth_of_isosurface_of_sea_water_potential_temperature" unit="m" axis_ref="iax_28C" /> - <field id="hc300" long_name="Heat content 0-300m" standard_name="integral_of_sea_water_potential_temperature_wrt_depth_expressed_as_heat_content" unit="J/m2" /> - <field id="hc700" long_name="Heat content 0-700m" standard_name="integral_of_sea_water_potential_temperature_wrt_depth_expressed_as_heat_content" unit="J/m2" /> - <field id="hc2000" long_name="Heat content 0-2000m" standard_name="integral_of_sea_water_potential_temperature_wrt_depth_expressed_as_heat_content" unit="J/m2" /> + <field id="mlddzt" long_name="Thermocline Depth (depth of max dT/dz)" standard_name="depth_at_maximum_upward_derivative_of_sea_water_potential_temperature" unit="m" grid_ref="grid_T_2D_inner" /> + <field id="mldr10_3" long_name="Mixed Layer Depth (dsigma = 0.03 wrt 10m)" standard_name="ocean_mixed_layer_thickness_defined_by_sigma_theta" unit="m" grid_ref="grid_T_2D_inner" /> + <field id="mldr0_1" long_name="Mixed Layer Depth (dsigma = 0.01 wrt sfc)" standard_name="ocean_mixed_layer_thickness_defined_by_sigma_theta" unit="m" grid_ref="grid_T_2D_inner" /> + <field id="mldr0_3" long_name="Mixed Layer Depth (dsigma = 0.03 wrt sfc)" standard_name="ocean_mixed_layer_thickness_defined_by_sigma_theta" unit="m" grid_ref="grid_T_2D_inner" /> + <field id="mld_dt02" long_name="Mixed Layer Depth (|dT| = 0.2 wrt 10m)" standard_name="ocean_mixed_layer_thickness_defined_by_temperature" unit="m" grid_ref="grid_T_2D_inner" /> + <field id="topthdep" long_name="Top of Thermocline Depth (dT = -0.2 wrt 10m)" standard_name="ocean_mixed_layer_thickness_defined_by_temperature" unit="m" grid_ref="grid_T_2D_inner" /> + <field id="pycndep" long_name="Pycnocline Depth (dsigma[dT=-0.2] wrt 10m)" standard_name="ocean_mixed_layer_thickness_defined_by_sigma_theta" unit="m" grid_ref="grid_T_2D_inner" /> + <field id="BLT" long_name="Barrier Layer Thickness" unit="m" grid_ref="grid_T_2D_inner" > topthdep - pycndep </field> + <field id="tinv" long_name="Max of vertical invertion of temperature" unit="degC" grid_ref="grid_T_2D_inner" /> + <field id="depti" long_name="Depth of max. vert. inv. of temperature" unit="m" grid_ref="grid_T_2D_inner" /> + <field id="20d" long_name="Depth of 20C isotherm" standard_name="depth_of_isosurface_of_sea_water_potential_temperature" unit="m" grid_ref="grid_T_2D_inner" axis_ref="iax_20C" /> + <field id="26d" long_name="Depth of 26C isotherm" standard_name="depth_of_isosurface_of_sea_water_potential_temperature" unit="m" grid_ref="grid_T_2D_inner" axis_ref="iax_26C" /> + <field id="28d" long_name="Depth of 28C isotherm" standard_name="depth_of_isosurface_of_sea_water_potential_temperature" unit="m" grid_ref="grid_T_2D_inner" axis_ref="iax_28C" /> + <field id="hc300" long_name="Heat content 0-300m" standard_name="integral_of_sea_water_potential_temperature_wrt_depth_expressed_as_heat_content" unit="J/m2" grid_ref="grid_T_2D_inner" /> + <field id="hc700" long_name="Heat content 0-700m" standard_name="integral_of_sea_water_potential_temperature_wrt_depth_expressed_as_heat_content" unit="J/m2" grid_ref="grid_T_2D_inner" /> + <field id="hc2000" long_name="Heat content 0-2000m" standard_name="integral_of_sea_water_potential_temperature_wrt_depth_expressed_as_heat_content" unit="J/m2" grid_ref="grid_T_2D_inner" /> <!-- variables available with diaar5 --> <field id="botpres" long_name="Sea Water Pressure at Sea Floor" standard_name="sea_water_pressure_at_sea_floor" unit="dbar" /> @@ -1018,7 +1018,7 @@ that are available in the tidal-forcing implementation (see </field_group> <!-- Total trends calculated every time step--> - <field_group id="trendT" grid_ref="grid_T_3D"> + <field_group id="trendT" grid_ref="grid_T_3D_inner"> <field id="ttrd_tot" long_name="temperature-trend: total model trend" unit="degC/s" /> <field id="strd_tot" long_name="salinity -trend: total model trend" unit="1e-3/s" /> <!-- Thickness weighted versions: --> @@ -1034,10 +1034,10 @@ that are available in the tidal-forcing implementation (see <field id="ketrd_spg" long_name="ke-trend: surface pressure gradient" unit="W/s^3" /> <field id="ketrd_spgexp" long_name="ke-trend: surface pressure gradient (explicit)" unit="W/s^3" /> <field id="ketrd_spgflt" long_name="ke-trend: surface pressure gradient (filter)" unit="W/s^3" /> - <field id="ssh_flt" long_name="filtered contribution to ssh (dynspg_flt)" unit="m" grid_ref="grid_T_2D" /> - <field id="w0" long_name="surface vertical velocity" unit="m/s" grid_ref="grid_T_2D" /> - <field id="pw0_exp" long_name="surface pressure flux due to ssh" unit="W/s^2" grid_ref="grid_T_2D" /> - <field id="pw0_flt" long_name="surface pressure flux due to filtered ssh" unit="W/s^2" grid_ref="grid_T_2D" /> + <field id="ssh_flt" long_name="filtered contribution to ssh (dynspg_flt)" unit="m" grid_ref="grid_T_2D_inner" /> + <field id="w0" long_name="surface vertical velocity" unit="m/s" grid_ref="grid_T_2D_inner" /> + <field id="pw0_exp" long_name="surface pressure flux due to ssh" unit="W/s^2" grid_ref="grid_T_2D_inner" /> + <field id="pw0_flt" long_name="surface pressure flux due to filtered ssh" unit="W/s^2" grid_ref="grid_T_2D_inner" /> <field id="ketrd_keg" long_name="ke-trend: KE gradient or hor. adv." unit="W/s^3" /> <field id="ketrd_rvo" long_name="ke-trend: relative vorticity or metric term" unit="W/s^3" /> <field id="ketrd_pvo" long_name="ke-trend: planetary vorticity" unit="W/s^3" /> @@ -1045,7 +1045,7 @@ that are available in the tidal-forcing implementation (see <field id="ketrd_udx" long_name="ke-trend: U.dx[U]" unit="W/s^3" /> <field id="ketrd_ldf" long_name="ke-trend: lateral diffusion" unit="W/s^3" /> <field id="ketrd_zdf" long_name="ke-trend: vertical diffusion" unit="W/s^3" /> - <field id="ketrd_tau" long_name="ke-trend: wind stress " unit="W/s^3" grid_ref="grid_T_2D" /> + <field id="ketrd_tau" long_name="ke-trend: wind stress " unit="W/s^3" grid_ref="grid_T_2D_inner" /> <field id="ketrd_bfr" long_name="ke-trend: bottom friction (explicit)" unit="W/s^3" /> <field id="ketrd_bfri" long_name="ke-trend: bottom friction (implicit)" unit="W/s^3" /> <field id="ketrd_atf" long_name="ke-trend: asselin time filter trend" unit="W/s^3" /> @@ -1062,7 +1062,7 @@ that are available in the tidal-forcing implementation (see <field id="petrd_xad" long_name="pe-trend: i-advection" unit="W/m^3" /> <field id="petrd_yad" long_name="pe-trend: j-advection" unit="W/m^3" /> <field id="petrd_zad" long_name="pe-trend: k-advection" unit="W/m^3" /> - <field id="petrd_sad" long_name="pe-trend: surface adv. (linssh true)" unit="W/m^3" grid_ref="grid_T_2D" /> + <field id="petrd_sad" long_name="pe-trend: surface adv. (linssh true)" unit="W/m^3" grid_ref="grid_T_2D_inner" /> <field id="petrd_ldf" long_name="pe-trend: lateral diffusion" unit="W/m^3" /> <field id="petrd_zdf" long_name="pe-trend: vertical diffusion" unit="W/m^3" /> <field id="petrd_zdfp" long_name="pe-trend: pure vert. diffusion" unit="W/m^3" /> @@ -1080,42 +1080,42 @@ that are available in the tidal-forcing implementation (see <field_group id="trendU" grid_ref="grid_U_3D"> <!-- variables available with ln_dyn_trd --> - <field id="utrd_hpg" long_name="i-trend: hydrostatic pressure gradient" unit="m/s^2" /> - <field id="utrd_spg" long_name="i-trend: surface pressure gradient" unit="m/s^2" /> - <field id="utrd_spgexp" long_name="i-trend: surface pressure gradient (explicit)" unit="m/s^2" /> - <field id="utrd_spgflt" long_name="i-trend: surface pressure gradient (filtered)" unit="m/s^2" /> - <field id="utrd_keg" long_name="i-trend: KE gradient or hor. adv." unit="m/s^2" /> - <field id="utrd_rvo" long_name="i-trend: relative vorticity or metric term" unit="m/s^2" /> - <field id="utrd_pvo" long_name="i-trend: planetary vorticity" unit="m/s^2" /> - <field id="utrd_zad" long_name="i-trend: vertical advection" unit="m/s^2" /> - <field id="utrd_udx" long_name="i-trend: U.dx[U]" unit="m/s^2" /> - <field id="utrd_ldf" long_name="i-trend: lateral diffusion" unit="m/s^2" /> - <field id="utrd_zdf" long_name="i-trend: vertical diffusion" unit="m/s^2" /> - <field id="utrd_tau" long_name="i-trend: wind stress " unit="m/s^2" grid_ref="grid_U_2D" /> - <field id="utrd_bfr" long_name="i-trend: bottom friction (explicit)" unit="m/s^2" /> - <field id="utrd_bfri" long_name="i-trend: bottom friction (implicit)" unit="m/s^2" /> - <field id="utrd_tot" long_name="i-trend: total momentum trend before atf" unit="m/s^2" /> - <field id="utrd_atf" long_name="i-trend: asselin time filter trend" unit="m/s^2" /> + <field id="utrd_hpg" long_name="i-trend: hydrostatic pressure gradient" unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="utrd_spg" long_name="i-trend: surface pressure gradient" unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="utrd_spgexp" long_name="i-trend: surface pressure gradient (explicit)" unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="utrd_spgflt" long_name="i-trend: surface pressure gradient (filtered)" unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="utrd_keg" long_name="i-trend: KE gradient or hor. adv." unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="utrd_rvo" long_name="i-trend: relative vorticity or metric term" unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="utrd_pvo" long_name="i-trend: planetary vorticity" unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="utrd_zad" long_name="i-trend: vertical advection" unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="utrd_udx" long_name="i-trend: U.dx[U]" unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="utrd_ldf" long_name="i-trend: lateral diffusion" unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="utrd_zdf" long_name="i-trend: vertical diffusion" unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="utrd_tau" long_name="i-trend: wind stress " unit="m/s^2" grid_ref="grid_U_2D_inner" /> + <field id="utrd_bfr" long_name="i-trend: bottom friction (explicit)" unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="utrd_bfri" long_name="i-trend: bottom friction (implicit)" unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="utrd_tot" long_name="i-trend: total momentum trend before atf" unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="utrd_atf" long_name="i-trend: asselin time filter trend" unit="m/s^2" grid_ref="grid_T_3D_inner" /> </field_group> <field_group id="trendV" grid_ref="grid_V_3D"> <!-- variables available with ln_dyn_trd --> - <field id="vtrd_hpg" long_name="j-trend: hydrostatic pressure gradient" unit="m/s^2" /> - <field id="vtrd_spg" long_name="j-trend: surface pressure gradient" unit="m/s^2" /> - <field id="vtrd_spgexp" long_name="j-trend: surface pressure gradient (explicit)" unit="m/s^2" /> - <field id="vtrd_spgflt" long_name="j-trend: surface pressure gradient (filtered)" unit="m/s^2" /> - <field id="vtrd_keg" long_name="j-trend: KE gradient or hor. adv." unit="m/s^2" /> - <field id="vtrd_rvo" long_name="j-trend: relative vorticity or metric term" unit="m/s^2" /> - <field id="vtrd_pvo" long_name="j-trend: planetary vorticity" unit="m/s^2" /> - <field id="vtrd_zad" long_name="j-trend: vertical advection" unit="m/s^2" /> - <field id="vtrd_vdy" long_name="i-trend: V.dx[V]" unit="m/s^2" /> - <field id="vtrd_ldf" long_name="j-trend: lateral diffusion" unit="m/s^2" /> - <field id="vtrd_zdf" long_name="j-trend: vertical diffusion" unit="m/s^2" /> - <field id="vtrd_tau" long_name="j-trend: wind stress " unit="m/s^2" grid_ref="grid_V_2D" /> - <field id="vtrd_bfr" long_name="j-trend: bottom friction (explicit)" unit="m/s^2" /> - <field id="vtrd_bfri" long_name="j-trend: bottom friction (implicit)" unit="m/s^2" /> - <field id="vtrd_tot" long_name="j-trend: total momentum trend before atf" unit="m/s^2" /> - <field id="vtrd_atf" long_name="j-trend: asselin time filter trend" unit="m/s^2" /> + <field id="vtrd_hpg" long_name="j-trend: hydrostatic pressure gradient" unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="vtrd_spg" long_name="j-trend: surface pressure gradient" unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="vtrd_spgexp" long_name="j-trend: surface pressure gradient (explicit)" unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="vtrd_spgflt" long_name="j-trend: surface pressure gradient (filtered)" unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="vtrd_keg" long_name="j-trend: KE gradient or hor. adv." unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="vtrd_rvo" long_name="j-trend: relative vorticity or metric term" unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="vtrd_pvo" long_name="j-trend: planetary vorticity" unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="vtrd_zad" long_name="j-trend: vertical advection" unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="vtrd_vdy" long_name="i-trend: V.dx[V]" unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="vtrd_ldf" long_name="j-trend: lateral diffusion" unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="vtrd_zdf" long_name="j-trend: vertical diffusion" unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="vtrd_tau" long_name="j-trend: wind stress " unit="m/s^2" grid_ref="grid_V_2D_inner" /> + <field id="vtrd_bfr" long_name="j-trend: bottom friction (explicit)" unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="vtrd_bfri" long_name="j-trend: bottom friction (implicit)" unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="vtrd_tot" long_name="j-trend: total momentum trend before atf" unit="m/s^2" grid_ref="grid_T_3D_inner" /> + <field id="vtrd_atf" long_name="j-trend: asselin time filter trend" unit="m/s^2" grid_ref="grid_T_3D_inner" /> </field_group> <!-- shared variables available with TOP interface --> diff --git a/cfgs/SHARED/grid_def_nemo.xml b/cfgs/SHARED/grid_def_nemo.xml index d66b314f3..374f67c65 100644 --- a/cfgs/SHARED/grid_def_nemo.xml +++ b/cfgs/SHARED/grid_def_nemo.xml @@ -336,41 +336,44 @@ <scalar /> </grid> <grid id="diamlr_grid_T_2D" > - <domain domain_ref="grid_T" /> + <domain domain_ref="grid_T_inner" /> <scalar /> </grid> + <grid id="diamlr_grid_T_2D_inner" > + <domain domain_ref="grid_T_inner" /> + </grid> <grid id="diamlr_grid_U_2D" > - <domain domain_ref="grid_U" /> + <domain domain_ref="grid_U_inner" /> <scalar /> </grid> <grid id="diamlr_grid_V_2D" > - <domain domain_ref="grid_V" /> + <domain domain_ref="grid_V_inner" /> <scalar /> </grid> <grid id="diamlr_grid_W_2D" > - <domain domain_ref="grid_W" /> + <domain domain_ref="grid_W_inner" /> <scalar /> </grid> <grid id="diamlr_grid_2D_to_grid_T_3D" > - <domain domain_ref="grid_T" /> + <domain domain_ref="grid_T_inner" /> <axis axis_ref="deptht"> <duplicate_scalar /> </axis> </grid> <grid id="diamlr_grid_2D_to_grid_U_3D" > - <domain domain_ref="grid_U" /> + <domain domain_ref="grid_U_inner" /> <axis axis_ref="depthu"> <duplicate_scalar /> </axis> </grid> <grid id="diamlr_grid_2D_to_grid_V_3D" > - <domain domain_ref="grid_V" /> + <domain domain_ref="grid_V_inner" /> <axis axis_ref="depthv"> <duplicate_scalar /> </axis> </grid> <grid id="diamlr_grid_2D_to_grid_W_3D" > - <domain domain_ref="grid_W" /> + <domain domain_ref="grid_W_inner" /> <axis axis_ref="depthw"> <duplicate_scalar /> </axis> @@ -383,37 +386,37 @@ </grid> <!-- grid definitions for the computation of daily detided model diagnostics (diadetide) --> <grid id="diadetide_grid_T_2D" > - <domain domain_ref="grid_T" /> + <domain domain_ref="grid_T_inner" /> <scalar /> </grid> <grid id="diadetide_grid_U_2D" > - <domain domain_ref="grid_U" /> + <domain domain_ref="grid_U_inner" /> <scalar /> </grid> <grid id="diadetide_grid_V_2D" > - <domain domain_ref="grid_V" /> + <domain domain_ref="grid_V_inner" /> <scalar /> </grid> <grid id="diadetide_grid_2D_to_grid_T_3D" > - <domain domain_ref="grid_T" /> + <domain domain_ref="grid_T_inner" /> <axis axis_ref="deptht"> <duplicate_scalar /> </axis> </grid> <grid id="diadetide_grid_2D_to_grid_U_3D" > - <domain domain_ref="grid_U" /> + <domain domain_ref="grid_U_inner" /> <axis axis_ref="depthu"> <duplicate_scalar /> </axis> </grid> <grid id="diadetide_grid_2D_to_grid_V_3D" > - <domain domain_ref="grid_V" /> + <domain domain_ref="grid_V_inner" /> <axis axis_ref="depthv"> <duplicate_scalar /> </axis> </grid> <grid id="diadetide_grid_2D_to_grid_W_3D" > - <domain domain_ref="grid_W" /> + <domain domain_ref="grid_W_inner" /> <axis axis_ref="depthw"> <duplicate_scalar /> </axis> diff --git a/src/OCE/DIA/diadetide.F90 b/src/OCE/DIA/diadetide.F90 index 9d6760039..43dcc23de 100644 --- a/src/OCE/DIA/diadetide.F90 +++ b/src/OCE/DIA/diadetide.F90 @@ -5,11 +5,11 @@ MODULE diadetide !!====================================================================== !! History : ! 2019 (S. Mueller) !!---------------------------------------------------------------------- - USE par_oce , ONLY : wp, jpi, jpj - USE in_out_manager , ONLY : lwp, numout - USE iom , ONLY : iom_put - USE dom_oce , ONLY : rn_Dt, nsec_day - USE phycst , ONLY : rpi + USE par_oce + USE in_out_manager + USE iom + USE dom_oce + USE phycst USE tide_mod #if defined key_xios USE xios @@ -24,6 +24,8 @@ MODULE diadetide PUBLIC :: dia_detide_init, dia_detide + !! * Substitutions +# include "do_loop_substitute.h90" !!---------------------------------------------------------------------- !! NEMO/OCE 4.0 , NEMO Consortium (2019) !! $Id$ @@ -90,9 +92,9 @@ CONTAINS !!---------------------------------------------------------------------- INTEGER, INTENT(in) :: kt - REAL(wp), DIMENSION(jpi,jpj) :: zwght_2D + REAL(wp), DIMENSION(A2D(0)) :: zwght_2D REAL(wp) :: zwght, ztmp - INTEGER :: jn + INTEGER :: ji, jj, jn ! Compute detiding weight at the current time-step; the daily total weight ! is one, and the daily summation of a diagnosed field multiplied by this @@ -104,7 +106,10 @@ CONTAINS zwght = zwght + 1.0_wp / REAL( ndiadetide, KIND=wp ) END IF END DO - zwght_2D(:,:) = zwght + + DO_2D( 0, 0, 0, 0 ) + zwght_2D(ji,jj) = zwght + END_2D CALL iom_put( "diadetide_weight", zwght_2D) END SUBROUTINE dia_detide diff --git a/src/OCE/DIA/diahsb.F90 b/src/OCE/DIA/diahsb.F90 index f9307c0af..1039d5c3c 100644 --- a/src/OCE/DIA/diahsb.F90 +++ b/src/OCE/DIA/diahsb.F90 @@ -50,6 +50,7 @@ MODULE diahsb REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: tmask_ini !! * Substitutions +# include "do_loop_substitute.h90" # include "domzgr_substitute.h90" !!---------------------------------------------------------------------- !! NEMO/OCE 4.0 , NEMO Consortium (2018) @@ -82,30 +83,61 @@ CONTAINS REAL(wp) :: z_frc_trd_v ! - - REAL(wp) :: z_wn_trd_t , z_wn_trd_s ! - - REAL(wp) :: z_ssh_hc , z_ssh_sc ! - - - REAL(wp), DIMENSION(jpi,jpj,13) :: ztmp - REAL(wp), DIMENSION(jpi,jpj,jpkm1,4) :: ztmpk + REAL(wp), DIMENSION(A2D(0),13) :: ztmp + REAL(wp), DIMENSION(A2D(0),jpkm1,4) :: ztmpk REAL(wp), DIMENSION(17) :: zbg !!--------------------------------------------------------------------------- IF( ln_timing ) CALL timing_start('dia_hsb') ! - ztmp (:,:,:) = 0._wp ! should be better coded - ztmpk(:,:,:,:) = 0._wp ! should be better coded - ! - ts(:,:,:,1,Kmm) = ts(:,:,:,1,Kmm) * tmask(:,:,:) ; ts(:,:,:,1,Kbb) = ts(:,:,:,1,Kbb) * tmask(:,:,:) ; - ts(:,:,:,2,Kmm) = ts(:,:,:,2,Kmm) * tmask(:,:,:) ; ts(:,:,:,2,Kbb) = ts(:,:,:,2,Kbb) * tmask(:,:,:) ; + DO_2D( 0, 0, 0, 0 ) + ztmp (ji,jj,:) = 0._wp ! should be better coded + ztmpk(ji,jj,:,:) = 0._wp ! should be better coded + ! + ts(ji,jj,:,1,Kmm) = ts(ji,jj,:,1,Kmm) * tmask(ji,jj,:) + ts(ji,jj,:,1,Kbb) = ts(ji,jj,:,1,Kbb) * tmask(ji,jj,:) + ! + ts(ji,jj,:,2,Kmm) = ts(ji,jj,:,2,Kmm) * tmask(ji,jj,:) + ts(ji,jj,:,2,Kbb) = ts(ji,jj,:,2,Kbb) * tmask(ji,jj,:) + END_2D ! ! ------------------------- ! ! 1 - Trends due to forcing ! ! ------------------------- ! ! prepare trends - ztmp(:,:,1) = - r1_rho0 * ( emp(:,:) - rnf(:,:) - fwfisf_cav(:,:) - fwfisf_par(:,:) ) * surf(:,:) ! volume - ztmp(:,:,2) = sbc_tsc(:,:,jp_tem) * surf(:,:) ! heat - ztmp(:,:,3) = sbc_tsc(:,:,jp_sal) * surf(:,:) ! salt - IF( ln_rnf ) ztmp(:,:,4) = rnf_tsc(:,:,jp_tem) * surf(:,:) ! runoff temp - IF( ln_rnf_sal ) ztmp(:,:,5) = rnf_tsc(:,:,jp_sal) * surf(:,:) ! runoff salt - IF( ln_isf ) ztmp(:,:,6) = ( risf_cav_tsc(:,:,jp_tem) + risf_par_tsc(:,:,jp_tem) ) * surf(:,:) ! isf temp - IF( ln_traqsr ) ztmp(:,:,7) = r1_rho0_rcp * qsr(:,:) * surf(:,:) ! penetrative solar radiation - IF( ln_trabbc ) ztmp(:,:,8) = qgh_trd0(:,:) * surf(:,:) ! geothermal heat + DO_2D( 0, 0, 0, 0 ) + ztmp(ji,jj,1) = - r1_rho0 * ( emp(ji,jj) & ! volume + & - rnf(ji,jj) & + & - fwfisf_cav(ji,jj) & + & - fwfisf_par(ji,jj) ) * surf(ji,jj) + ztmp(ji,jj,2) = sbc_tsc(ji,jj,jp_tem) * surf(ji,jj) ! heat + ztmp(ji,jj,3) = sbc_tsc(ji,jj,jp_sal) * surf(ji,jj) ! salt + END_2D + IF( ln_rnf ) THEN + DO_2D( 0, 0, 0, 0 ) + ztmp(ji,jj,4) = rnf_tsc(ji,jj,jp_tem) * surf(ji,jj) ! runoff temp + END_2D + END IF + IF( ln_rnf_sal ) THEN + DO_2D( 0, 0, 0, 0 ) + ztmp(ji,jj,5) = rnf_tsc(ji,jj,jp_sal) * surf(ji,jj) ! runoff salt + END_2D + END IF + IF( ln_isf ) THEN + DO_2D( 0, 0, 0, 0 ) + ztmp(ji,jj,6) = ( risf_cav_tsc(ji,jj,jp_tem) & + & + risf_par_tsc(ji,jj,jp_tem) ) * surf(ji,jj) ! isf temp + END_2D + END IF + IF( ln_traqsr ) THEN + DO_2D( 0, 0, 0, 0 ) + ztmp(ji,jj,7) = r1_rho0_rcp * qsr(ji,jj) * surf(ji,jj) ! penetrative solar radiation + END_2D + END IF + IF( ln_trabbc ) THEN + DO_2D( 0, 0, 0, 0 ) + ztmp(ji,jj,8) = qgh_trd0(ji,jj) * surf(ji,jj) ! geothermal heat + END_2D + END IF ! IF( ln_linssh ) THEN ! Advection flux through fixed surface (z=0) IF( ln_isfcav ) THEN @@ -116,8 +148,10 @@ CONTAINS END DO END DO ELSE - ztmp(:,:,9 ) = - surf(:,:) * ww(:,:,1) * ts(:,:,1,jp_tem,Kbb) - ztmp(:,:,10) = - surf(:,:) * ww(:,:,1) * ts(:,:,1,jp_sal,Kbb) + DO_2D( 0, 0, 0, 0 ) + ztmp(ji,jj,9 ) = - surf(ji,jj) * ww(ji,jj,1) * ts(ji,jj,1,jp_tem,Kbb) + ztmp(ji,jj,10) = - surf(ji,jj) * ww(ji,jj,1) * ts(ji,jj,1,jp_sal,Kbb) + END_2D END IF ENDIF @@ -152,7 +186,9 @@ CONTAINS ! glob_sum is needed because you keep only the interior domain to compute the sum (iscpl) ! ! ! volume variation (calculated with ssh) - ztmp(:,:,11) = surf(:,:)*ssh(:,:,Kmm) - surf_ini(:,:)*ssh_ini(:,:) + DO_2D( 0, 0, 0, 0 ) + ztmp(ji,jj,11) = surf(ji,jj)*ssh(ji,jj,Kmm) - surf_ini(ji,jj)*ssh_ini(ji,jj) + END_2D ! ! heat & salt content variation (associated with ssh) IF( ln_linssh ) THEN ! linear free surface case @@ -164,8 +200,10 @@ CONTAINS END DO END DO ELSE ! no under ice-shelf seas - ztmp(:,:,12) = surf(:,:) * ( ts(:,:,1,jp_tem,Kmm) * ssh(:,:,Kmm) - ssh_hc_loc_ini(:,:) ) - ztmp(:,:,13) = surf(:,:) * ( ts(:,:,1,jp_sal,Kmm) * ssh(:,:,Kmm) - ssh_sc_loc_ini(:,:) ) + DO_2D( 0, 0, 0, 0 ) + ztmp(ji,jj,12) = surf(ji,jj) * ( ts(ji,jj,1,jp_tem,Kmm) * ssh(ji,jj,Kmm) - ssh_hc_loc_ini(ji,jj) ) + ztmp(ji,jj,13) = surf(ji,jj) * ( ts(ji,jj,1,jp_sal,Kmm) * ssh(ji,jj,Kmm) - ssh_sc_loc_ini(ji,jj) ) + END_2D END IF ENDIF @@ -185,19 +223,27 @@ CONTAINS ! glob_sum is needed because you keep only the interior domain to compute the sum (iscpl) ! DO jk = 1, jpkm1 ! volume - ztmpk(:,:,jk,1) = surf (:,:) * e3t(:,:,jk,Kmm)*tmask(:,:,jk) & - & - surf_ini(:,:) * e3t_ini(:,:,jk )*tmask_ini(:,:,jk) + DO_2D( 0, 0, 0, 0 ) + ztmpk(ji,jj,jk,1) = surf (ji,jj) * e3t(ji,jj,jk,Kmm)*tmask(ji,jj,jk) & + & - surf_ini(ji,jj) * e3t_ini(ji,jj,jk )*tmask_ini(ji,jj,jk) + END_2D END DO DO jk = 1, jpkm1 ! heat - ztmpk(:,:,jk,2) = ( surf (:,:) * e3t(:,:,jk,Kmm)*ts(:,:,jk,jp_tem,Kmm) & - & - surf_ini(:,:) * hc_loc_ini(:,:,jk) ) + DO_2D( 0, 0, 0, 0 ) + ztmpk(ji,jj,jk,2) = ( surf (ji,jj) * e3t(ji,jj,jk,Kmm)*ts(ji,jj,jk,jp_tem,Kmm) & + & - surf_ini(ji,jj) * hc_loc_ini(ji,jj,jk) ) + END_2D END DO DO jk = 1, jpkm1 ! salt - ztmpk(:,:,jk,3) = ( surf (:,:) * e3t(:,:,jk,Kmm)*ts(:,:,jk,jp_sal,Kmm) & - & - surf_ini(:,:) * sc_loc_ini(:,:,jk) ) + DO_2D( 0, 0, 0, 0 ) + ztmpk(ji,jj,jk,3) = ( surf (ji,jj) * e3t(ji,jj,jk,Kmm)*ts(ji,jj,jk,jp_sal,Kmm) & + & - surf_ini(ji,jj) * sc_loc_ini(ji,jj,jk) ) + END_2D END DO DO jk = 1, jpkm1 ! total ocean volume - ztmpk(:,:,jk,4) = surf(:,:) * e3t(:,:,jk,Kmm) * tmask(:,:,jk) + DO_2D( 0, 0, 0, 0 ) + ztmpk(ji,jj,jk,4) = surf(ji,jj) * e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk) + END_2D END DO ! global sum @@ -315,14 +361,18 @@ CONTAINS IF(lwp) WRITE(numout,*) IF(lwp) WRITE(numout,*) ' dia_hsb_rst : initialise hsb at initial state ' IF(lwp) WRITE(numout,*) - surf_ini(:,:) = e1e2t(:,:) * tmask_i(:,:) ! initial ocean surface - ssh_ini(:,:) = ssh(:,:,Kmm) ! initial ssh + DO_2D( 0, 0, 0, 0 ) + surf_ini(ji,jj) = e1e2t(ji,jj) * tmask_i(ji,jj) ! initial ocean surface + ssh_ini(ji,jj) = ssh(ji,jj,Kmm) ! initial ssh + END_2D DO jk = 1, jpk ! if ice sheet/oceqn coupling, need to mask ini variables here (mask could change at the next NEMO instance). - e3t_ini (:,:,jk) = e3t(:,:,jk,Kmm) * tmask(:,:,jk) ! initial vertical scale factors - tmask_ini (:,:,jk) = tmask(:,:,jk) ! initial mask - hc_loc_ini(:,:,jk) = ts(:,:,jk,jp_tem,Kmm) * e3t(:,:,jk,Kmm) * tmask(:,:,jk) ! initial heat content - sc_loc_ini(:,:,jk) = ts(:,:,jk,jp_sal,Kmm) * e3t(:,:,jk,Kmm) * tmask(:,:,jk) ! initial salt content + DO_2D( 0, 0, 0, 0 ) + e3t_ini (ji,jj,jk) = e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk) ! initial vertical scale factors + tmask_ini (ji,jj,jk) = tmask(ji,jj,jk) ! initial mask + hc_loc_ini(ji,jj,jk) = ts(ji,jj,jk,jp_tem,Kmm) * e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk) ! initial heat content + sc_loc_ini(ji,jj,jk) = ts(ji,jj,jk,jp_sal,Kmm) * e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk) ! initial salt content + END_2D END DO frc_v = 0._wp ! volume trend due to forcing frc_t = 0._wp ! heat content - - - - @@ -334,13 +384,15 @@ CONTAINS ssh_hc_loc_ini(ji,jj) = ts(ji,jj,mikt(ji,jj),jp_tem,Kmm) * ssh(ji,jj,Kmm) ! initial heat content in ssh ssh_sc_loc_ini(ji,jj) = ts(ji,jj,mikt(ji,jj),jp_sal,Kmm) * ssh(ji,jj,Kmm) ! initial salt content in ssh END DO - END DO - ELSE - ssh_hc_loc_ini(:,:) = ts(:,:,1,jp_tem,Kmm) * ssh(:,:,Kmm) ! initial heat content in ssh - ssh_sc_loc_ini(:,:) = ts(:,:,1,jp_sal,Kmm) * ssh(:,:,Kmm) ! initial salt content in ssh + END DO + ELSE + DO_2D( 0, 0, 0, 0 ) + ssh_hc_loc_ini(ji,jj) = ts(ji,jj,1,jp_tem,Kmm) * ssh(ji,jj,Kmm) ! initial heat content in ssh + ssh_sc_loc_ini(ji,jj) = ts(ji,jj,1,jp_sal,Kmm) * ssh(ji,jj,Kmm) ! initial salt content in ssh + END_2D END IF - frc_wn_t = 0._wp ! initial heat content misfit due to free surface - frc_wn_s = 0._wp ! initial salt content misfit due to free surface + frc_wn_t = 0._wp ! initial heat content misfit due to free surface + frc_wn_s = 0._wp ! initial salt content misfit due to free surface ENDIF ENDIF ! @@ -388,6 +440,7 @@ CONTAINS INTEGER, INTENT(in) :: Kmm ! time level index ! INTEGER :: ierror, ios ! local integer + INTEGER :: ji, jj ! loop index !! NAMELIST/namhsb/ ln_diahsb !!---------------------------------------------------------------------- @@ -427,7 +480,10 @@ CONTAINS ! ----------------------------------------------- ! ! 2 - Time independant variables and file opening ! ! ----------------------------------------------- ! - surf(:,:) = e1e2t(:,:) * tmask_i(:,:) ! masked surface grid cell area + + DO_2D( 0, 0, 0, 0 ) + surf(ji,jj) = e1e2t(ji,jj) * tmask_i(ji,jj) ! masked surface grid cell area + END_2D surf_tot = glob_sum( 'diahsb', surf(:,:) ) ! total ocean surface area IF( ln_bdy ) CALL ctl_warn( 'dia_hsb_init: heat/salt budget does not consider open boundary fluxes' ) diff --git a/src/OCE/DIA/diahth.F90 b/src/OCE/DIA/diahth.F90 index 0b8650acf..dcec1291f 100644 --- a/src/OCE/DIA/diahth.F90 +++ b/src/OCE/DIA/diahth.F90 @@ -86,22 +86,22 @@ CONTAINS INTEGER, INTENT( in ) :: kt ! ocean time-step index INTEGER, INTENT( in ) :: Kmm ! ocean time level index !! - INTEGER :: ji, jj, jk ! dummy loop arguments - REAL(wp) :: zrho3 = 0.03_wp ! density criterion for mixed layer depth - REAL(wp) :: zrho1 = 0.01_wp ! density criterion for mixed layer depth - REAL(wp) :: ztem2 = 0.2_wp ! temperature criterion for mixed layer depth - REAL(wp) :: zztmp, zzdep ! temporary scalars inside do loop - REAL(wp) :: zu, zv, zw, zut, zvt ! temporary workspace - REAL(wp), DIMENSION(jpi,jpj) :: zabs2 ! MLD: abs( tn - tn(10m) ) = ztem2 - REAL(wp), DIMENSION(jpi,jpj) :: ztm2 ! Top of thermocline: tn = tn(10m) - ztem2 - REAL(wp), DIMENSION(jpi,jpj) :: zrho10_3 ! MLD: rho = rho10m + zrho3 - REAL(wp), DIMENSION(jpi,jpj) :: zpycn ! pycnocline: rho = rho10m + (dr/dT)(T,S,10m)*(-0.2 degC) - REAL(wp), DIMENSION(jpi,jpj) :: ztinv ! max of temperature inversion - REAL(wp), DIMENSION(jpi,jpj) :: zdepinv ! depth of temperature inversion - REAL(wp), DIMENSION(jpi,jpj) :: zrho0_3 ! MLD rho = rho(surf) = 0.03 - REAL(wp), DIMENSION(jpi,jpj) :: zrho0_1 ! MLD rho = rho(surf) = 0.01 - REAL(wp), DIMENSION(jpi,jpj) :: zmaxdzT ! max of dT/dz - REAL(wp), DIMENSION(jpi,jpj) :: zdelr ! delta rho equivalent to deltaT = 0.2 + INTEGER :: ji, jj, jk ! dummy loop arguments + REAL(wp) :: zrho3 = 0.03_wp ! density criterion for mixed layer depth + REAL(wp) :: zrho1 = 0.01_wp ! density criterion for mixed layer depth + REAL(wp) :: ztem2 = 0.2_wp ! temperature criterion for mixed layer depth + REAL(wp) :: zztmp, zzdep ! temporary scalars inside do loop + REAL(wp) :: zu, zv, zw, zut, zvt ! temporary workspace + REAL(wp), DIMENSION(A2D(0)) :: zabs2 ! MLD: abs( tn - tn(10m) ) = ztem2 + REAL(wp), DIMENSION(A2D(0)) :: ztm2 ! Top of thermocline: tn = tn(10m) - ztem2 + REAL(wp), DIMENSION(A2D(0)) :: zrho10_3 ! MLD: rho = rho10m + zrho3 + REAL(wp), DIMENSION(A2D(0)) :: zpycn ! pycnocline: rho = rho10m + (dr/dT)(T,S,10m)*(-0.2 degC) + REAL(wp), DIMENSION(A2D(0)) :: ztinv ! max of temperature inversion + REAL(wp), DIMENSION(A2D(0)) :: zdepinv ! depth of temperature inversion + REAL(wp), DIMENSION(A2D(0)) :: zrho0_3 ! MLD rho = rho(surf) = 0.03 + REAL(wp), DIMENSION(A2D(0)) :: zrho0_1 ! MLD rho = rho(surf) = 0.01 + REAL(wp), DIMENSION(A2D(0)) :: zmaxdzT ! max of dT/dz + REAL(wp), DIMENSION(A2D(0)) :: zdelr ! delta rho equivalent to deltaT = 0.2 !!---------------------------------------------------------------------- IF( ln_timing ) CALL timing_start('dia_hth') @@ -131,7 +131,7 @@ CONTAINS IF( iom_use( 'mlddzt' ) ) zmaxdzT(:,:) = 0._wp IF( iom_use( 'mlddzt' ) .OR. iom_use( 'mld_dt02' ) .OR. iom_use( 'topthdep' ) & & .OR. iom_use( 'mldr10_3' ) .OR. iom_use( 'pycndep' ) ) THEN - DO_2D( 1, 1, 1, 1 ) + DO_2D( 0, 0, 0, 0 ) zztmp = gdepw(ji,jj,mbkt(ji,jj)+1,Kmm) hth (ji,jj) = zztmp zabs2 (ji,jj) = zztmp @@ -142,7 +142,7 @@ CONTAINS ENDIF IF( iom_use( 'mldr0_3' ) .OR. iom_use( 'mldr0_1' ) ) THEN IF( nla10 > 1 ) THEN - DO_2D( 1, 1, 1, 1 ) + DO_2D( 0, 0, 0, 0 ) zztmp = gdepw(ji,jj,mbkt(ji,jj)+1,Kmm) zrho0_3(ji,jj) = zztmp zrho0_1(ji,jj) = zztmp @@ -157,7 +157,7 @@ CONTAINS ! MLD: rho = rho(1) + zrho3 ! ! MLD: rho = rho(1) + zrho1 ! ! ------------------------------------------------------------- ! - DO_3DS( 1, 1, 1, 1, jpkm1, 2, -1 ) ! loop from bottom to 2 + DO_3DS( 0, 0, 0, 0, jpkm1, 2, -1 ) ! loop from bottom to 2 ! zzdep = gdepw(ji,jj,jk,Kmm) zztmp = ( ts(ji,jj,jk-1,jp_tem,Kmm) - ts(ji,jj,jk,jp_tem,Kmm) ) & @@ -189,7 +189,7 @@ CONTAINS ! ! Preliminary computation ! computation of zdelr = (dr/dT)(T,S,10m)*(-0.2 degC) - DO_2D( 1, 1, 1, 1 ) + DO_2D( 0, 0, 0, 0 ) IF( tmask(ji,jj,nla10) == 1. ) THEN zu = 1779.50 + 11.250 * ts(ji,jj,nla10,jp_tem,Kmm) - 3.80 * ts(ji,jj,nla10,jp_sal,Kmm) & & - 0.0745 * ts(ji,jj,nla10,jp_tem,Kmm) * ts(ji,jj,nla10,jp_tem,Kmm) & @@ -213,7 +213,7 @@ CONTAINS ! temperature inversion: max( 0, max of tn - tn(10m) ) ! ! depth of temperature inversion ! ! ------------------------------------------------------------- ! - DO_3DS( 1, 1, 1, 1, jpkm1, nlb10, -1 ) ! loop from bottom to nlb10 + DO_3DS( 0, 0, 0, 0, jpkm1, nlb10, -1 ) ! loop from bottom to nlb10 ! zzdep = gdepw(ji,jj,jk,Kmm) * tmask(ji,jj,1) ! @@ -305,13 +305,16 @@ CONTAINS ! INTEGER :: ji, jj, jk, iid REAL(wp) :: zztmp, zzdep - INTEGER, DIMENSION(jpi,jpj) :: iktem + INTEGER, DIMENSION(A2D(0)) :: iktem ! --------------------------------------- ! ! search deepest level above ptem ! ! --------------------------------------- ! - iktem(:,:) = 1 - DO_3D( 1, 1, 1, 1, 1, jpkm1 ) ! beware temperature is not always decreasing with depth => loop from top to bottom + DO_2D( 0, 0, 0, 0 ) + iktem(ji,jj) = 1 + END_2D + + DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! beware temperature is not always decreasing with depth => loop from top to bottom zztmp = ts(ji,jj,jk,jp_tem,Kmm) IF( zztmp >= ptem ) iktem(ji,jj) = jk END_3D @@ -319,7 +322,7 @@ CONTAINS ! ------------------------------- ! ! Depth of ptem isotherm ! ! ------------------------------- ! - DO_2D( 1, 1, 1, 1 ) + DO_2D( 0, 0, 0, 0 ) ! zzdep = gdepw(ji,jj,mbkt(ji,jj)+1,Kmm) ! depth of the ocean bottom ! @@ -346,18 +349,29 @@ CONTAINS REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: phtc ! INTEGER :: ji, jj, jk, ik - REAL(wp), DIMENSION(jpi,jpj) :: zthick - INTEGER , DIMENSION(jpi,jpj) :: ilevel + REAL(wp), DIMENSION(A2D(0)) :: zthick + INTEGER , DIMENSION(A2D(0)) :: ilevel ! surface boundary condition - IF( .NOT. ln_linssh ) THEN ; zthick(:,:) = 0._wp ; phtc(:,:) = 0._wp - ELSE ; zthick(:,:) = ssh(:,:,Kmm) ; phtc(:,:) = pt(:,:,1) * ssh(:,:,Kmm) * tmask(:,:,1) + IF( .NOT. ln_linssh ) THEN + DO_2D( 0, 0, 0, 0 ) + zthick(ji,jj) = 0._wp + phtc (ji,jj) = 0._wp + END_2D + ELSE + DO_2D( 0, 0, 0, 0 ) + zthick(ji,jj) = ssh(ji,jj,Kmm) + phtc (ji,jj) = pt(ji,jj,1) * ssh(ji,jj,Kmm) * tmask(ji,jj,1) + END_2D ENDIF ! - ilevel(:,:) = 1 - DO_3D( 1, 1, 1, 1, 1, jpkm1 ) + DO_2D( 0, 0, 0, 0 ) + ilevel(ji,jj) = 1 + END_2D + ! + DO_3D( 0, 0, 0, 0, 1, jpkm1 ) IF( ( gdepw(ji,jj,jk+1,Kmm) < pdep ) .AND. ( tmask(ji,jj,jk) == 1 ) ) THEN ilevel(ji,jj) = jk+1 zthick(ji,jj) = zthick(ji,jj) + e3t(ji,jj,jk,Kmm) @@ -365,7 +379,7 @@ CONTAINS ENDIF END_3D ! - DO_2D( 1, 1, 1, 1 ) + DO_2D( 0, 0, 0, 0 ) ik = ilevel(ji,jj) IF( tmask(ji,jj,ik) == 1 ) THEN zthick(ji,jj) = MIN ( gdepw(ji,jj,ik+1,Kmm), pdep ) - zthick(ji,jj) ! remaining thickness to reach dephw pdep diff --git a/src/OCE/DIA/diamlr.F90 b/src/OCE/DIA/diamlr.F90 index edd8e19f9..dbcf5bbea 100644 --- a/src/OCE/DIA/diamlr.F90 +++ b/src/OCE/DIA/diamlr.F90 @@ -6,7 +6,7 @@ MODULE diamlr !! History : 4.0 ! 2019 (S. Mueller) Original code !!---------------------------------------------------------------------- - USE par_oce , ONLY : wp, jpi, jpj + USE par_oce , ONLY : wp, jpi, jpj, ntsi, ntei, ntsj, ntej USE phycst , ONLY : rpi USE dom_oce , ONLY : adatrj USE tide_mod @@ -407,8 +407,9 @@ CONTAINS !! ** Purpose : update time used in multiple-linear-regression analysis !! !!---------------------------------------------------------------------- - REAL(wp), DIMENSION(jpi,jpj) :: zadatrj2d + REAL(wp), DIMENSION(A2D(0)) :: zadatrj2d !!---------------------------------------------------------------------- + INTEGER :: ji, jj IF( ln_timing ) CALL timing_start('dia_mlr') @@ -417,7 +418,9 @@ CONTAINS ! ! A 2-dimensional field of constant value is sent, and subsequently used directly ! or transformed to a scalar or a constant 3-dimensional field as required. - zadatrj2d(:,:) = adatrj*86400.0_wp + DO_2D( 0, 0, 0, 0 ) + zadatrj2d(ji,jj) = adatrj*86400.0_wp + END_2D IF ( iom_use('diamlr_time') ) CALL iom_put('diamlr_time', zadatrj2d) ! IF( ln_timing ) CALL timing_stop('dia_mlr') diff --git a/src/OCE/DIA/diaptr.F90 b/src/OCE/DIA/diaptr.F90 index 8962bd371..6c4f48a3d 100644 --- a/src/OCE/DIA/diaptr.F90 +++ b/src/OCE/DIA/diaptr.F90 @@ -78,7 +78,7 @@ CONTAINS !!---------------------------------------------------------------------- INTEGER , INTENT(in) :: kt ! ocean time-step index INTEGER , INTENT(in) :: Kmm ! time level index - REAL(wp), DIMENSION(A2D(nn_hls),jpk) , INTENT(in), OPTIONAL :: pvtr ! j-effective transport + REAL(wp), DIMENSION(A2D(0),jpk) , INTENT(in), OPTIONAL :: pvtr ! j-effective transport !!---------------------------------------------------------------------- ! IF( ln_timing ) CALL timing_start('dia_ptr') @@ -110,13 +110,13 @@ CONTAINS !!---------------------------------------------------------------------- !! ** Purpose : Calculate diagnostics and send to XIOS !!---------------------------------------------------------------------- - INTEGER , INTENT(in) :: kt ! ocean time-step index - INTEGER , INTENT(in) :: Kmm ! time level index - REAL(wp), DIMENSION(A2D(nn_hls),jpk) , INTENT(in), OPTIONAL :: pvtr ! j-effective transport + INTEGER , INTENT(in) :: kt ! ocean time-step index + INTEGER , INTENT(in) :: Kmm ! time level index + REAL(wp), DIMENSION(:,:,:), INTENT(in), OPTIONAL :: pvtr ! j-effective transport (used only by PRESENT) ! INTEGER :: ji, jj, jk, jn ! dummy loop indices - REAL(wp), DIMENSION(jpi,jpj) :: z2d ! 2D workspace - REAL(wp), DIMENSION(jpj) :: zvsum, ztsum, zssum ! 1D workspace + REAL(wp), DIMENSION(Nis0:Nie0,Njs0:Nje0) :: z2d ! 2D workspace + REAL(wp), DIMENSION( Njs0:Nje0) :: zvsum, ztsum, zssum ! 1D workspace ! !overturning calculation REAL(wp), DIMENSION(:,:,: ), ALLOCATABLE :: sjk, r1_sjk, v_msf ! i-mean i-k-surface and its inverse @@ -126,19 +126,19 @@ CONTAINS REAL(wp), DIMENSION(:,:,: ), ALLOCATABLE :: z3dtr !!---------------------------------------------------------------------- ! - ALLOCATE( z3dtr(jpi,jpj,nbasin) ) + ALLOCATE( z3dtr(Nis0:Nie0,Njs0:Nje0,nbasin) ) IF( PRESENT( pvtr ) ) THEN IF( iom_use( 'zomsf' ) ) THEN ! effective MSF - ALLOCATE( z4d1(jpi,jpj,jpk,nbasin) ) + ALLOCATE( z4d1(Nis0:Nie0,Njs0:Nje0,jpk,nbasin) ) ! DO jn = 1, nbasin ! by sub-basins - z4d1(1,:,:,jn) = pvtr_int(:,:,jp_vtr,jn) ! zonal cumulative effective transport excluding closed seas + z4d1(Nis0,:,:,jn) = pvtr_int(:,:,jp_vtr,jn) ! zonal cumulative effective transport excluding closed seas DO jk = jpkm1, 1, -1 - z4d1(1,:,jk,jn) = z4d1(1,:,jk+1,jn) - z4d1(1,:,jk,jn) ! effective j-Stream-Function (MSF) + z4d1(Nis0,:,jk,jn) = z4d1(Nis0,:,jk+1,jn) - z4d1(Nis0,:,jk,jn) ! effective j-Stream-Function (MSF) END DO - DO ji = 2, jpi - z4d1(ji,:,:,jn) = z4d1(1,:,:,jn) + DO ji = Nis0+1, Nie0 + z4d1(ji,:,:,jn) = z4d1(Nis0,:,:,jn) ENDDO END DO CALL iom_put( 'zomsf', z4d1 * rc_sv ) @@ -146,8 +146,8 @@ CONTAINS DEALLOCATE( z4d1 ) ENDIF IF( iom_use( 'sopstove' ) .OR. iom_use( 'sophtove' ) ) THEN - ALLOCATE( sjk(jpj,jpk,nbasin), r1_sjk(jpj,jpk,nbasin), v_msf(jpj,jpk,nbasin), & - & zt_jk(jpj,jpk,nbasin), zs_jk(jpj,jpk,nbasin) ) + ALLOCATE( sjk( Njs0:Nje0,jpk,nbasin), r1_sjk(Njs0:Nje0,jpk,nbasin), v_msf(Njs0:Nje0,jpk,nbasin), & + & zt_jk(Njs0:Nje0,jpk,nbasin), zs_jk( Njs0:Nje0,jpk,nbasin) ) ! DO jn = 1, nbasin sjk(:,:,jn) = pvtr_int(:,:,jp_msk,jn) @@ -162,16 +162,16 @@ CONTAINS ! ENDDO DO jn = 1, nbasin - z3dtr(1,:,jn) = hstr_ove(:,jp_tem,jn) * rc_pwatt ! (conversion in PW) - DO ji = 2, jpi - z3dtr(ji,:,jn) = z3dtr(1,:,jn) + z3dtr(Nis0,:,jn) = hstr_ove(:,jp_tem,jn) * rc_pwatt ! (conversion in PW) + DO ji = Nis0+1, Nie0 + z3dtr(ji,:,jn) = z3dtr(Nis0,:,jn) ENDDO ENDDO CALL iom_put( 'sophtove', z3dtr ) DO jn = 1, nbasin - z3dtr(1,:,jn) = hstr_ove(:,jp_sal,jn) * rc_ggram ! (conversion in Gg) - DO ji = 2, jpi - z3dtr(ji,:,jn) = z3dtr(1,:,jn) + z3dtr(Nis0,:,jn) = hstr_ove(:,jp_sal,jn) * rc_ggram ! (conversion in Gg) + DO ji = Nis0+1, Nie0 + z3dtr(ji,:,jn) = z3dtr(Nis0,:,jn) ENDDO ENDDO CALL iom_put( 'sopstove', z3dtr ) @@ -181,7 +181,7 @@ CONTAINS IF( iom_use( 'sopstbtr' ) .OR. iom_use( 'sophtbtr' ) ) THEN ! Calculate barotropic heat and salt transport here - ALLOCATE( sjk(jpj,1,nbasin), r1_sjk(jpj,1,nbasin) ) + ALLOCATE( sjk(A1Dj(0),1,nbasin), r1_sjk(A1Dj(0),1,nbasin) ) ! DO jn = 1, nbasin sjk(:,1,jn) = SUM( pvtr_int(:,:,jp_msk,jn), 2 ) @@ -196,16 +196,16 @@ CONTAINS ! ENDDO DO jn = 1, nbasin - z3dtr(1,:,jn) = hstr_btr(:,jp_tem,jn) * rc_pwatt ! (conversion in PW) - DO ji = 2, jpi - z3dtr(ji,:,jn) = z3dtr(1,:,jn) + z3dtr(Nis0,:,jn) = hstr_btr(:,jp_tem,jn) * rc_pwatt ! (conversion in PW) + DO ji = Nis0+1, Nie0 + z3dtr(ji,:,jn) = z3dtr(Nis0,:,jn) ENDDO ENDDO CALL iom_put( 'sophtbtr', z3dtr ) DO jn = 1, nbasin - z3dtr(1,:,jn) = hstr_btr(:,jp_sal,jn) * rc_ggram ! (conversion in Gg) - DO ji = 2, jpi - z3dtr(ji,:,jn) = z3dtr(1,:,jn) + z3dtr(Nis0,:,jn) = hstr_btr(:,jp_sal,jn) * rc_ggram ! (conversion in Gg) + DO ji = Nis0+1, Nie0 + z3dtr(ji,:,jn) = z3dtr(Nis0,:,jn) ENDDO ENDDO CALL iom_put( 'sopstbtr', z3dtr ) @@ -218,28 +218,28 @@ CONTAINS pvtr_int(:,:,:,:) = 0._wp ELSE IF( iom_use( 'zotem' ) .OR. iom_use( 'zosal' ) .OR. iom_use( 'zosrf' ) ) THEN ! i-mean i-k-surface - ALLOCATE( z4d1(jpi,jpj,jpk,nbasin), z4d2(jpi,jpj,jpk,nbasin) ) + ALLOCATE( z4d1(Nis0:Nie0,Njs0:Nje0,jpk,nbasin), z4d2(Nis0:Nie0,Njs0:Nje0,jpk,nbasin) ) ! DO jn = 1, nbasin - z4d1(1,:,:,jn) = pzon_int(:,:,jp_msk,jn) - DO ji = 2, jpi - z4d1(ji,:,:,jn) = z4d1(1,:,:,jn) + z4d1(Nis0,:,:,jn) = pzon_int(:,:,jp_msk,jn) + DO ji = Nis0+1, Nie0 + z4d1(ji,:,:,jn) = z4d1(Nis0,:,:,jn) ENDDO ENDDO CALL iom_put( 'zosrf', z4d1 ) ! DO jn = 1, nbasin - z4d2(1,:,:,jn) = pzon_int(:,:,jp_tem,jn) / MAX( z4d1(1,:,:,jn), 10.e-15 ) - DO ji = 2, jpi - z4d2(ji,:,:,jn) = z4d2(1,:,:,jn) + z4d2(Nis0,:,:,jn) = pzon_int(:,:,jp_tem,jn) / MAX( z4d1(Nis0,:,:,jn), 10.e-15 ) + DO ji = Nis0+1, Nie0 + z4d2(ji,:,:,jn) = z4d2(Nis0,:,:,jn) ENDDO ENDDO CALL iom_put( 'zotem', z4d2 ) ! DO jn = 1, nbasin - z4d2(1,:,:,jn) = pzon_int(:,:,jp_sal,jn) / MAX( z4d1(1,:,:,jn), 10.e-15 ) - DO ji = 2, jpi - z4d2(ji,:,:,jn) = z4d2(1,:,:,jn) + z4d2(Nis0,:,:,jn) = pzon_int(:,:,jp_sal,jn) / MAX( z4d1(Nis0,:,:,jn), 10.e-15 ) + DO ji = Nis0+1, Nie0 + z4d2(ji,:,:,jn) = z4d2(Nis0,:,:,jn) ENDDO ENDDO CALL iom_put( 'zosal', z4d2 ) @@ -251,16 +251,16 @@ CONTAINS IF( iom_use( 'sophtadv' ) .OR. iom_use( 'sopstadv' ) ) THEN ! DO jn = 1, nbasin - z3dtr(1,:,jn) = hstr_adv(:,jp_tem,jn) * rc_pwatt ! (conversion in PW) - DO ji = 2, jpi - z3dtr(ji,:,jn) = z3dtr(1,:,jn) + z3dtr(Nis0,:,jn) = hstr_adv(:,jp_tem,jn) * rc_pwatt ! (conversion in PW) + DO ji = Nis0+1, Nie0 + z3dtr(ji,:,jn) = z3dtr(Nis0,:,jn) ENDDO ENDDO CALL iom_put( 'sophtadv', z3dtr ) DO jn = 1, nbasin - z3dtr(1,:,jn) = hstr_adv(:,jp_sal,jn) * rc_ggram ! (conversion in Gg) - DO ji = 2, jpi - z3dtr(ji,:,jn) = z3dtr(1,:,jn) + z3dtr(Nis0,:,jn) = hstr_adv(:,jp_sal,jn) * rc_ggram ! (conversion in Gg) + DO ji = Nis0+1, Nie0 + z3dtr(ji,:,jn) = z3dtr(Nis0,:,jn) ENDDO ENDDO CALL iom_put( 'sopstadv', z3dtr ) @@ -269,16 +269,16 @@ CONTAINS IF( iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) ) THEN ! DO jn = 1, nbasin - z3dtr(1,:,jn) = hstr_ldf(:,jp_tem,jn) * rc_pwatt ! (conversion in PW) - DO ji = 2, jpi - z3dtr(ji,:,jn) = z3dtr(1,:,jn) + z3dtr(Nis0,:,jn) = hstr_ldf(:,jp_tem,jn) * rc_pwatt ! (conversion in PW) + DO ji = Nis0+1, Nie0 + z3dtr(ji,:,jn) = z3dtr(Nis0,:,jn) ENDDO ENDDO CALL iom_put( 'sophtldf', z3dtr ) DO jn = 1, nbasin - z3dtr(1,:,jn) = hstr_ldf(:,jp_sal,jn) * rc_ggram ! (conversion in Gg) - DO ji = 2, jpi - z3dtr(ji,:,jn) = z3dtr(1,:,jn) + z3dtr(Nis0,:,jn) = hstr_ldf(:,jp_sal,jn) * rc_ggram ! (conversion in Gg) + DO ji = Nis0+1, Nie0 + z3dtr(ji,:,jn) = z3dtr(Nis0,:,jn) ENDDO ENDDO CALL iom_put( 'sopstldf', z3dtr ) @@ -287,16 +287,16 @@ CONTAINS IF( iom_use( 'sophteiv' ) .OR. iom_use( 'sopsteiv' ) ) THEN ! DO jn = 1, nbasin - z3dtr(1,:,jn) = hstr_eiv(:,jp_tem,jn) * rc_pwatt ! (conversion in PW) - DO ji = 2, jpi - z3dtr(ji,:,jn) = z3dtr(1,:,jn) + z3dtr(Nis0,:,jn) = hstr_eiv(:,jp_tem,jn) * rc_pwatt ! (conversion in PW) + DO ji = Nis0+1, Nie0 + z3dtr(ji,:,jn) = z3dtr(Nis0,:,jn) ENDDO ENDDO CALL iom_put( 'sophteiv', z3dtr ) DO jn = 1, nbasin - z3dtr(1,:,jn) = hstr_eiv(:,jp_sal,jn) * rc_ggram ! (conversion in Gg) - DO ji = 2, jpi - z3dtr(ji,:,jn) = z3dtr(1,:,jn) + z3dtr(Nis0,:,jn) = hstr_eiv(:,jp_sal,jn) * rc_ggram ! (conversion in Gg) + DO ji = Nis0+1, Nie0 + z3dtr(ji,:,jn) = z3dtr(Nis0,:,jn) ENDDO ENDDO CALL iom_put( 'sopsteiv', z3dtr ) @@ -304,16 +304,16 @@ CONTAINS ! IF( iom_use( 'sopstvtr' ) .OR. iom_use( 'sophtvtr' ) ) THEN DO jn = 1, nbasin - z3dtr(1,:,jn) = hstr_vtr(:,jp_tem,jn) * rc_pwatt ! (conversion in PW) - DO ji = 2, jpi - z3dtr(ji,:,jn) = z3dtr(1,:,jn) + z3dtr(Nis0,:,jn) = hstr_vtr(:,jp_tem,jn) * rc_pwatt ! (conversion in PW) + DO ji = Nis0+1, Nie0 + z3dtr(ji,:,jn) = z3dtr(Nis0,:,jn) ENDDO ENDDO CALL iom_put( 'sophtvtr', z3dtr ) DO jn = 1, nbasin - z3dtr(1,:,jn) = hstr_vtr(:,jp_sal,jn) * rc_ggram ! (conversion in Gg) - DO ji = 2, jpi - z3dtr(ji,:,jn) = z3dtr(1,:,jn) + z3dtr(Nis0,:,jn) = hstr_vtr(:,jp_sal,jn) * rc_ggram ! (conversion in Gg) + DO ji = Nis0+1, Nie0 + z3dtr(ji,:,jn) = z3dtr(Nis0,:,jn) ENDDO ENDDO CALL iom_put( 'sopstvtr', z3dtr ) @@ -349,8 +349,8 @@ CONTAINS !! ** Action : pvtr_int - terms for volume streamfunction, heat/salt transport barotropic/overturning terms !! pzon_int - terms for i mean temperature/salinity !!---------------------------------------------------------------------- - INTEGER , INTENT(in) :: Kmm ! time level index - REAL(wp), DIMENSION(A2D(nn_hls),jpk), INTENT(in), OPTIONAL :: pvtr ! j-effective transport + INTEGER , INTENT(in) :: Kmm ! time level index + REAL(wp), DIMENSION(A2D(0),jpk), INTENT(in), OPTIONAL :: pvtr ! j-effective transport REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: zmask ! 3D workspace REAL(wp), DIMENSION(:,:,:,:), ALLOCATABLE :: zts ! 4D workspace REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: sjk, v_msf ! Zonal sum: i-k surface area, j-effective transport @@ -362,7 +362,7 @@ CONTAINS IF( PRESENT( pvtr ) ) THEN ! i sum of effective j transport excluding closed seas IF( iom_use( 'zomsf' ) .OR. iom_use( 'sopstove' ) .OR. iom_use( 'sophtove' ) ) THEN - ALLOCATE( v_msf(A1Dj(nn_hls),jpk,nbasin) ) + ALLOCATE( v_msf(A1Dj(0),jpk,nbasin) ) DO jn = 1, nbasin v_msf(:,:,jn) = ptr_sjk( pvtr(:,:,:), btmsk34(:,:,jn) ) @@ -374,16 +374,16 @@ CONTAINS ENDIF ! i sum of j surface area, j surface area - temperature/salinity product on V grid - IF( iom_use( 'sopstove' ) .OR. iom_use( 'sophtove' ) .OR. & + IF( iom_use( 'sopstove' ) .OR. iom_use( 'sophtove' ) .OR. & & iom_use( 'sopstbtr' ) .OR. iom_use( 'sophtbtr' ) ) THEN - ALLOCATE( zmask(A2D(nn_hls),jpk), zts(A2D(nn_hls),jpk,jpts), & - & sjk(A1Dj(nn_hls),jpk,nbasin), & - & zt_jk(A1Dj(nn_hls),jpk,nbasin), zs_jk(A1Dj(nn_hls),jpk,nbasin) ) + ALLOCATE( zmask( A2D(0),jpk ), zts( A2D(0),jpk,jpts ), & + & sjk( A1Dj(0),jpk,nbasin), & + & zt_jk(A1Dj(0),jpk,nbasin), zs_jk(A1Dj(0),jpk,nbasin) ) zmask(:,:,:) = 0._wp zts(:,:,:,:) = 0._wp - DO_3D( 1, 1, 1, 0, 1, jpkm1 ) + DO_3D( 0, 0, 0, 0, 1, jpkm1 ) zvfc = e1v(ji,jj) * e3v(ji,jj,jk,Kmm) zmask(ji,jj,jk) = vmask(ji,jj,jk) * zvfc zts(ji,jj,jk,jp_tem) = (ts(ji,jj,jk,jp_tem,Kmm)+ts(ji,jj+1,jk,jp_tem,Kmm)) * 0.5 * zvfc !Tracers averaged onto V grid @@ -405,14 +405,14 @@ CONTAINS ELSE ! i sum of j surface area - temperature/salinity product on T grid IF( iom_use( 'zotem' ) .OR. iom_use( 'zosal' ) .OR. iom_use( 'zosrf' ) ) THEN - ALLOCATE( zmask(A2D(nn_hls),jpk), zts(A2D(nn_hls),jpk,jpts), & - & sjk(A1Dj(nn_hls),jpk,nbasin), & - & zt_jk(A1Dj(nn_hls),jpk,nbasin), zs_jk(A1Dj(nn_hls),jpk,nbasin) ) + ALLOCATE( zmask( A2D(0),jpk ), zts( A2D(0),jpk,jpts ), & + & sjk( A1Dj(0),jpk,nbasin), & + & zt_jk(A1Dj(0),jpk,nbasin), zs_jk(A1Dj(0),jpk,nbasin) ) zmask(:,:,:) = 0._wp zts(:,:,:,:) = 0._wp - DO_3D( 1, 1, 1, 1, 1, jpkm1 ) + DO_3D( 0, 0, 0, 0, 1, jpkm1 ) zsfc = e1t(ji,jj) * e3t(ji,jj,jk,Kmm) zmask(ji,jj,jk) = tmask(ji,jj,jk) * zsfc zts(ji,jj,jk,jp_tem) = ts(ji,jj,jk,jp_tem,Kmm) * zsfc @@ -434,11 +434,12 @@ CONTAINS ! i-k sum of j surface area - temperature/salinity product on V grid IF( iom_use( 'sopstvtr' ) .OR. iom_use( 'sophtvtr' ) ) THEN + ! TODO: Can be A2D(0) once all dia_ptr_hst calls have arguments with consistent declarations ALLOCATE( zts(A2D(nn_hls),jpk,jpts) ) zts(:,:,:,:) = 0._wp - DO_3D( 1, 1, 1, 0, 1, jpkm1 ) + DO_3D( 0, 0, 0, 0, 1, jpkm1 ) zvfc = e1v(ji,jj) * e3v(ji,jj,jk,Kmm) zts(ji,jj,jk,jp_tem) = (ts(ji,jj,jk,jp_tem,Kmm)+ts(ji,jj+1,jk,jp_tem,Kmm)) * 0.5 * zvfc !Tracers averaged onto V grid zts(ji,jj,jk,jp_sal) = (ts(ji,jj,jk,jp_sal,Kmm)+ts(ji,jj+1,jk,jp_sal,Kmm)) * 0.5 * zvfc @@ -538,11 +539,12 @@ CONTAINS !! Wrapper for heat and salt transport calculations to calculate them for each basin !! Called from all advection and/or diffusion routines !!---------------------------------------------------------------------- - INTEGER , INTENT(in ) :: ktra ! tracer index - CHARACTER(len=3) , INTENT(in) :: cptr ! transport type 'adv'/'ldf'/'eiv' - REAL(wp), DIMENSION(A2D(nn_hls),jpk) , INTENT(in) :: pvflx ! 3D input array of advection/diffusion - REAL(wp), DIMENSION(A1Dj(nn_hls),nbasin) :: zsj ! - INTEGER :: jn ! + INTEGER , INTENT(in) :: ktra ! tracer index + CHARACTER(len=3) , INTENT(in) :: cptr ! transport type 'adv'/'ldf'/'eiv' + ! TODO: Can be A2D(0) once all dia_ptr_hst calls have arguments with consistent declarations + REAL(wp), DIMENSION(A2D(nn_hls),jpk), INTENT(in) :: pvflx ! 3D input array of advection/diffusion + REAL(wp), DIMENSION(A1Dj(0),nbasin) :: zsj ! + INTEGER :: jn ! DO jn = 1, nbasin zsj(:,jn) = ptr_sj( pvflx(:,:,:), btmsk(:,:,jn) ) @@ -576,13 +578,13 @@ CONTAINS !! !! ** Action : phstr !!---------------------------------------------------------------------- - REAL(wp), DIMENSION(jpj,nbasin) , INTENT(inout) :: phstr ! - REAL(wp), DIMENSION(A1Dj(nn_hls),nbasin), INTENT(in) :: pva ! + REAL(wp), DIMENSION(Njs0:Nje0,nbasin), INTENT(inout) :: phstr ! + REAL(wp), DIMENSION(A1Dj(0) ,nbasin), INTENT(in ) :: pva ! INTEGER :: jj #if ! defined key_mpi_off - INTEGER, DIMENSION(1) :: ish1d - INTEGER, DIMENSION(2) :: ish2d - REAL(wp), DIMENSION(jpj*nbasin) :: zwork + INTEGER, DIMENSION(1) :: ish1d + INTEGER, DIMENSION(2) :: ish2d + REAL(wp), DIMENSION(:), ALLOCATABLE :: zwork #endif DO jj = ntsj, ntej @@ -591,11 +593,13 @@ CONTAINS #if ! defined key_mpi_off IF( .NOT. l_istiled .OR. ntile == nijtile ) THEN - ish1d(1) = jpj*nbasin - ish2d(1) = jpj ; ish2d(2) = nbasin + ALLOCATE( zwork(Nj_0*nbasin) ) + ish1d(1) = Nj_0*nbasin + ish2d(1) = Nj_0 ; ish2d(2) = nbasin zwork(:) = RESHAPE( phstr(:,:), ish1d ) CALL mpp_sum( 'diaptr', zwork, ish1d(1), ncomm_znl ) phstr(:,:) = RESHAPE( zwork, ish2d ) + DEALLOCATE( zwork ) ENDIF #endif END SUBROUTINE ptr_sum_2d @@ -612,13 +616,13 @@ CONTAINS !! !! ** Action : phstr !!---------------------------------------------------------------------- - REAL(wp), DIMENSION(jpj,jpk,nbasin) , INTENT(inout) :: phstr ! - REAL(wp), DIMENSION(A1Dj(nn_hls),jpk,nbasin), INTENT(in) :: pva ! - INTEGER :: jj, jk + REAL(wp), DIMENSION(Njs0:Nje0,jpk,nbasin), INTENT(inout) :: phstr ! + REAL(wp), DIMENSION(A1Dj(0) ,jpk,nbasin), INTENT(in ) :: pva ! + INTEGER :: jj, jk #if ! defined key_mpi_off - INTEGER, DIMENSION(1) :: ish1d - INTEGER, DIMENSION(3) :: ish3d - REAL(wp), DIMENSION(jpj*jpk*nbasin) :: zwork + INTEGER, DIMENSION(1) :: ish1d + INTEGER, DIMENSION(3) :: ish3d + REAL(wp), DIMENSION(:), ALLOCATABLE :: zwork #endif DO jk = 1, jpk @@ -629,11 +633,13 @@ CONTAINS #if ! defined key_mpi_off IF( .NOT. l_istiled .OR. ntile == nijtile ) THEN - ish1d(1) = jpj*jpk*nbasin - ish3d(1) = jpj ; ish3d(2) = jpk ; ish3d(3) = nbasin + ALLOCATE( zwork(Nj_0*jpk*nbasin) ) + ish1d(1) = Nj_0*jpk*nbasin + ish3d(1) = Nj_0 ; ish3d(2) = jpk ; ish3d(3) = nbasin zwork(:) = RESHAPE( phstr(:,:,:), ish1d ) CALL mpp_sum( 'diaptr', zwork, ish1d(1), ncomm_znl ) phstr(:,:,:) = RESHAPE( zwork, ish3d ) + DEALLOCATE( zwork ) ENDIF #endif END SUBROUTINE ptr_sum_3d @@ -651,13 +657,13 @@ CONTAINS ! nbasin has been initialized in iom_init to define the axis "basin" ! IF( .NOT. ALLOCATED( btmsk ) ) THEN - ALLOCATE( btmsk(jpi,jpj,nbasin) , btmsk34(jpi,jpj,nbasin), & - & hstr_adv(jpj,jpts,nbasin), hstr_eiv(jpj,jpts,nbasin), & - & hstr_ove(jpj,jpts,nbasin), hstr_btr(jpj,jpts,nbasin), & - & hstr_ldf(jpj,jpts,nbasin), hstr_vtr(jpj,jpts,nbasin), STAT=ierr(1) ) + ALLOCATE( btmsk(jpi,jpj,nbasin) , btmsk34(jpi,jpj,nbasin), & + & hstr_adv(Njs0:Nje0,jpts,nbasin), hstr_eiv(Njs0:Nje0,jpts,nbasin), & + & hstr_ove(Njs0:Nje0,jpts,nbasin), hstr_btr(Njs0:Nje0,jpts,nbasin), & + & hstr_ldf(Njs0:Nje0,jpts,nbasin), hstr_vtr(Njs0:Nje0,jpts,nbasin), STAT=ierr(1) ) ! - ALLOCATE( pvtr_int(jpj,jpk,jpts+2,nbasin), & - & pzon_int(jpj,jpk,jpts+1,nbasin), STAT=ierr(2) ) + ALLOCATE( pvtr_int(Njs0:Nje0,jpk,jpts+2,nbasin), & + & pzon_int(Njs0:Nje0,jpk,jpts+1,nbasin), STAT=ierr(2) ) ! dia_ptr_alloc = MAXVAL( ierr ) CALL mpp_sum( 'diaptr', dia_ptr_alloc ) @@ -677,11 +683,12 @@ CONTAINS !! !! ** Action : - p_fval: i-k-mean poleward flux of pvflx !!---------------------------------------------------------------------- + ! TODO: Can be A2D(0) once all dia_ptr_hst calls have arguments with consistent declarations REAL(wp), INTENT(in), DIMENSION(A2D(nn_hls),jpk) :: pvflx ! mask flux array at V-point - REAL(wp), INTENT(in), DIMENSION(jpi,jpj) :: pmsk ! Optional 2D basin mask + REAL(wp), INTENT(in), DIMENSION(jpi,jpj ) :: pmsk ! Optional 2D basin mask ! - INTEGER :: ji, jj, jk ! dummy loop arguments - REAL(wp), DIMENSION(A1Dj(nn_hls)) :: p_fval ! function value + INTEGER :: ji, jj, jk ! dummy loop arguments + REAL(wp), DIMENSION(A1Dj(0)) :: p_fval ! function value !!-------------------------------------------------------------------- ! p_fval(:) = 0._wp @@ -702,11 +709,12 @@ CONTAINS !! !! ** Action : - p_fval: i-k-mean poleward flux of pvflx !!---------------------------------------------------------------------- - REAL(wp) , INTENT(in), DIMENSION(A2D(nn_hls)) :: pvflx ! mask flux array at V-point - REAL(wp) , INTENT(in), DIMENSION(jpi,jpj) :: pmsk ! Optional 2D basin mask + ! TODO: Can be A2D(0) once all dia_ptr_hst calls have arguments with consistent declarations + REAL(wp) , INTENT(in), DIMENSION(A2D(nn_hls)) :: pvflx ! mask flux array at V-point + REAL(wp) , INTENT(in), DIMENSION(jpi,jpj ) :: pmsk ! Optional 2D basin mask ! - INTEGER :: ji,jj ! dummy loop arguments - REAL(wp), DIMENSION(A1Dj(nn_hls)) :: p_fval ! function value + INTEGER :: ji, jj ! dummy loop arguments + REAL(wp), DIMENSION(A1Dj(0)) :: p_fval ! function value !!-------------------------------------------------------------------- ! p_fval(:) = 0._wp @@ -725,14 +733,13 @@ CONTAINS !! !! ** Action : - p_fval: j-cumulated sum of pva !!---------------------------------------------------------------------- - REAL(wp) , INTENT(in), DIMENSION(jpi,jpj) :: pva ! mask flux array at V-point + REAL(wp) , INTENT(in), DIMENSION(A2D(0)) :: pva ! mask flux array at V-point ! - INTEGER :: ji,jj,jc ! dummy loop arguments - INTEGER :: ijpj ! ??? - REAL(wp), DIMENSION(jpi,jpj) :: p_fval ! function value + INTEGER :: ji,jj,jc ! dummy loop arguments + INTEGER :: ijpj ! ??? + REAL(wp), DIMENSION(A2D(0)) :: p_fval ! function value !!-------------------------------------------------------------------- ! - ijpj = jpj ! ??? p_fval(:,:) = 0._wp DO jc = 1, jpnj ! looping over all processors in j axis DO_2D( 0, 0, 0, 0 ) @@ -756,11 +763,11 @@ CONTAINS !!---------------------------------------------------------------------- !! IMPLICIT none - REAL(wp) , INTENT(in), DIMENSION(A2D(nn_hls),jpk) :: pta ! mask flux array at V-point - REAL(wp) , INTENT(in), DIMENSION(jpi,jpj) :: pmsk ! Optional 2D basin mask + REAL(wp) , INTENT(in), DIMENSION(A2D(0) ,jpk) :: pta ! mask flux array at V-point + REAL(wp) , INTENT(in), DIMENSION(jpi,jpj ) :: pmsk ! Optional 2D basin mask !! INTEGER :: ji, jj, jk ! dummy loop arguments - REAL(wp), DIMENSION(A1Dj(nn_hls),jpk) :: p_fval ! return function value + REAL(wp), DIMENSION(A1Dj(0),jpk) :: p_fval ! return function value !!-------------------------------------------------------------------- ! p_fval(:,:) = 0._wp diff --git a/src/OCE/DIA/diawri.F90 b/src/OCE/DIA/diawri.F90 index 0deeb3370..b6dd96cda 100644 --- a/src/OCE/DIA/diawri.F90 +++ b/src/OCE/DIA/diawri.F90 @@ -135,8 +135,8 @@ CONTAINS ENDIF ! initialize arrays - z2d(:,:) = 0._wp - z3d(:,:,:) = 0._wp + z2d(A2D(0)) = 0._wp + z3d(A2D(0),:) = 0._wp ! Output of initial vertical scale factor CALL iom_put("e3t_0", e3t_0(:,:,:) ) diff --git a/src/OCE/IOM/iom.F90 b/src/OCE/IOM/iom.F90 index d66b13d19..728c75ebc 100644 --- a/src/OCE/IOM/iom.F90 +++ b/src/OCE/IOM/iom.F90 @@ -1327,14 +1327,21 @@ CONTAINS IF( irankpv == 1 ) ishape(1:1) = SHAPE(pv_r1d) IF( irankpv == 2 ) ishape(1:2) = SHAPE(pv_r2d) IF( irankpv == 3 ) ishape(1:3) = SHAPE(pv_r3d) + ix1 = 1 ; ix2 = icnt(1) ; iy1 = 1 ; iy2 = icnt(2) ! index of the array to be read ctmp1 = 'd' ELSE - IF( irankpv == 2 ) THEN - ishape(1:2) = SHAPE(pv_r2d(Nis0:Nie0,Njs0:Nje0 )) ; ctmp1 = 'd(Nis0:Nie0,Njs0:Nje0)' - ENDIF - IF( irankpv == 3 ) THEN - ishape(1:3) = SHAPE(pv_r3d(Nis0:Nie0,Njs0:Nje0,:)) ; ctmp1 = 'd(Nis0:Nie0,Njs0:Nje0,:)' + IF( irankpv == 2 ) ishape(1:2) = SHAPE(pv_r2d) + IF( irankpv == 3 ) ishape(1:3) = SHAPE(pv_r3d) + IF( ishape(1) == jpi .AND. ishape(2) == jpj ) THEN + ishape(1:2) = (/ Ni_0, Nj_0 /) + ix1 = Nis0 ; ix2 = Nie0 ; iy1 = Njs0 ; iy2 = Nje0 ! index of the array to be read + ctmp1 = 'd(Nis0:Nie0,Njs0:Nje0' + ELSE + ix1 = 1 ; ix2 = icnt(1) ; iy1 = 1 ; iy2 = icnt(2) ! index of the array to be read + ctmp1 = 'd(:,:' ENDIF + IF( irankpv == 3 ) ctmp1 = TRIM(ctmp1)//',:' + ctmp1 = TRIM(ctmp1)//')' ENDIF DO jl = 1, irankpv WRITE( ctmp2, FMT="(', ', i1,'): ', i5,' /= icnt(', i1,'):', i5)" ) jl, ishape(jl), jl, icnt(jl) @@ -1347,11 +1354,6 @@ CONTAINS !- IF( idvar > 0 .AND. istop == nstop ) THEN ! no additional errors until this point... ! - ! find the right index of the array to be read - IF( idom /= jpdom_unknown ) THEN ; ix1 = Nis0 ; ix2 = Nie0 ; iy1 = Njs0 ; iy2 = Nje0 - ELSE ; ix1 = 1 ; ix2 = icnt(1) ; iy1 = 1 ; iy2 = icnt(2) - ENDIF - CALL iom_nf90_get( kiomid, idvar, inbdim, istart, icnt, ix1, ix2, iy1, iy2, pv_r1d, pv_r2d, pv_r3d ) IF( istop == nstop ) THEN ! no additional errors until this point... @@ -1362,9 +1364,11 @@ CONTAINS zsgn = 1._wp IF( PRESENT(psgn ) ) zsgn = psgn !--- overlap areas and extra hallows (mpp) - IF( PRESENT(pv_r2d) .AND. idom /= jpdom_unknown .AND. cl_type /= 'Z' ) THEN + llok = idom /= jpdom_unknown .AND. cl_type /= 'Z' & + & .AND. ix1 == Nis0 .AND. ix2 == Nie0 .AND. iy1 == Njs0 .AND. iy2 == Nje0 + IF( PRESENT(pv_r2d) .AND. llok ) THEN CALL lbc_lnk( 'iom', pv_r2d, cl_type, zsgn, kfillmode = kfill ) - ELSEIF( PRESENT(pv_r3d) .AND. idom /= jpdom_unknown .AND. cl_type /= 'Z' ) THEN + ELSEIF( PRESENT(pv_r3d) .AND. llok ) THEN CALL lbc_lnk( 'iom', pv_r3d, cl_type, zsgn, kfillmode = kfill ) ENDIF ! @@ -2336,14 +2340,14 @@ CONTAINS idb(jn) = -nn_hls ! Tile data offset (halo size) END DO - ! Tile_[ij]begin are defined with respect to the processor data domain, so data_[ij]begin is added CALL iom_set_domain_attr("grid_"//cdgrd, ntiles=nijtile, & - & tile_ibegin=ntsi_a(1:nijtile) + idb(:) - 1, tile_jbegin=ntsj_a(1:nijtile) + idb(:) - 1, & + & tile_ibegin=ntsi_a(1:nijtile) - nn_hls - 1, tile_jbegin=ntsj_a(1:nijtile) - nn_hls - 1, & & tile_ni=ini(:), tile_nj=inj(:), & & tile_data_ibegin=idb(:), tile_data_jbegin=idb(:), & & tile_data_ni=ini(:) - 2 * idb(:), tile_data_nj=inj(:) - 2 * idb(:)) + idb(:) = 0 CALL iom_set_domain_attr("grid_"//cdgrd//"_inner", ntiles=nijtile, & - & tile_ibegin=ntsi_a(1:nijtile) + idb(:) - 1, tile_jbegin=ntsj_a(1:nijtile) + idb(:) - 1, & + & tile_ibegin=ntsi_a(1:nijtile) - nn_hls - 1, tile_jbegin=ntsj_a(1:nijtile) - nn_hls - 1, & & tile_ni=ini(:), tile_nj=inj(:), & & tile_data_ibegin=idb(:), tile_data_jbegin=idb(:), & & tile_data_ni=ini(:) - 2 * idb(:), tile_data_nj=inj(:) - 2 * idb(:)) @@ -2453,7 +2457,7 @@ CONTAINS ! CALL dom_ngb( -168.53_wp, 65.03_wp, ix, iy, 'T' ) ! i-line that passes through Bering Strait: Reference latitude (used in plots) CALL dom_ngb( 180.0_wp, 90.0_wp, ix, iy, 'T' ) ! i-line that passes near the North Pole : Reference latitude (used in plots) CALL iom_set_domain_attr("gznl", ni_glo=Ni0glo, nj_glo=Nj0glo, ibegin=mig0(Nis0)-1, jbegin=mjg0(Njs0)-1, ni=Ni_0, nj=Nj_0) - CALL iom_set_domain_attr("gznl", data_dim=2, data_ibegin = -nn_hls, data_ni = jpi, data_jbegin = -nn_hls, data_nj = jpj) + CALL iom_set_domain_attr("gznl", data_dim=2, data_ibegin=0, data_ni=Ni_0, data_jbegin=0, data_nj=Nj_0) CALL iom_set_domain_attr("gznl", lonvalue = real(zlon, dp), & & latvalue = real(RESHAPE(plat(Nis0:Nie0, Njs0:Nje0),(/ Ni_0*Nj_0 /)),dp)) CALL iom_set_zoom_domain_attr("ptr", ibegin=ix-1, jbegin=0, ni=1, nj=Nj0glo) diff --git a/src/OCE/IOM/iom_nf90.F90 b/src/OCE/IOM/iom_nf90.F90 index c49d9538e..287d943a9 100644 --- a/src/OCE/IOM/iom_nf90.F90 +++ b/src/OCE/IOM/iom_nf90.F90 @@ -40,6 +40,8 @@ MODULE iom_nf90 MODULE PROCEDURE iom_nf90_rp0123d_dp END INTERFACE + !! * Substitutions +# include "do_loop_substitute.h90" !!---------------------------------------------------------------------- !! NEMO/OCE 4.0 , NEMO Consortium (2018) !! $Id: iom_nf90.F90 14433 2021-02-11 08:06:49Z smasson $ @@ -544,7 +546,7 @@ CONTAINS INTEGER :: idvar ! variable id INTEGER :: jd ! dimension loop counter INTEGER :: ix1, ix2, iy1, iy2 ! subdomain indexes - INTEGER, DIMENSION(4) :: idimsz ! dimensions size + INTEGER, DIMENSION(3) :: ishape ! dimensions size INTEGER, DIMENSION(4) :: idimid ! dimensions id CHARACTER(LEN=256) :: clinfo ! info character INTEGER :: if90id ! nf90 file identifier @@ -627,11 +629,9 @@ CONTAINS itype = NF90_DOUBLE ENDIF IF( PRESENT(pv_r0d) ) THEN - CALL iom_nf90_check(NF90_DEF_VAR( if90id, TRIM(cdvar), itype, & - & iom_file(kiomid)%nvid(idvar) ), clinfo ) + CALL iom_nf90_check(NF90_DEF_VAR( if90id, TRIM(cdvar), itype, iom_file(kiomid)%nvid(idvar) ), clinfo ) ELSE - CALL iom_nf90_check(NF90_DEF_VAR( if90id, TRIM(cdvar), itype, idimid(1:idims), & - & iom_file(kiomid)%nvid(idvar) ), clinfo ) + CALL iom_nf90_check(NF90_DEF_VAR( if90id, TRIM(cdvar), itype, idimid(1:idims), iom_file(kiomid)%nvid(idvar) ), clinfo ) ENDIF lchunk = .false. IF( snc4set%luse .AND. idims == 4 ) lchunk = .true. @@ -673,23 +673,13 @@ CONTAINS ENDIF ! on what kind of domain must the data be written? IF( PRESENT(pv_r2d) .OR. PRESENT(pv_r3d) ) THEN - idimsz(1:2) = iom_file(kiomid)%dimsz(1:2,idvar) - IF( idimsz(1) == Ni_0 .AND. idimsz(2) == Nj_0 ) THEN - ix1 = Nis0 ; ix2 = Nie0 ; iy1 = Njs0 ; iy2 = Nje0 - ELSEIF( idimsz(1) == jpi .AND. idimsz(2) == jpj ) THEN - ix1 = 1 ; ix2 = jpi ; iy1 = 1 ; iy2 = jpj - ELSEIF( idimsz(1) == jpi .AND. idimsz(2) == jpj ) THEN - ix1 = 1 ; ix2 = jpi ; iy1 = 1 ; iy2 = jpj - ELSE - CALL ctl_stop( 'iom_nf90_rp0123d: should have been an impossible case...' ) - ENDIF ! write dimension variables if it is not already done ! ============= ! trick: is defined to 0 => dimension variable are defined but not yet written IF( iom_file(kiomid)%dimsz(1, 4) == 0 ) THEN ! time_counter = 0 - CALL iom_nf90_check( NF90_PUT_VAR( if90id, 1, glamt(ix1:ix2, iy1:iy2) ), clinfo ) - CALL iom_nf90_check( NF90_PUT_VAR( if90id, 2, gphit(ix1:ix2, iy1:iy2) ), clinfo ) + CALL iom_nf90_check( NF90_PUT_VAR( if90id, 1, glamt(A2D(0)) ), clinfo ) + CALL iom_nf90_check( NF90_PUT_VAR( if90id, 2, gphit(A2D(0)) ), clinfo ) SELECT CASE (iom_file(kiomid)%comp) CASE ('OCE') CALL iom_nf90_check( NF90_PUT_VAR( if90id, 3, gdept_1d ), clinfo ) @@ -704,6 +694,17 @@ CONTAINS iom_file(kiomid)%dimsz(1, 4) = 1 ! so we don't enter this IF case any more... IF(lwp) WRITE(numout,*) TRIM(clinfo)//' write dimension variables done' ENDIF + + IF( PRESENT(pv_r2d) ) ishape(1:2) = SHAPE(pv_r2d) + IF( PRESENT(pv_r3d) ) ishape(1:3) = SHAPE(pv_r3d) + IF( ishape(1) == Ni_0 .AND. ishape(2) == Nj_0 ) THEN + ix1 = 1 ; ix2 = Ni_0 ; iy1 = 1 ; iy2 = Nj_0 + ELSEIF( ishape(1) == jpi .AND. ishape(2) == jpj ) THEN + ix1 = Nis0 ; ix2 = Nie0 ; iy1 = Njs0 ; iy2 = Nje0 + ELSE + CALL ctl_stop( 'iom_nf90_rp0123d: should have been an impossible case...' ) + ENDIF + ENDIF ! write the data @@ -712,7 +713,7 @@ CONTAINS CALL iom_nf90_check( NF90_PUT_VAR( if90id, idvar, pv_r0d ), clinfo ) ELSEIF( PRESENT(pv_r1d) ) THEN CALL iom_nf90_check( NF90_PUT_VAR( if90id, idvar, pv_r1d(:) ), clinfo ) - ELSEIF( PRESENT(pv_r2d) ) THEN + ELSEIF( PRESENT(pv_r2d) ) THEN CALL iom_nf90_check( NF90_PUT_VAR( if90id, idvar, pv_r2d(ix1:ix2,iy1:iy2) ), clinfo ) ELSEIF( PRESENT(pv_r3d) ) THEN CALL iom_nf90_check( NF90_PUT_VAR( if90id, idvar, pv_r3d(ix1:ix2,iy1:iy2,:) ), clinfo ) diff --git a/src/OCE/IOM/prtctl.F90 b/src/OCE/IOM/prtctl.F90 index fae7c0d80..6fcca35f9 100644 --- a/src/OCE/IOM/prtctl.F90 +++ b/src/OCE/IOM/prtctl.F90 @@ -137,9 +137,14 @@ CONTAINS INTEGER :: jn, jl, kdir INTEGER :: iis, iie, jjs, jje INTEGER :: itra, inum + INTEGER, DIMENSION(4) :: ishape REAL(2*wp) :: zsum1, zsum2, zvctl1, zvctl2 !!---------------------------------------------------------------------- ! + IF( ( ktab2d_1 * ktab3d_1 * ktab4d_1 * ktab2d_2 * ktab3d_2 ) /= 0 ) THEN + CALL ctl_stop( 'prt_ctl is not working with tiles' ) + ENDIF + ! Arrays, scalars initialization cl1 = '' cl2 = '' @@ -157,12 +162,19 @@ CONTAINS ! Loop over each sub-domain, i.e. the total number of processors ijsplt DO jl = 1, SIZE(nall_ictls) - - ! define shoter names... - iis = MAX( nall_ictls(jl), ntsi ) - iie = MIN( nall_ictle(jl), ntei ) - jjs = MAX( nall_jctls(jl), ntsj ) - jje = MIN( nall_jctle(jl), ntej ) + + IF( PRESENT(tab2d_1) ) ishape(1:2) = SHAPE(tab2d_1) + IF( PRESENT(tab3d_1) ) ishape(1:3) = SHAPE(tab3d_1) + IF( PRESENT(tab4d_1) ) ishape(1:4) = SHAPE(tab4d_1) + IF( ishape(1) == jpi .AND. ishape(2) == jpj ) THEN + iis = Nis0 ; iie = Nie0 ; jjs = Njs0 ; jje = Nje0 + ELSE + iis = 1 ; iie = ishape(1) ; jjs = 1 ; jje = ishape(2) + ENDIF + iis = MAX( nall_ictls(jl), iis ) + iie = MIN( nall_ictle(jl), iie ) + jjs = MAX( nall_jctls(jl), jjs ) + jje = MIN( nall_jctle(jl), jje ) IF( PRESENT(clinfo) ) THEN ; inum = numprt_top(jl) ELSE ; inum = numprt_oce(jl) @@ -188,32 +200,32 @@ CONTAINS ! 2D arrays IF( PRESENT(tab2d_1) ) THEN - IF( PRESENT(mask1) ) THEN ; zsum1 = SUM( tab2d_1(iis:iie,jjs:jje) * mask1(iis:iie,jjs:jje,1) ) - ELSE ; zsum1 = SUM( tab2d_1(iis:iie,jjs:jje) ) + IF( PRESENT(mask1) ) THEN ; zsum1 = SUM( tab2d_1(iis:iie,jjs:jje) * mask1(A2D(0),1) ) + ELSE ; zsum1 = SUM( tab2d_1(iis:iie,jjs:jje) ) ENDIF ENDIF IF( PRESENT(tab2d_2) ) THEN - IF( PRESENT(mask2) ) THEN ; zsum2 = SUM( tab2d_2(iis:iie,jjs:jje) * mask2(iis:iie,jjs:jje,1) ) - ELSE ; zsum2 = SUM( tab2d_2(iis:iie,jjs:jje) ) + IF( PRESENT(mask2) ) THEN ; zsum2 = SUM( tab2d_2(iis:iie,jjs:jje) * mask2(A2D(0),1) ) + ELSE ; zsum2 = SUM( tab2d_2(iis:iie,jjs:jje) ) ENDIF ENDIF ! 3D arrays IF( PRESENT(tab3d_1) ) THEN - IF( PRESENT(mask1) ) THEN ; zsum1 = SUM( tab3d_1(iis:iie,jjs:jje,1:kdir) * mask1(iis:iie,jjs:jje,1:kdir) ) - ELSE ; zsum1 = SUM( tab3d_1(iis:iie,jjs:jje,1:kdir) ) + IF( PRESENT(mask1) ) THEN ; zsum1 = SUM( tab3d_1(iis:iie,jjs:jje,1:kdir) * mask1(A2D(0),1:kdir) ) + ELSE ; zsum1 = SUM( tab3d_1(iis:iie,jjs:jje,1:kdir) ) ENDIF ENDIF IF( PRESENT(tab3d_2) ) THEN - IF( PRESENT(mask2) ) THEN ; zsum2 = SUM( tab3d_2(iis:iie,jjs:jje,1:kdir) * mask2(iis:iie,jjs:jje,1:kdir) ) - ELSE ; zsum2 = SUM( tab3d_2(iis:iie,jjs:jje,1:kdir) ) + IF( PRESENT(mask2) ) THEN ; zsum2 = SUM( tab3d_2(iis:iie,jjs:jje,1:kdir) * mask2(A2D(0),1:kdir) ) + ELSE ; zsum2 = SUM( tab3d_2(iis:iie,jjs:jje,1:kdir) ) ENDIF ENDIF ! 4D arrays IF( PRESENT(tab4d_1) ) THEN - IF( PRESENT(mask1) ) THEN ; zsum1 = SUM( tab4d_1(iis:iie,jjs:jje,1:kdir,jn) * mask1(iis:iie,jjs:jje,1:kdir) ) - ELSE ; zsum1 = SUM( tab4d_1(iis:iie,jjs:jje,1:kdir,jn) ) + IF( PRESENT(mask1) ) THEN ; zsum1 = SUM( tab4d_1(iis:iie,jjs:jje,1:kdir,jn) * mask1(A2D(0),1:kdir) ) + ELSE ; zsum1 = SUM( tab4d_1(iis:iie,jjs:jje,1:kdir,jn) ) ENDIF ENDIF diff --git a/src/OCE/TRA/eosbn2.F90 b/src/OCE/TRA/eosbn2.F90 index 81fda25e0..1f808dfb9 100644 --- a/src/OCE/TRA/eosbn2.F90 +++ b/src/OCE/TRA/eosbn2.F90 @@ -1336,10 +1336,10 @@ CONTAINS !! pab_pe(:,:,:,jp_tem) is alpha_pe !! pab_pe(:,:,:,jp_sal) is beta_pe !!---------------------------------------------------------------------- - INTEGER , INTENT(in ) :: Kmm ! time level index - REAL(wp), DIMENSION(jpi,jpj,jpk,jpts), INTENT(in ) :: pts ! pot. temperature & salinity - REAL(wp), DIMENSION(jpi,jpj,jpk,jpts), INTENT( out) :: pab_pe ! alpha_pe and beta_pe - REAL(wp), DIMENSION(jpi,jpj,jpk) , INTENT( out) :: ppen ! potential energy anomaly + INTEGER , INTENT(in ) :: Kmm ! time level index + REAL(wp), DIMENSION(:,:,:,:), INTENT(in ) :: pts ! pot. temperature & salinity + REAL(wp), DIMENSION(:,:,:,:), INTENT( out) :: pab_pe ! alpha_pe and beta_pe + REAL(wp), DIMENSION(:,:,:) , INTENT( out) :: ppen ! potential energy anomaly ! INTEGER :: ji, jj, jk ! dummy loop indices REAL(wp) :: zt , zh , zs , ztm ! local scalars @@ -1352,7 +1352,7 @@ CONTAINS ! CASE( np_teos10, np_eos80 ) !== polynomial TEOS-10 / EOS-80 ==! ! - DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1 ) + DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! zh = gdept(ji,jj,jk,Kmm) * r1_Z0 ! depth zt = pts (ji,jj,jk,jp_tem) * r1_T0 ! temperature @@ -1411,9 +1411,9 @@ CONTAINS ! CASE( np_seos ) !== Vallis (2006) simplified EOS ==! ! - DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1 ) + DO_3D( 0, 0, 0, 0, 1, jpkm1 ) zt = pts(ji,jj,jk,jp_tem) - 10._wp ! temperature anomaly (t-T0) - zs = pts (ji,jj,jk,jp_sal) - 35._wp ! abs. salinity anomaly (s-S0) + zs = pts (ji,jj,jk,jp_sal) - 35._wp ! abs. salinity anomaly (s-S0) zh = gdept(ji,jj,jk,Kmm) ! depth in meters at t-point ztm = tmask(ji,jj,jk) ! tmask zn = 0.5_wp * zh * r1_rho0 * ztm diff --git a/src/OCE/TRA/traadv.F90 b/src/OCE/TRA/traadv.F90 index a55ffff58..190a6a7ad 100644 --- a/src/OCE/TRA/traadv.F90 +++ b/src/OCE/TRA/traadv.F90 @@ -178,7 +178,7 @@ CONTAINS ! !!gm ??? ! TEMP: [tiling] This copy-in not necessary after all lbc_lnks removed in the nn_hls = 2 case in tra_adv_fct - IF( l_diaptr ) CALL dia_ptr( kt, Kmm, zvv(A2D(nn_hls),:) ) ! diagnose the effective MSF + IF( l_diaptr ) CALL dia_ptr( kt, Kmm, zvv(A2D(0),:) ) ! diagnose the effective MSF !!gm ??? ! diff --git a/src/OCE/TRD/trddyn.F90 b/src/OCE/TRD/trddyn.F90 index 32e95a865..dd2daf7ef 100644 --- a/src/OCE/TRD/trddyn.F90 +++ b/src/OCE/TRD/trddyn.F90 @@ -56,10 +56,13 @@ CONTAINS INTEGER , INTENT(in ) :: ktrd ! trend index INTEGER , INTENT(in ) :: kt ! time step INTEGER , INTENT(in ) :: Kmm ! time level index + INTEGER :: ji, jj, jk ! lopp indices !!---------------------------------------------------------------------- ! - putrd(:,:,:) = putrd(:,:,:) * umask(:,:,:) ! mask the trends - pvtrd(:,:,:) = pvtrd(:,:,:) * vmask(:,:,:) + DO_3D( 0, 0, 0, 0, 1, jpkm1 ) + putrd(ji,jj,jk) = putrd(ji,jj,jk) * umask(ji,jj,jk) ! mask the trends + pvtrd(ji,jj,jk) = pvtrd(ji,jj,jk) * vmask(ji,jj,jk) + END_3D ! !!gm NB : here a lbc_lnk should probably be added @@ -120,10 +123,10 @@ CONTAINS CALL iom_put( "vtrd_rvo", pvtrd ) CASE( jpdyn_keg ) ; CALL iom_put( "utrd_keg", putrd ) ! Kinetic Energy gradient (or had) CALL iom_put( "vtrd_keg", pvtrd ) - ALLOCATE( z3dx(jpi,jpj,jpk) , z3dy(jpi,jpj,jpk) ) - z3dx(:,:,:) = 0._wp ! U.dxU & V.dyV (approximation) - z3dy(:,:,:) = 0._wp - DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! no mask as un,vn are masked + ALLOCATE( z3dx(A2D(0),jpk) , z3dy(A2D(0),jpk) ) ! U.dxU & V.dyV (approximation) + z3dx(A2D(0),jpk) = 0._wp + z3dy(A2D(0),jpk) = 0._wp + DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! no mask as un,vn are masked z3dx(ji,jj,jk) = uu(ji,jj,jk,Kmm) * ( uu(ji+1,jj,jk,Kmm) - uu(ji-1,jj,jk,Kmm) ) / ( 2._wp * e1u(ji,jj) ) z3dy(ji,jj,jk) = vv(ji,jj,jk,Kmm) * ( vv(ji,jj+1,jk,Kmm) - vv(ji,jj-1,jk,Kmm) ) / ( 2._wp * e2v(ji,jj) ) END_3D @@ -139,7 +142,7 @@ CONTAINS CALL iom_put( "vtrd_zdf", pvtrd ) ! ! ! wind stress trends - ALLOCATE( z2dx(jpi,jpj) , z2dy(jpi,jpj) ) + ALLOCATE( z2dx(A2D(0)) , z2dy(A2D(0)) ) DO_2D( 0, 0, 0, 0 ) z2dx(ji,jj) = ( utau_b(ji,jj) + utauU(ji,jj) ) / ( e3u(ji,jj,1,Kmm) * rho0 ) z2dy(ji,jj) = ( vtau_b(ji,jj) + vtauV(ji,jj) ) / ( e3v(ji,jj,1,Kmm) * rho0 ) diff --git a/src/OCE/TRD/trdglo.F90 b/src/OCE/TRD/trdglo.F90 index a19770158..40c7b05b9 100644 --- a/src/OCE/TRD/trdglo.F90 +++ b/src/OCE/TRD/trdglo.F90 @@ -42,6 +42,7 @@ MODULE trdglo REAL(wp) :: tvolv ! volume of the whole ocean computed at v-points REAL(wp) :: rpktrd ! potential to kinetic energy conversion REAL(wp) :: peke ! conversion potential energy - kinetic energy trend + REAL(wp) :: xcof ! !!! domain averaged trends REAL(wp), DIMENSION(jptot_tra) :: tmo, smo ! temperature and salinity trends @@ -77,44 +78,47 @@ CONTAINS INTEGER :: ji, jj, jk ! dummy loop indices INTEGER :: ikbu, ikbv ! local integers REAL(wp):: zvm, zvt, zvs, z1_2rho0 ! local scalars - REAL(wp), DIMENSION(jpi,jpj) :: ztswu, ztswv, z2dx, z2dy ! 2D workspace + REAL(wp), DIMENSION(A2D(0)) :: ztswu, ztswv, z2dx, z2dy ! 2D workspace !!---------------------------------------------------------------------- ! IF( MOD(kt,nn_trd) == 0 .OR. kt == nit000 .OR. kt == nitend ) THEN ! - SELECT CASE( ctype ) - ! - CASE( 'TRA' ) !== Tracers (T & S) ==! - DO_3D( 1, 1, 1, 1, 1, jpkm1 ) ! global sum of mask volume trend and trend*T (including interior mask) - zvm = e1e2t(ji,jj) * e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk) * tmask_i(ji,jj) - zvt = ptrdx(ji,jj,jk) * zvm - zvs = ptrdy(ji,jj,jk) * zvm - tmo(ktrd) = tmo(ktrd) + zvt - smo(ktrd) = smo(ktrd) + zvs - t2 (ktrd) = t2(ktrd) + zvt * ts(ji,jj,jk,jp_tem,Kmm) - s2 (ktrd) = s2(ktrd) + zvs * ts(ji,jj,jk,jp_sal,Kmm) - END_3D + SELECT CASE( ctype ) + ! + CASE( 'TRA' ) !== Tracers (T & S) ==! + DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! global sum of mask volume trend and trend*T (including interior mask) + zvm = e1e2t(ji,jj) * e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk) * tmask_i(ji,jj) + zvt = ptrdx(ji,jj,jk) * zvm + zvs = ptrdy(ji,jj,jk) * zvm + tmo(ktrd) = tmo(ktrd) + zvt + smo(ktrd) = smo(ktrd) + zvs + t2 (ktrd) = t2(ktrd) + zvt * ts(ji,jj,jk,jp_tem,Kmm) + s2 (ktrd) = s2(ktrd) + zvs * ts(ji,jj,jk,jp_sal,Kmm) + END_3D ! ! linear free surface: diagnose advective flux trough the fixed k=1 w-surface - IF( ln_linssh .AND. ktrd == jptra_zad ) THEN - tmo(jptra_sad) = SUM( ww(:,:,1) * ts(:,:,1,jp_tem,Kmm) * e1e2t(:,:) * tmask_i(:,:) ) - smo(jptra_sad) = SUM( ww(:,:,1) * ts(:,:,1,jp_sal,Kmm) * e1e2t(:,:) * tmask_i(:,:) ) - t2 (jptra_sad) = SUM( ww(:,:,1) * ts(:,:,1,jp_tem,Kmm) * ts(:,:,1,jp_tem,Kmm) * e1e2t(:,:) * tmask_i(:,:) ) - s2 (jptra_sad) = SUM( ww(:,:,1) * ts(:,:,1,jp_sal,Kmm) * ts(:,:,1,jp_sal,Kmm) * e1e2t(:,:) * tmask_i(:,:) ) - ENDIF + IF( ln_linssh .AND. ktrd == jptra_zad ) THEN + DO_2D( 0, 0, 0, 0 ) ! global sum of mask volume trend and trend*T (including interior mask) + zvm = ww(ji,jj,1) * e1e2t(ji,jj) * tmask_i(ji,jj) + tmo(jptra_sad) = tmo(jptra_sad) + ts(ji,jj,1,jp_tem,Kmm) * zvm + smo(jptra_sad) = smo(jptra_sad) + ts(ji,jj,1,jp_sal,Kmm) * zvm + t2 (jptra_sad) = t2 (jptra_sad) + ts(ji,jj,1,jp_tem,Kmm) * ts(ji,jj,1,jp_tem,Kmm) * zvm + s2 (jptra_sad) = s2 (jptra_sad) + ts(ji,jj,1,jp_sal,Kmm) * ts(ji,jj,1,jp_sal,Kmm) * zvm + END_2D + ENDIF ! - IF( ktrd == jptra_atf ) THEN ! last trend (asselin time filter) - ! - CALL glo_tra_wri( kt ) ! print the results in ocean.output - ! - tmo(:) = 0._wp ! prepare the next time step (domain averaged array reset to zero) - smo(:) = 0._wp - t2 (:) = 0._wp - s2 (:) = 0._wp - ! - ENDIF + IF( ktrd == jptra_atf ) THEN ! last trend (asselin time filter) + ! + CALL glo_tra_wri( kt ) ! print the results in ocean.output + ! + tmo(:) = 0._wp ! prepare the next time step (domain averaged array reset to zero) + smo(:) = 0._wp + t2 (:) = 0._wp + s2 (:) = 0._wp + ! + ENDIF ! CASE( 'DYN' ) !== Momentum and KE ==! - DO_3D( 1, 0, 1, 0, 1, jpkm1 ) + DO_3D( 0, 0, 0, 0, 1, jpkm1 ) zvt = ptrdx(ji,jj,jk) * tmask_i(ji+1,jj) * tmask_i(ji,jj) * umask(ji,jj,jk) & & * e1e2u (ji,jj) * e3u(ji,jj,jk,Kmm) zvs = ptrdy(ji,jj,jk) * tmask_i(ji,jj+1) * tmask_i(ji,jj) * vmask(ji,jj,jk) & @@ -126,7 +130,7 @@ CONTAINS ! IF( ktrd == jpdyn_zdf ) THEN ! zdf trend: compute separately the surface forcing trend z1_2rho0 = 0.5_wp / rho0 - DO_2D( 1, 0, 1, 0 ) + DO_2D( 0, 0, 0, 0 ) zvt = ( utau_b(ji,jj) + utauU(ji,jj) ) * tmask_i(ji+1,jj) * tmask_i(ji,jj) * umask(ji,jj,jk) & & * z1_2rho0 * e1e2u(ji,jj) zvs = ( vtau_b(ji,jj) + vtauV(ji,jj) ) * tmask_i(ji,jj+1) * tmask_i(ji,jj) * vmask(ji,jj,jk) & @@ -184,8 +188,7 @@ CONTAINS INTEGER, INTENT(in) :: Kmm ! time level index ! INTEGER :: ji, jj, jk ! dummy loop indices - REAL(wp) :: zcof ! local scalar - REAL(wp), DIMENSION(jpi,jpj,jpk) :: zkx, zky, zkz, zkepe + REAL(wp), DIMENSION(A2D(0),jpk) :: zkx, zky, zkz, zkepe !!---------------------------------------------------------------------- ! I. Momentum trends @@ -196,24 +199,24 @@ CONTAINS ! I.1 Conversion potential energy - kinetic energy ! -------------------------------------------------- ! c a u t i o n here, trends are computed at kt+1 (now , but after the swap) - zkx (:,:,:) = 0._wp - zky (:,:,:) = 0._wp - zkz (:,:,:) = 0._wp - zkepe(:,:,:) = 0._wp + DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! loop from top to bottom + zkx (ji,jj,jk) = 0._wp + zky (ji,jj,jk) = 0._wp + zkz (ji,jj,jk) = 0._wp + zkepe(ji,jj,jk) = 0._wp + END_3D CALL eos( ts(:,:,:,:,Kmm), rhd, rhop ) ! now potential density - zcof = 0.5_wp / rho0 ! Density flux at w-point - zkz(:,:,1) = 0._wp - DO jk = 2, jpk - zkz(:,:,jk) = zcof * e1e2t(:,:) * ww(:,:,jk) * ( rhop(:,:,jk) + rhop(:,:,jk-1) ) * tmask_i(:,:) - END DO + zkz(A2D(0),1) = 0._wp + DO_3D( 0, 0, 0, 0, 1, jpk ) ! loop from top to bottom + zkz(ji,jj,jk) = xcof * e1e2t(ji,jj) * ww(ji,jj,jk) * ( rhop(ji,jj,jk) + rhop(ji,jj,jk-1) ) * tmask_i(ji,jj) + END_3D - zcof = 0.5_wp / rho0 ! Density flux at u and v-points - DO_3D( 1, 0, 1, 0, 1, jpkm1 ) - zkx(ji,jj,jk) = zcof * e2u(ji,jj) * e3u(ji,jj,jk,Kmm) & + DO_3D( 0, 0, 0, 0, 1, jpkm1 ) + zkx(ji,jj,jk) = xcof * e2u(ji,jj) * e3u(ji,jj,jk,Kmm) & & * uu(ji,jj,jk,Kmm) * ( rhop(ji,jj,jk) + rhop(ji+1,jj,jk) ) - zky(ji,jj,jk) = zcof * e1v(ji,jj) * e3v(ji,jj,jk,Kmm) & + zky(ji,jj,jk) = xcof * e1v(ji,jj) * e3v(ji,jj,jk,Kmm) & & * vv(ji,jj,jk,Kmm) * ( rhop(ji,jj,jk) + rhop(ji,jj+1,jk) ) END_3D @@ -227,10 +230,9 @@ CONTAINS ! I.2 Basin averaged kinetic energy trend ! ---------------------------------------- peke = 0._wp - DO jk = 1, jpkm1 - peke = peke + SUM( zkepe(:,:,jk) * gdept(:,:,jk,Kmm) * e1e2t(:,:) & - & * e3t(:,:,jk,Kmm) ) - END DO + DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! loop from top to bottom + peke = peke + zkepe(ji,jj,jk) * gdept(ji,jj,jk,Kmm) * e1e2t(ji,jj) * e3t(ji,jj,jk,Kmm) + END_3D peke = grav * peke ! I.3 Sums over the global domain @@ -509,11 +511,14 @@ CONTAINS WRITE(numout,*) '~~~~~~~~~~~~~' ENDIF + xcof = 0.5_wp / rho0 + ! Total volume at t-points: tvolt = 0._wp - DO jk = 1, jpkm1 - tvolt = tvolt + SUM( e1e2t(:,:) * e3t(:,:,jk,Kmm) * tmask(:,:,jk) * tmask_i(:,:) ) - END DO + + DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! Density flux divergence at t-point + tvolt = tvolt + e1e2t(ji,jj) * e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk) * tmask_i(ji,jj) + END_3D CALL mpp_sum( 'trdglo', tvolt ) ! sum over the global domain IF(lwp) WRITE(numout,*) ' total ocean volume at T-point tvolt = ',tvolt diff --git a/src/OCE/TRD/trdken.F90 b/src/OCE/TRD/trdken.F90 index 1920e3e72..d397a4d15 100644 --- a/src/OCE/TRD/trdken.F90 +++ b/src/OCE/TRD/trdken.F90 @@ -30,6 +30,10 @@ MODULE trdken IMPLICIT NONE PRIVATE + !! * Substitutions +# include "do_loop_substitute.h90" +# include "domzgr_substitute.h90" + PUBLIC trd_ken ! called by trddyn module PUBLIC trd_ken_init ! called by trdini module @@ -38,9 +42,6 @@ MODULE trdken REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: bu, bv ! volume of u- and v-boxes REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: r1_bt ! inverse of t-box volume - !! * Substitutions -# include "do_loop_substitute.h90" -# include "domzgr_substitute.h90" !!---------------------------------------------------------------------- !! NEMO/OCE 4.0 , NEMO Consortium (2018) !! $Id: trdken.F90 15104 2021-07-07 14:36:00Z clem $ @@ -52,7 +53,7 @@ CONTAINS !!--------------------------------------------------------------------- !! *** FUNCTION trd_ken_alloc *** !!--------------------------------------------------------------------- - ALLOCATE( bu(jpi,jpj,jpk) , bv(jpi,jpj,jpk) , r1_bt(jpi,jpj,jpk) , STAT= trd_ken_alloc ) + ALLOCATE( bu(A2D(0),jpk) , bv(A2D(0),jpk) , r1_bt(A2D(0),jpk) , STAT= trd_ken_alloc ) ! CALL mpp_sum ( 'trdken', trd_ken_alloc ) IF( trd_ken_alloc /= 0 ) CALL ctl_stop( 'STOP', 'trd_ken_alloc: failed to allocate arrays' ) @@ -77,31 +78,32 @@ CONTAINS ! ! !!---------------------------------------------------------------------- - REAL(wp), DIMENSION(:,:,:), INTENT(inout) :: putrd, pvtrd ! U and V masked trends - INTEGER , INTENT(in ) :: ktrd ! trend index - INTEGER , INTENT(in ) :: kt ! time step - INTEGER , INTENT(in ) :: Kmm ! time level index + REAL(wp), DIMENSION(:,:,:) , INTENT(inout) :: putrd, pvtrd ! U and V masked trends + INTEGER , INTENT(in ) :: ktrd ! trend index + INTEGER , INTENT(in ) :: kt ! time step + INTEGER , INTENT(in ) :: Kmm ! time level index ! INTEGER :: ji, jj, jk ! dummy loop indices INTEGER :: ikbu , ikbv ! local integers INTEGER :: ikbum1, ikbvm1 ! - - - REAL(wp), DIMENSION(:,:), ALLOCATABLE :: z2dx, z2dy, zke2d ! 2D workspace - REAL(wp), DIMENSION(jpi,jpj,jpk) :: zke ! 3D workspace + REAL(wp), DIMENSION(:,:), ALLOCATABLE :: zke2d ! 2D workspace + REAL(wp) :: z2dx, z2dy, z2dxm1, z2dym1 ! 2D workspace + REAL(wp), DIMENSION(A2D(0),jpk) :: zke ! 3D workspace !!---------------------------------------------------------------------- ! CALL lbc_lnk( 'trdken', putrd, 'U', -1.0_wp , pvtrd, 'V', -1.0_wp ) ! lateral boundary conditions ! nkstp = kt DO jk = 1, jpkm1 - bu (:,:,jk) = e1e2u(:,:) * e3u(:,:,jk,Kmm) - bv (:,:,jk) = e1e2v(:,:) * e3v(:,:,jk,Kmm) - r1_bt(:,:,jk) = r1_e1e2t(:,:) / e3t(:,:,jk,Kmm) * tmask(:,:,jk) + DO_2D( 0, 0, 0, 0 ) + bu (ji,jj,jk) = e1e2u(ji,jj) * e3u(ji,jj,jk,Kmm) + bv (ji,jj,jk) = e1e2v(ji,jj) * e3v(ji,jj,jk,Kmm) + r1_bt(ji,jj,jk) = r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk) + END_2D END DO ! - zke(:,:,jpk) = 0._wp - zke(1:nn_hls,:, : ) = 0._wp - zke(:,1:nn_hls, : ) = 0._wp - DO_3D( 0, nn_hls, 0, nn_hls, 1, jpkm1 ) + zke(A2D(0),:) = 0._wp + DO_3D( 0, 0, 0, 0, 1, jpkm1 ) zke(ji,jj,jk) = 0.5_wp * rho0 *( uu(ji ,jj,jk,Kmm) * putrd(ji ,jj,jk) * bu(ji ,jj,jk) & & + uu(ji-1,jj,jk,Kmm) * putrd(ji-1,jj,jk) * bu(ji-1,jj,jk) & & + vv(ji,jj ,jk,Kmm) * pvtrd(ji,jj ,jk) * bv(ji,jj ,jk) & @@ -118,37 +120,31 @@ CONTAINS CASE( jpdyn_ldf ) ; CALL iom_put( "ketrd_ldf" , zke ) ! lateral diffusion CASE( jpdyn_zdf ) ; CALL iom_put( "ketrd_zdf" , zke ) ! vertical diffusion ! ! ! wind stress trends - ALLOCATE( z2dx(jpi,jpj) , z2dy(jpi,jpj) , zke2d(jpi,jpj) ) - DO_2D( 1, 0, 1, 0 ) - z2dx(ji,jj) = uu(ji,jj,1,Kmm) * ( utau_b(ji,jj) + utauU(ji,jj) ) * e1e2u(ji,jj) * umask(ji,jj,1) - z2dy(ji,jj) = vv(ji,jj,1,Kmm) * ( vtau_b(ji,jj) + vtauV(ji,jj) ) * e1e2v(ji,jj) * vmask(ji,jj,1) - END_2D - DO_2D( 0, 0, 0, 0 ) - zke2d(ji,jj) = r1_rho0 * 0.5_wp * ( z2dx(ji,jj) + z2dx(ji-1,jj) & - & + z2dy(ji,jj) + z2dy(ji,jj-1) ) * r1_bt(ji,jj,1) - END_2D - CALL iom_put( "ketrd_tau" , zke2d ) ! - DEALLOCATE( z2dx , z2dy , zke2d ) - ! + ALLOCATE( zke2d(A2D(0)) ) ; zke2d(A2D(0)) = 0._wp + DO_2D( 0, 0, 0, 0 ) + z2dx = uu(ji,jj,1,Kmm) * ( utau_b(ji,jj) + utauU(ji,jj) ) * e1e2u(ji,jj) * umask(ji,jj,1) + z2dy = vv(ji,jj,1,Kmm) * ( vtau_b(ji,jj) + vtauV(ji,jj) ) * e1e2v(ji,jj) * vmask(ji,jj,1) + z2dxm1 = uu(ji-1,jj,1,Kmm) * ( utau_b(ji-1,jj) + utauU(ji-1,jj) ) * e1e2u(ji-1,jj) * umask(ji-1,jj,1) + z2dym1 = vv(ji,jj-1,1,Kmm) * ( vtau_b(ji,jj-1) + vtauV(ji,jj-1) ) * e1e2v(ji,jj-1) * vmask(ji,jj-1,1) + zke2d(ji,jj) = r1_rho0 * 0.5_wp * ( z2dx + z2dxm1 + z2dy + z2dym1 ) * r1_bt(ji,jj,1) + END_2D + CALL iom_put( "ketrd_tau" , zke2d ) ! + DEALLOCATE( zke2d ) CASE( jpdyn_bfr ) ; CALL iom_put( "ketrd_bfr" , zke ) ! bottom friction (explicit case) !!gm TO BE DONE properly !!gm only valid if ln_drgimp=F otherwise the bottom stress as to be recomputed at the end of the computation.... ! IF(.NOT. ln_drgimp) THEN -! DO jj = 1, jpj ! -! DO ji = 1, jpi -! ikbu = mbku(ji,jj) ! deepest ocean u- & v-levels -! ikbv = mbkv(ji,jj) -! z2dx(ji,jj) = uu(ji,jj,ikbu,Kmm) * bfrua(ji,jj) * uu(ji,jj,ikbu,Kmm) -! z2dy(ji,jj) = vv(ji,jj,ikbu,Kmm) * bfrva(ji,jj) * vv(ji,jj,ikbv,Kmm) -! END DO -! END DO +! zke2d(A2D(0)) = 0._wp +! DO_2D( 0, 0, 0, 0 ) +! ikbu = mbku(ji,jj) ! deepest ocean u- & v-levels +! ikbv = mbkv(ji,jj) +! z2dx = uu(ji,jj,ikbu,Kmm) * bfrua(ji,jj) * uu(ji,jj,ikbu,Kmm) +! z2dy = vv(ji,jj,ikbu,Kmm) * bfrva(ji,jj) * vv(ji,jj,ikbv,Kmm) +! z2dxm1 = uu(ji-1,jj,ikbu,Kmm) * bfrua(ji-1,jj) * uu(ji-1,jj,ikbu,Kmm) +! z2dym1 = vv(ji,jj-1,ikbu,Kmm) * bfrva(ji,jj-1) * vv(ji,jj-1,ikbv,Kmm) +! zke2d(ji,jj) = 0.5_wp * ( z2dx + z2dxm1 + z2dy + z2dym1 ) * r1_bt(ji,jj,1) +! END_2D ! zke2d(1,:) = 0._wp ; zke2d(:,1) = 0._wp -! DO jj = 2, jpj -! DO ji = 2, jpi -! zke2d(ji,jj) = 0.5_wp * ( z2dx(ji,jj) + z2dx(ji-1,jj) & -! & + z2dy(ji,jj) + z2dy(ji,jj-1) ) * r1_bt(ji,jj, BEURK!!! -! END DO -! END DO ! CALL iom_put( "ketrd_bfr" , zke2d ) ! bottom friction (explicit case) ! ENDIF !!gm end @@ -177,11 +173,11 @@ CONTAINS CASE( jpdyn_ken ) ; ! kinetic energy ! called in dynnxt.F90 before asselin time filter ! with putrd=uu(Krhs) and pvtrd=vv(Krhs) - zke(:,:,:) = 0.5_wp * zke(:,:,:) + zke(A2D(0),:) = 0.5_wp * zke(A2D(0),:) CALL iom_put( "KE", zke ) ! CALL ken_p2k( kt , zke, Kmm ) - CALL iom_put( "ketrd_convP2K", zke ) ! conversion -rau*g*w + CALL iom_put( "ketrd_convP2K", zke ) ! conversion -rau*g*w ! END SELECT ! @@ -205,22 +201,23 @@ CONTAINS INTEGER :: ji, jj, jk ! dummy loop indices INTEGER :: iku, ikv ! local integers REAL(wp) :: zcoef ! local scalars - REAL(wp), DIMENSION(jpi,jpj,jpk) :: zconv ! 3D workspace + REAL(wp), DIMENSION(A2D(0),jpk) :: zconv ! 3D workspace !!---------------------------------------------------------------------- ! ! Local constant initialization zcoef = - rho0 * grav * 0.5_wp ! Surface value (also valid in partial step case) - zconv(:,:,1) = zcoef * ( 2._wp * rhd(:,:,1) ) * ww(:,:,1) * e3w(:,:,1,Kmm) - + DO_2D( 0, 0, 0, 0 ) + zconv(ji,jj,1) = zcoef * ( 2._wp * rhd(ji,jj,1) ) * ww(ji,jj,1) * e3w(ji,jj,1,Kmm) + END_2D ! interior value (2=<jk=<jpkm1) - DO jk = 2, jpk - zconv(:,:,jk) = zcoef * ( rhd(:,:,jk) + rhd(:,:,jk-1) ) * ww(:,:,jk) * e3w(:,:,jk,Kmm) - END DO + DO_3D( 0, 0, 0, 0 , 2, jpk ) + zconv(ji,jj,jk) = zcoef * ( rhd(ji,jj,jk) + rhd(ji,jj,jk-1) ) * ww(ji,jj,jk) * e3w(ji,jj,jk,Kmm) + END_3D ! conv value on T-point - DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1 ) + DO_3D( 0, 0, 0, 0, 1, jpkm1 ) zcoef = 0.5_wp / e3t(ji,jj,jk,Kmm) pconv(ji,jj,jk) = zcoef * ( zconv(ji,jj,jk) + zconv(ji,jj,jk+1) ) * tmask(ji,jj,jk) END_3D diff --git a/src/OCE/TRD/trdmxl.F90 b/src/OCE/TRD/trdmxl.F90 index bce0dd536..8edb0c54b 100644 --- a/src/OCE/TRD/trdmxl.F90 +++ b/src/OCE/TRD/trdmxl.F90 @@ -45,20 +45,7 @@ MODULE trdmxl PUBLIC trd_mxl_init ! routine called by opa.F90 PUBLIC trd_mxl_zint ! routine called by tracers routines - INTEGER :: nkstp ! current time step - -!!gm to be moved from trdmxl_oce -! REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: hml ! ML depth (sum of e3t over nmln-1 levels) [m] -! REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: tml , sml ! now ML averaged T & S -! REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: tmlb_nf, smlb_nf ! not filtered before ML averaged T & S -! -! -! REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: hmlb, hmln ! before, now, and after Mixed Layer depths [m] -! -! REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: tb_mlb, tb_mln ! before (not filtered) tracer averaged over before and now ML -! -! REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: tn_mln ! now tracer averaged over now ML -!!gm end + INTEGER :: nkstp ! current time step CHARACTER (LEN=40) :: clhstnam ! name of the trends NetCDF file INTEGER :: nh_t, nmoymltrd @@ -111,44 +98,41 @@ CONTAINS IF ( kt /= nkstp ) THEN != 1st call at kt time step =! ! !==============================! nkstp = kt - - ! !== reset trend arrays to zero ==! - tmltrd(:,:,:) = 0._wp ; smltrd(:,:,:) = 0._wp - - + DO_3D( 0, 0, 0, 0, 1, jpktrd ) + tmltrd(ji,jj,jk) = 0._wp + smltrd(ji,jj,jk) = 0._wp + END_3D ! - wkx(:,:,:) = 0._wp !== now ML weights for vertical averaging ==! - DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpktrd ) ! initialize wkx with vertical scale factor in mixed-layer + ! !== now ML weights for vertical averaging ==! + DO_3D( 0, 0, 0, 0, 1, jpktrd ) ! initialize wkx with vertical scale factor in mixed-layer IF( jk - kmxln(ji,jj) < 0 ) THEN wkx(ji,jj,jk) = e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk) + ELSE + wkx(ji,jj,jk) = 0._wp ENDIF END_3D + ! hmxl(:,:) = 0._wp ! NOW mixed-layer depth - DO jk = 1, jpktrd - hmxl(:,:) = hmxl(:,:) + wkx(:,:,jk) - END DO - DO jk = 1, jpktrd ! integration weights - wkx(:,:,jk) = wkx(:,:,jk) / MAX( 1.e-20_wp, hmxl(:,:) ) * tmask(:,:,1) - END DO - - + DO_3D( 0, 0, 0, 0, 1, jpktrd ) + hmxl(ji,jj) = hmxl(ji,jj) + wkx(ji,jj,jk) + wkx (ji,jj,jk) = wkx (ji,jj,jk) / MAX( 1.e-20_wp, hmxl(ji,jj) ) * tmask(ji,jj,1) + END_3D ! ! !== Vertically averaged T and S ==! tml(:,:) = 0._wp ; sml(:,:) = 0._wp - DO jk = 1, jpktrd - tml(:,:) = tml(:,:) + wkx(:,:,jk) * ts(:,:,jk,jp_tem,Kmm) - sml(:,:) = sml(:,:) + wkx(:,:,jk) * ts(:,:,jk,jp_sal,Kmm) - END DO + DO_3D( 0, 0, 0, 0, 1, jpktrd ) + tml(ji,jj) = tml(ji,jj) + wkx(ji,jj,jk) * ts(ji,jj,jk,jp_tem,Kmm) + sml(ji,jj) = sml(ji,jj) + wkx(ji,jj,jk) * ts(ji,jj,jk,jp_sal,Kmm) + END_3D ! ENDIF - - ! mean now trends over the now ML - tmltrd(:,:,ktrd) = tmltrd(:,:,ktrd) + ptrdx(:,:,jk) * wkx(:,:,jk) ! temperature - smltrd(:,:,ktrd) = smltrd(:,:,ktrd) + ptrdy(:,:,jk) * wkx(:,:,jk) ! salinity - + DO_2D( 0, 0, 0, 0 ) + tmltrd(ji,jj,ktrd) = tmltrd(ji,jj,ktrd) + ptrdx(ji,jj,jk) * wkx(ji,jj,jk) ! temperature + smltrd(ji,jj,ktrd) = smltrd(ji,jj,ktrd) + ptrdy(ji,jj,jk) * wkx(ji,jj,jk) ! salinity + END_2D !!gm to be put juste before the output ! @@ -249,11 +233,11 @@ CONTAINS !!---------------------------------------------------------------------- INTEGER , INTENT( in ) :: ktrd ! ocean trend index CHARACTER(len=2) , INTENT( in ) :: ctype ! 2D surface/bottom or 3D interior physics - REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( in ) :: pttrdmxl ! temperature trend - REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( in ) :: pstrdmxl ! salinity trend + REAL(wp), DIMENSION(A2D(0),jpk), INTENT( in ) :: pttrdmxl ! temperature trend + REAL(wp), DIMENSION(A2D(0),jpk), INTENT( in ) :: pstrdmxl ! salinity trend ! INTEGER :: ji, jj, jk, isum - REAL(wp), DIMENSION(jpi,jpj) :: zvlmsk + REAL(wp), DIMENSION(A2D(0)) :: zvlmsk !!---------------------------------------------------------------------- ! I. Definition of control surface and associated fields @@ -268,11 +252,13 @@ CONTAINS ! ... Set nmxl(ji,jj) = index of first T point below control surf. or outside mixed-layer - IF( nn_ctls == 0 ) THEN ! * control surface = mixed-layer with density criterion - nmxl(:,:) = nmln(:,:) ! array nmln computed in zdfmxl.F90 - ELSEIF( nn_ctls == 1 ) THEN ! * control surface = read index from file + IF( nn_ctls == 0 ) THEN ! * control surface = mixed-layer with density criterion + DO_2D( 0, 0, 0, 0 ) + nmxl(ji,jj) = nmln(ji,jj) ! array nmln computed in zdfmxl.F90 + END_2D + ELSEIF( nn_ctls == 1 ) THEN ! * control surface = read index from file nmxl(:,:) = nbol(:,:) - ELSEIF( nn_ctls >= 2 ) THEN ! * control surface = model level + ELSEIF( nn_ctls >= 2 ) THEN ! * control surface = model level nn_ctls = MIN( nn_ctls, jpktrd - 1 ) nmxl(:,:) = nn_ctls + 1 ENDIF @@ -335,9 +321,9 @@ CONTAINS REAL(wp) :: zavt, zfn, zfn2 ! ! z(ts)mltot : dT/dt over the anlysis window (including Asselin) ! ! z(ts)mlres : residual = dh/dt entrainment term - REAL(wp), DIMENSION(jpi,jpj ) :: ztmltot , zsmltot , ztmlres , zsmlres , ztmlatf , zsmlatf - REAL(wp), DIMENSION(jpi,jpj ) :: ztmltot2, zsmltot2, ztmlres2, zsmlres2, ztmlatf2, zsmlatf2, ztmltrdm2, zsmltrdm2 - REAL(wp), DIMENSION(jpi,jpj,jpk) :: ztmltrd2, zsmltrd2 ! only needed for mean diagnostics + REAL(wp), DIMENSION(A2D(0) ) :: ztmltot , zsmltot , ztmlres , zsmlres , ztmlatf , zsmlatf + REAL(wp), DIMENSION(A2D(0) ) :: ztmltot2, zsmltot2, ztmlres2, zsmlres2, ztmlatf2, zsmlatf2, ztmltrdm2, zsmltrdm2 + REAL(wp), DIMENSION(A2D(0),jpk) :: ztmltrd2, zsmltrd2 ! only needed for mean diagnostics !!---------------------------------------------------------------------- ! ====================================================================== @@ -446,12 +432,7 @@ CONTAINS itmod = kt - nit000 + 1 MODULO_NTRD : IF( MOD( itmod, nn_trd ) == 0 ) THEN ! nitend MUST be multiple of nn_trd - ! - ztmltot (:,:) = 0.e0 ; zsmltot (:,:) = 0.e0 ! reset arrays to zero - ztmlres (:,:) = 0.e0 ; zsmlres (:,:) = 0.e0 - ztmltot2(:,:) = 0.e0 ; zsmltot2(:,:) = 0.e0 - ztmlres2(:,:) = 0.e0 ; zsmlres2(:,:) = 0.e0 - + ! zfn = REAL( nmoymltrd, wp ) ; zfn2 = zfn * zfn ! III.1 Prepare fields for output ("instantaneous" diagnostics) @@ -469,25 +450,18 @@ CONTAINS ztmlatf(:,:) = tmlatfm(:,:) - tmlatfn(:,:) + tmlatfb(:,:) zsmlatf(:,:) = smlatfm(:,:) - smlatfn(:,:) + smlatfb(:,:) - !-- Lateral boundary conditions - ! ... temperature ... ... salinity ... - CALL lbc_lnk( 'trdmxl', ztmltot , 'T', 1.0_wp, zsmltot , 'T', 1.0_wp, & - & ztmlres , 'T', 1.0_wp, zsmlres , 'T', 1.0_wp, & - & ztmlatf , 'T', 1.0_wp, zsmlatf , 'T', 1.0_wp ) - - ! III.2 Prepare fields for output ("mean" diagnostics) ! ---------------------------------------------------- !-- Update the ML depth time sum (to build the Leap-Frog time mean) - hmxl_sum(:,:) = hmxlbn(:,:) + 2 * ( hmxl_sum(:,:) - hmxl(:,:) ) + hmxl(:,:) + hmxl_sum(:,:) = hmxlbn(:,:) + 2._wp * ( hmxl_sum(:,:) - hmxl(:,:) ) + hmxl(:,:) !-- Compute temperature total trends - tml_sum (:,:) = tmlbn(:,:) + 2 * ( tml_sum(:,:) - tml(:,:) ) + tml(:,:) + tml_sum (:,:) = tmlbn(:,:) + 2._wp * ( tml_sum(:,:) - tml(:,:) ) + tml(:,:) ztmltot2(:,:) = ( tml_sum(:,:) - tml_sumb(:,:) ) / p2dt ! now in degC/s !-- Compute salinity total trends - sml_sum (:,:) = smlbn(:,:) + 2 * ( sml_sum(:,:) - sml(:,:) ) + sml(:,:) + sml_sum (:,:) = smlbn(:,:) + 2._wp * ( sml_sum(:,:) - sml(:,:) ) + sml(:,:) zsmltot2(:,:) = ( sml_sum(:,:) - sml_sumb(:,:) ) / p2dt ! now in psu/s !-- Compute temperature residuals @@ -495,7 +469,7 @@ CONTAINS ztmltrd2(:,:,jl) = tmltrd_csum_ub(:,:,jl) + tmltrd_csum_ln(:,:,jl) END DO - ztmltrdm2(:,:) = 0.e0 + ztmltrdm2(:,:) = 0._wp DO jl = 1, jpltrd ztmltrdm2(:,:) = ztmltrdm2(:,:) + ztmltrd2(:,:,jl) END DO @@ -508,7 +482,7 @@ CONTAINS zsmltrd2(:,:,jl) = smltrd_csum_ub(:,:,jl) + smltrd_csum_ln(:,:,jl) END DO - zsmltrdm2(:,:) = 0. + zsmltrdm2(:,:) = 0._wp DO jl = 1, jpltrd zsmltrdm2(:,:) = zsmltrdm2(:,:) + zsmltrd2(:,:,jl) END DO @@ -519,13 +493,6 @@ CONTAINS !-- Diagnose Asselin trend over the analysis window ztmlatf2(:,:) = ztmltrd2(:,:,jpmxl_atf) - tmltrd_sum(:,:,jpmxl_atf) + tmltrd_atf_sumb(:,:) zsmlatf2(:,:) = zsmltrd2(:,:,jpmxl_atf) - smltrd_sum(:,:,jpmxl_atf) + smltrd_atf_sumb(:,:) - - !-- Lateral boundary conditions - ! ... temperature ... ... salinity ... - CALL lbc_lnk( 'trdmxl', ztmltot2, 'T', 1.0_wp, zsmltot2, 'T', 1.0_wp, & - & ztmlres2, 'T', 1.0_wp, zsmlres2, 'T', 1.0_wp ) - ! - CALL lbc_lnk( 'trdmxl', ztmltrd2(:,:,:), 'T', 1.0_wp, zsmltrd2(:,:,:), 'T', 1.0_wp ) ! / in the NetCDF trends file ! III.3 Time evolution array swap ! ------------------------------- @@ -622,6 +589,8 @@ CONTAINS ! ====================================================================== !-- Write the trends for T/S instantaneous diagnostics + + ! clem => these fields do not exist in field_def IF( ln_trdmxl_instant ) THEN diff --git a/src/OCE/TRD/trdmxl_oce.F90 b/src/OCE/TRD/trdmxl_oce.F90 index c2d4d4be9..b1266b5ce 100644 --- a/src/OCE/TRD/trdmxl_oce.F90 +++ b/src/OCE/TRD/trdmxl_oce.F90 @@ -80,6 +80,8 @@ MODULE trdmxl_oce smltrd_csum_ln, & !: ( idem for salinity ) smltrd_csum_ub !: + !! * Substitutions +# include "do_loop_substitute.h90" !!---------------------------------------------------------------------- !! NEMO/OCE 4.0 , NEMO Consortium (2018) !! $Id: trdmxl_oce.F90 10425 2018-12-19 21:54:16Z smasson $ @@ -101,29 +103,29 @@ CONTAINS ierr(:) = 0 - ALLOCATE( nmxl (jpi,jpj) , nbol (jpi,jpj), & - & wkx (jpi,jpj,jpk), hmxl (jpi,jpj), & - & tml (jpi,jpj) , sml (jpi,jpj), & - & tmlb (jpi,jpj) , smlb (jpi,jpj), & - & tmlbb(jpi,jpj) , smlbb(jpi,jpj), STAT = ierr(1) ) - - ALLOCATE( tmlbn(jpi,jpj) , smlbn(jpi,jpj), & - & tmltrdm(jpi,jpj), smltrdm(jpi,jpj), & - & tml_sum(jpi,jpj), tml_sumb(jpi,jpj),& - & tmltrd_atf_sumb(jpi,jpj) , STAT=ierr(2) ) - - ALLOCATE( sml_sum(jpi,jpj), sml_sumb(jpi,jpj), & - & smltrd_atf_sumb(jpi,jpj), & - & hmxl_sum(jpi,jpj), hmxlbn(jpi,jpj), & - & tmlatfb(jpi,jpj), tmlatfn(jpi,jpj), STAT = ierr(3) ) - - ALLOCATE( smlatfb(jpi,jpj), smlatfn(jpi,jpj), & - & tmlatfm(jpi,jpj), smlatfm(jpi,jpj), & - & tmltrd(jpi,jpj,jpltrd), smltrd(jpi,jpj,jpltrd), STAT=ierr(4)) - - ALLOCATE( tmltrd_sum(jpi,jpj,jpltrd),tmltrd_csum_ln(jpi,jpj,jpltrd), & - & tmltrd_csum_ub(jpi,jpj,jpltrd), smltrd_sum(jpi,jpj,jpltrd), & - & smltrd_csum_ln(jpi,jpj,jpltrd), smltrd_csum_ub(jpi,jpj,jpltrd), STAT=ierr(5) ) + ALLOCATE( nmxl (A2D(0)) , nbol (A2D(0)), & + & wkx (A2D(0),jpk), hmxl (A2D(0)), & + & tml (A2D(0)) , sml (A2D(0)), & + & tmlb (A2D(0)) , smlb (A2D(0)), & + & tmlbb(A2D(0)) , smlbb(A2D(0)), STAT = ierr(1) ) + + ALLOCATE( tmlbn(A2D(0)) , smlbn(A2D(0)), & + & tmltrdm(A2D(0)), smltrdm(A2D(0)), & + & tml_sum(A2D(0)), tml_sumb(A2D(0)),& + & tmltrd_atf_sumb(A2D(0)) , STAT=ierr(2) ) + + ALLOCATE( sml_sum(A2D(0)), sml_sumb(A2D(0)), & + & smltrd_atf_sumb(A2D(0)), & + & hmxl_sum(A2D(0)), hmxlbn(A2D(0)), & + & tmlatfb(A2D(0)), tmlatfn(A2D(0)), STAT = ierr(3) ) + + ALLOCATE( smlatfb(A2D(0)), smlatfn(A2D(0)), & + & tmlatfm(A2D(0)), smlatfm(A2D(0)), & + & tmltrd(A2D(0),jpltrd), smltrd(A2D(0),jpltrd), STAT=ierr(4)) + + ALLOCATE( tmltrd_sum(A2D(0),jpltrd),tmltrd_csum_ln(A2D(0),jpltrd), & + & tmltrd_csum_ub(A2D(0),jpltrd), smltrd_sum(A2D(0),jpltrd), & + & smltrd_csum_ln(A2D(0),jpltrd), smltrd_csum_ub(A2D(0),jpltrd), STAT=ierr(5) ) ! trdmxl_oce_alloc = MAXVAL( ierr ) CALL mpp_sum ( 'trdmxl_oce', trdmxl_oce_alloc ) diff --git a/src/OCE/TRD/trdpen.F90 b/src/OCE/TRD/trdpen.F90 index 0c356f841..1a1fdec15 100644 --- a/src/OCE/TRD/trdpen.F90 +++ b/src/OCE/TRD/trdpen.F90 @@ -35,6 +35,7 @@ MODULE trdpen REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: rab_pe ! partial derivatives of PE anomaly with respect to T and S !! * Substitutions +# include "do_loop_substitute.h90" # include "domzgr_substitute.h90" !!---------------------------------------------------------------------- !! NEMO/OCE 4.0 , NEMO Consortium (2018) @@ -48,7 +49,7 @@ CONTAINS !!--------------------------------------------------------------------- !! *** FUNCTION trd_tra_alloc *** !!--------------------------------------------------------------------- - ALLOCATE( rab_pe(jpi,jpj,jpk,jpts) , STAT= trd_pen_alloc ) + ALLOCATE( rab_pe(A2D(0),jpk,jpts) , STAT= trd_pen_alloc ) ! CALL mpp_sum ( 'trdpen', trd_pen_alloc ) IF( trd_pen_alloc /= 0 ) CALL ctl_stop( 'STOP', 'trd_pen_alloc: failed to allocate arrays' ) @@ -69,37 +70,40 @@ CONTAINS INTEGER , INTENT(in) :: Kmm ! time level index REAL(wp) , INTENT(in) :: pdt ! time step [s] ! - INTEGER :: jk ! dummy loop indices - REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: z2d ! 2D workspace - REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpe ! 3D workspace + INTEGER :: ji, jj, jk ! dummy loop indices + REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: z2d ! 2D workspace + REAL(wp), DIMENSION(A2D(0),jpk) :: zpe ! 3D workspace !!---------------------------------------------------------------------- ! - zpe(:,:,:) = 0._wp + zpe(A2D(0),:) = 0._wp ! IF( kt /= nkstp ) THEN ! full eos: set partial derivatives at the 1st call of kt time step nkstp = kt - CALL eos_pen( ts(:,:,:,:,Kmm), rab_PE, zpe, Kmm ) + CALL eos_pen( ts(A2D(0),:,:,Kmm), rab_pe, zpe, Kmm ) CALL iom_put( "alphaPE", rab_pe(:,:,:,jp_tem) ) CALL iom_put( "betaPE" , rab_pe(:,:,:,jp_sal) ) CALL iom_put( "PEanom" , zpe ) ENDIF ! - zpe(:,:,jpk) = 0._wp - DO jk = 1, jpkm1 - zpe(:,:,jk) = ( - ( rab_n(:,:,jk,jp_tem) + rab_pe(:,:,jk,jp_tem) ) * ptrdx(:,:,jk) & - & + ( rab_n(:,:,jk,jp_sal) + rab_pe(:,:,jk,jp_sal) ) * ptrdy(:,:,jk) ) - END DO + zpe(A2D(0),jpk) = 0._wp + ! + DO_3D( 0, 0, 0, 0, 1, jpkm1 ) + zpe(ji,jj,jk) = ( - ( rab_n(ji,jj,jk,jp_tem) + rab_pe(ji,jj,jk,jp_tem) ) * ptrdx(ji,jj,jk) & + & + ( rab_n(ji,jj,jk,jp_sal) + rab_pe(ji,jj,jk,jp_sal) ) * ptrdy(ji,jj,jk) ) + END_3D SELECT CASE ( ktrd ) CASE ( jptra_xad ) ; CALL iom_put( "petrd_xad", zpe ) ! zonal advection CASE ( jptra_yad ) ; CALL iom_put( "petrd_yad", zpe ) ! merid. advection CASE ( jptra_zad ) ; CALL iom_put( "petrd_zad", zpe ) ! vertical advection IF( ln_linssh ) THEN ! cst volume : adv flux through z=0 surface - ALLOCATE( z2d(jpi,jpj) ) - z2d(:,:) = ww(:,:,1) * ( & - & - ( rab_n(:,:,1,jp_tem) + rab_pe(:,:,1,jp_tem) ) * ts(:,:,1,jp_tem,Kmm) & - & + ( rab_n(:,:,1,jp_sal) + rab_pe(:,:,1,jp_sal) ) * ts(:,:,1,jp_sal,Kmm) & - & ) / e3t(:,:,1,Kmm) + ALLOCATE( z2d(A2D(0)) ) + DO_2D( 0, 0, 0, 0 ) + z2d(ji,jj) = ww(ji,jj,1) * ( & + & - ( rab_n(ji,jj,1,jp_tem) + rab_pe(ji,jj,1,jp_tem) ) * ts(ji,jj,1,jp_tem,Kmm) & + & + ( rab_n(ji,jj,1,jp_sal) + rab_pe(ji,jj,1,jp_sal) ) * ts(ji,jj,1,jp_sal,Kmm) & + & ) / e3t(ji,jj,1,Kmm) + END_2D CALL iom_put( "petrd_sad" , z2d ) DEALLOCATE( z2d ) ENDIF diff --git a/src/OCE/TRD/trdtra.F90 b/src/OCE/TRD/trdtra.F90 index a33778632..6c978ddd1 100644 --- a/src/OCE/TRD/trdtra.F90 +++ b/src/OCE/TRD/trdtra.F90 @@ -53,7 +53,7 @@ CONTAINS !!--------------------------------------------------------------------- !! *** FUNCTION trd_tra_alloc *** !!--------------------------------------------------------------------- - ALLOCATE( trdtx(jpi,jpj,jpk) , trdty(jpi,jpj,jpk) , trdt(jpi,jpj,jpk) , avt_evd(jpi,jpj,jpk), STAT= trd_tra_alloc ) + ALLOCATE( trdtx(A2D(0),jpk), trdty(A2D(0),jpk), trdt(A2D(0),jpk), avt_evd(A2D(0),jpk), STAT= trd_tra_alloc ) ! CALL mpp_sum ( 'trdtra', trd_tra_alloc ) IF( trd_tra_alloc /= 0 ) CALL ctl_stop( 'STOP', 'trd_tra_alloc: failed to allocate arrays' ) @@ -73,24 +73,25 @@ CONTAINS !! - 'TRA' case : regroup T & S trends !! - send the trends to trd_tra_mng (trdtrc) for further processing !!---------------------------------------------------------------------- - INTEGER , INTENT(in) :: kt ! time step - CHARACTER(len=3) , INTENT(in) :: ctype ! tracers trends type 'TRA'/'TRC' - INTEGER , INTENT(in) :: ktra ! tracer index - INTEGER , INTENT(in) :: ktrd ! tracer trend index - INTEGER , INTENT(in) :: Kmm, Krhs ! time level indices - REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in) :: ptrd ! tracer trend or flux - REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in), OPTIONAL :: pu ! now velocity - REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in), OPTIONAL :: ptra ! now tracer variable + INTEGER , INTENT(in) :: kt ! time step + CHARACTER(len=3) , INTENT(in) :: ctype ! tracers trends type 'TRA'/'TRC' + INTEGER , INTENT(in) :: ktra ! tracer index + INTEGER , INTENT(in) :: ktrd ! tracer trend index + INTEGER , INTENT(in) :: Kmm, Krhs ! time level indices + REAL(wp), DIMENSION(:,:,:), INTENT(in) :: ptrd ! tracer trend or flux + REAL(wp), DIMENSION(:,:,:), INTENT(in), OPTIONAL :: pu ! now velocity + REAL(wp), DIMENSION(:,:,:), INTENT(in), OPTIONAL :: ptra ! now tracer variable ! - INTEGER :: jk ! loop indices - INTEGER :: i01 ! 0 or 1 - REAL(wp), DIMENSION(jpi,jpj,jpk) :: ztrds ! 3D workspace + INTEGER :: ji,jj,jk ! loop indices + INTEGER :: i01 ! 0 or 1 + REAL(wp) :: z1e3w + REAL(wp), DIMENSION(A2D(0),jpk) :: ztrds ! 3D workspace REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zwt, zws, ztrdt ! 3D workspace !!---------------------------------------------------------------------- ! IF( .NOT. ALLOCATED( trdtx ) ) THEN ! allocate trdtra arrays IF( trd_tra_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'trd_tra : unable to allocate arrays' ) - avt_evd(:,:,:) = 0._wp + avt_evd(A2D(0),:) = 0._wp ENDIF ! i01 = COUNT( (/ PRESENT(pu) .OR. ( ktrd /= jptra_xad .AND. ktrd /= jptra_yad .AND. ktrd /= jptra_zad ) /) ) @@ -103,13 +104,13 @@ CONTAINS CASE( jptra_yad ) ; CALL trd_tra_adv( ptrd, pu, ptra, 'Y', trdty, Kmm ) CASE( jptra_zad ) ; CALL trd_tra_adv( ptrd, pu, ptra, 'Z', trdt, Kmm ) CASE( jptra_bbc, & ! qsr, bbc: on temperature only, send to trd_tra_mng - & jptra_qsr ) ; trdt(:,:,:) = ptrd(:,:,:) * tmask(:,:,:) - ztrds(:,:,:) = 0._wp + & jptra_qsr ) ; trdt(A2D(0),:) = ptrd(A2D(0),:) * tmask(A2D(0),:) + ztrds(A2D(0),:) = 0._wp CALL trd_tra_mng( trdt, ztrds, ktrd, kt, Kmm ) !!gm Gurvan, verify the jptra_evd trend please ! - CASE( jptra_evd ) ; avt_evd(:,:,:) = ptrd(:,:,:) * tmask(:,:,:) + CASE( jptra_evd ) ; avt_evd(A2D(0),:) = ptrd(A2D(0),:) * tmask(A2D(0),:) CASE DEFAULT ! other trends: masked trends - trdt(:,:,:) = ptrd(:,:,:) * tmask(:,:,:) ! mask & store + trdt(A2D(0),:) = ptrd(A2D(0),:) * tmask(A2D(0),:) ! mask & store END SELECT ! ENDIF @@ -127,44 +128,44 @@ CONTAINS CALL trd_tra_mng( trdt , ztrds, ktrd, kt, Kmm ) CASE( jptra_zdfp ) ! diagnose the "PURE" Kz trend (here: just before the swap) ! ! iso-neutral diffusion case otherwise jptra_zdf is "PURE" - ALLOCATE( zwt(jpi,jpj,jpk), zws(jpi,jpj,jpk), ztrdt(jpi,jpj,jpk) ) + ALLOCATE( zwt(A2D(0),jpk), zws(A2D(0),jpk), ztrdt(A2D(0),jpk) ) ! - zwt(:,:, 1 ) = 0._wp ; zws(:,:, 1 ) = 0._wp ! vertical diffusive fluxes - zwt(:,:,jpk) = 0._wp ; zws(:,:,jpk) = 0._wp - DO jk = 2, jpk - zwt(:,:,jk) = avt(:,:,jk) * ( ts(:,:,jk-1,jp_tem,Krhs) - ts(:,:,jk,jp_tem,Krhs) ) & - & / e3w(:,:,jk,Kmm) * tmask(:,:,jk) - zws(:,:,jk) = avs(:,:,jk) * ( ts(:,:,jk-1,jp_sal,Krhs) - ts(:,:,jk,jp_sal,Krhs) ) & - & / e3w(:,:,jk,Kmm) * tmask(:,:,jk) - END DO + zwt(A2D(0), 1 ) = 0._wp ; zws(A2D(0), 1 ) = 0._wp ! vertical diffusive fluxes + zwt(A2D(0),jpk) = 0._wp ; zws(A2D(0),jpk) = 0._wp + DO_3D( 0, 0, 0, 0, 2, jpk ) + z1e3w = 1._wp / ( e3w(ji,jj,jk,Kmm) * tmask(ji,jj,jk) ) + zwt(ji,jj,jk) = avt(ji,jj,jk) * ( ts(ji,jj,jk-1,jp_tem,Krhs) - ts(ji,jj,jk,jp_tem,Krhs) ) * z1e3w + zws(ji,jj,jk) = avs(ji,jj,jk) * ( ts(ji,jj,jk-1,jp_sal,Krhs) - ts(ji,jj,jk,jp_sal,Krhs) ) * z1e3w + END_3D ! - ztrdt(:,:,jpk) = 0._wp ; ztrds(:,:,jpk) = 0._wp - DO jk = 1, jpkm1 - ztrdt(:,:,jk) = ( zwt(:,:,jk) - zwt(:,:,jk+1) ) / e3t(:,:,jk,Kmm) - ztrds(:,:,jk) = ( zws(:,:,jk) - zws(:,:,jk+1) ) / e3t(:,:,jk,Kmm) - END DO + ztrdt(A2D(0),jpk) = 0._wp ; ztrds(A2D(0),jpk) = 0._wp + DO_3D( 0, 0, 0, 0, 1, jpkm1 ) + z1e3w = 1._wp / e3w(ji,jj,jk,Kmm) + ztrdt(ji,jj,jk) = ( zwt(ji,jj,jk) - zwt(ji,jj,jk+1) ) * z1e3w + ztrds(ji,jj,jk) = ( zws(ji,jj,jk) - zws(ji,jj,jk+1) ) * z1e3w + END_3D CALL trd_tra_mng( ztrdt, ztrds, jptra_zdfp, kt, Kmm ) ! ! ! Also calculate EVD trend at this point. - zwt(:,:,:) = 0._wp ; zws(:,:,:) = 0._wp ! vertical diffusive fluxes - DO jk = 2, jpk - zwt(:,:,jk) = avt_evd(:,:,jk) * ( ts(:,:,jk-1,jp_tem,Krhs) - ts(:,:,jk,jp_tem,Krhs) ) & - & / e3w(:,:,jk,Kmm) * tmask(:,:,jk) - zws(:,:,jk) = avt_evd(:,:,jk) * ( ts(:,:,jk-1,jp_sal,Krhs) - ts(:,:,jk,jp_sal,Krhs) ) & - & / e3w(:,:,jk,Kmm) * tmask(:,:,jk) - END DO + zwt(A2D(0), :) = 0._wp ; zws(A2D(0), :) = 0._wp ! vertical diffusive fluxes + DO_3D( 0, 0, 0, 0, 2, jpk ) + z1e3w = 1._wp / ( e3w(ji,jj,jk,Kmm) * tmask(ji,jj,jk) ) + zwt(ji,jj,jk) = avt_evd(ji,jj,jk) * ( ts(ji,jj,jk-1,jp_tem,Krhs) - ts(ji,jj,jk,jp_tem,Krhs) ) * z1e3w + zws(ji,jj,jk) = avt_evd(ji,jj,jk) * ( ts(ji,jj,jk-1,jp_sal,Krhs) - ts(ji,jj,jk,jp_sal,Krhs) ) * z1e3w + END_3D ! - ztrdt(:,:,jpk) = 0._wp ; ztrds(:,:,jpk) = 0._wp - DO jk = 1, jpkm1 - ztrdt(:,:,jk) = ( zwt(:,:,jk) - zwt(:,:,jk+1) ) / e3t(:,:,jk,Kmm) - ztrds(:,:,jk) = ( zws(:,:,jk) - zws(:,:,jk+1) ) / e3t(:,:,jk,Kmm) - END DO + ztrdt(A2D(0),jpk) = 0._wp ; ztrds(A2D(0),jpk) = 0._wp + DO_3D( 0, 0, 0, 0, 1, jpkm1 ) + z1e3w = 1._wp / e3w(ji,jj,jk,Kmm) + ztrdt(ji,jj,jk) = ( zwt(ji,jj,jk) - zwt(ji,jj,jk+1) ) * z1e3w + ztrds(ji,jj,jk) = ( zws(ji,jj,jk) - zws(ji,jj,jk+1) ) * z1e3w + END_3D CALL trd_tra_mng( ztrdt, ztrds, jptra_evd, kt, Kmm ) ! DEALLOCATE( zwt, zws, ztrdt ) ! CASE DEFAULT ! other trends: mask and send T & S trends to trd_tra_mng - ztrds(:,:,:) = ptrd(:,:,:) * tmask(:,:,:) + ztrds(A2D(0),:) = ptrd(A2D(0),:) * tmask(A2D(0),:) CALL trd_tra_mng( trdt, ztrds, ktrd, kt, Kmm ) END SELECT ENDIF @@ -177,7 +178,7 @@ CONTAINS CASE( jptra_yad ) ; CALL trd_tra_adv( ptrd , pu , ptra, 'Y', ztrds, Kmm ) CASE( jptra_zad ) ; CALL trd_tra_adv( ptrd , pu , ptra, 'Z', ztrds, Kmm ) CASE DEFAULT ! other trends: just masked - ztrds(:,:,:) = ptrd(:,:,:) * tmask(:,:,:) + ztrds(A2D(0),:) = ptrd(A2D(0),:) * tmask(A2D(0),:) END SELECT ! ! send trend to trd_trc CALL trd_trc( ztrds, ktra, ktrd, kt, Kmm ) @@ -199,12 +200,12 @@ CONTAINS !! k-advective trends = -un. di+1[T] = -( dk+1[fi] - tn dk+1[un] ) !! where fi is the incoming advective flux. !!---------------------------------------------------------------------- - REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in ) :: pf ! advective flux in one direction - REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in ) :: pu ! now velocity in one direction - REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in ) :: pt ! now or before tracer - CHARACTER(len=1) , INTENT(in ) :: cdir ! X/Y/Z direction - REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( out) :: ptrd ! advective trend in one direction - INTEGER, INTENT(in) :: Kmm ! time level index + REAL(wp), DIMENSION(:,:,:), INTENT(in ) :: pf ! advective flux in one direction + REAL(wp), DIMENSION(:,:,:), INTENT(in ) :: pu ! now velocity in one direction + REAL(wp), DIMENSION(:,:,:), INTENT(in ) :: pt ! now or before tracer + CHARACTER(len=1) , INTENT(in ) :: cdir ! X/Y/Z direction + REAL(wp), DIMENSION(:,:,:), INTENT( out) :: ptrd ! advective trend in one direction + INTEGER, INTENT(in) :: Kmm ! time level index ! INTEGER :: ji, jj, jk ! dummy loop indices INTEGER :: ii, ij, ik ! index shift as function of the direction @@ -217,9 +218,7 @@ CONTAINS END SELECT ! ! ! set to zero uncomputed values - ptrd(jpi,:,:) = 0._wp ; ptrd(1,:,:) = 0._wp - ptrd(:,jpj,:) = 0._wp ; ptrd(:,1,:) = 0._wp - ptrd(:,:,jpk) = 0._wp + ptrd(:,:,:) = 0._wp ! DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! advective trend ptrd(ji,jj,jk) = - ( pf (ji,jj,jk) - pf (ji-ii,jj-ij,jk-ik) & @@ -248,7 +247,7 @@ CONTAINS ! ! 3D output of tracers trends using IOM interface IF( ln_tra_trd ) CALL trd_tra_iom ( ptrdx, ptrdy, ktrd, kt, Kmm ) - ! ! Integral Constraints Properties for tracers trends !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< + ! ! Integral Constraints Properties for tracers trends IF( ln_glo_trd ) CALL trd_glo( ptrdx, ptrdy, ktrd, 'TRA', kt, Kmm ) ! ! Potential ENergy trends @@ -305,8 +304,8 @@ CONTAINS INTEGER , INTENT(in ) :: kt ! time step INTEGER , INTENT(in ) :: Kmm ! time level index !! - INTEGER :: ji, jj, jk ! dummy loop indices - INTEGER :: ikbu, ikbv ! local integers + INTEGER :: ji, jj + REAL(wp) :: z1e3 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: z2dx, z2dy ! 2D workspace !!---------------------------------------------------------------------- ! @@ -329,9 +328,12 @@ CONTAINS CASE( jptra_zad ) ; CALL iom_put( "ttrd_zad" , ptrdx ) ! z- vertical advection CALL iom_put( "strd_zad" , ptrdy ) IF( ln_linssh ) THEN ! cst volume : adv flux through z=0 surface - ALLOCATE( z2dx(jpi,jpj), z2dy(jpi,jpj) ) - z2dx(:,:) = ww(:,:,1) * ts(:,:,1,jp_tem,Kmm) / e3t(:,:,1,Kmm) - z2dy(:,:) = ww(:,:,1) * ts(:,:,1,jp_sal,Kmm) / e3t(:,:,1,Kmm) + ALLOCATE( z2dx(A2D(0)), z2dy(A2D(0)) ) + DO_2D( 0, 0, 0, 0 ) + z1e3 = ww(ji,jj,1) / e3t(ji,jj,1,Kmm) + z2dx(ji,jj) = ts(ji,jj,1,jp_tem,Kmm) * z1e3 + z2dy(ji,jj) = ts(ji,jj,1,jp_sal,Kmm) * z1e3 + END_2D CALL iom_put( "ttrd_sad", z2dx ) CALL iom_put( "strd_sad", z2dy ) DEALLOCATE( z2dx, z2dy ) diff --git a/src/OCE/TRD/trdvor.F90 b/src/OCE/TRD/trdvor.F90 index 8515d666d..b01062697 100644 --- a/src/OCE/TRD/trdvor.F90 +++ b/src/OCE/TRD/trdvor.F90 @@ -68,9 +68,9 @@ CONTAINS !!---------------------------------------------------------------------------- !! *** ROUTINE trd_vor_alloc *** !!---------------------------------------------------------------------------- - ALLOCATE( vor_avr (jpi,jpj) , vor_avrb(jpi,jpj) , vor_avrbb (jpi,jpj) , & - & vor_avrbn (jpi,jpj) , rotot (jpi,jpj) , vor_avrtot(jpi,jpj) , & - & vor_avrres(jpi,jpj) , vortrd (jpi,jpj,jpltot_vor) , & + ALLOCATE( vor_avr (A2D(0)) , vor_avrb(A2D(0)) , vor_avrbb (A2D(0)) , & + & vor_avrbn (A2D(0)) , rotot (A2D(0)) , vor_avrtot(A2D(0)) , & + & vor_avrres(A2D(0)) , vortrd (A2D(0),jpltot_vor) , & & ndexvor1 (jpi*jpj) , STAT= trd_vor_alloc ) ! CALL mpp_sum ( 'trdvor', trd_vor_alloc ) @@ -375,10 +375,6 @@ CONTAINS zmean = 1._wp / REAL( nmoydpvor, wp ) vor_avrres(:,:) = vor_avrtot(:,:) - rotot(:,:) / zmean - ! Boundary conditions - CALL lbc_lnk( 'trdvor', vor_avrtot, 'F', 1.0_wp , vor_avrres, 'F', 1.0_wp ) - - ! III.3 time evolution array swap ! ------------------------------ vor_avrbb(:,:) = vor_avrb(:,:) diff --git a/src/OCE/lib_fortran_generic.h90 b/src/OCE/lib_fortran_generic.h90 index 5201ad6d0..d09c2cd02 100644 --- a/src/OCE/lib_fortran_generic.h90 +++ b/src/OCE/lib_fortran_generic.h90 @@ -51,8 +51,8 @@ iis = Nis0 ; iie = Nie0 ijs = Njs0 ; ije = Nje0 ELSE - iis = 1 ; iie = jpi - ijs = 1 ; ije = jpj + iis = 1 ; iie = ipi + ijs = 1 ; ije = ipj ENDIF ! ctmp = CMPLX( 0.e0, 0.e0, dp ) ! warning ctmp is cumulated -- GitLab