From 6cc152dfc2914d4fc0e5ae93a1a0ed19260b0ade Mon Sep 17 00:00:00 2001 From: Sebastien Masson <sebastien.masson@locean.ipsl.fr> Date: Tue, 28 Jun 2022 07:41:20 +0000 Subject: [PATCH] Resolve "stress on T-points" --- cfgs/AGRIF_DEMO/EXPREF/file_def_nemo-oce.xml | 4 +- cfgs/AMM12/EXPREF/file_def_nemo-oce.xml | 4 +- cfgs/GYRE_PISCES/EXPREF/file_def_nemo.xml | 4 +- .../EXPREF/file_def_nemo-oce.xml | 4 +- .../EXPREF/file_def_nemo-oce.xml | 4 +- .../EXPREF/file_def_nemo-oce.xml | 4 +- cfgs/SHARED/field_def_nemo-oce.xml | 8 +- cfgs/SHARED/namelist_ref | 29 +- cfgs/SPITZ12/EXPREF/file_def_nemo-oce.xml | 4 +- cfgs/WED025/EXPREF/file_def_nemo-oce.xml | 4 +- sette/BATCH_TEMPLATE/batch-X64_JEANZAY | 1 + src/ABL/ablmod.F90 | 43 +- src/ABL/sbcabl.F90 | 4 +- src/ICE/icectl.F90 | 8 +- src/ICE/icedyn_rhg_eap.F90 | 8 +- src/ICE/icedyn_rhg_evp.F90 | 8 +- src/ICE/icedyn_rhg_vp.F90 | 10 +- src/ICE/icesbc.F90 | 10 +- src/ICE/icethd_do.F90 | 4 +- src/ICE/iceupdate.F90 | 33 +- src/OCE/CRS/crsfld.F90 | 4 +- src/OCE/DIA/diawri.F90 | 20 +- src/OCE/DYN/dynatf.F90 | 4 +- src/OCE/DYN/dynatf_qco.F90 | 4 +- src/OCE/DYN/dynspg_ts.F90 | 8 +- src/OCE/DYN/dynzdf.F90 | 8 +- src/OCE/SBC/cpl_oasis3.F90 | 1 + src/OCE/SBC/sbc_ice.F90 | 4 +- src/OCE/SBC/sbc_oce.F90 | 62 +-- src/OCE/SBC/sbcblk.F90 | 232 +++++------ src/OCE/SBC/sbccpl.F90 | 385 ++++++------------ src/OCE/SBC/sbcflx.F90 | 20 +- src/OCE/SBC/sbcmod.F90 | 108 ++--- src/OCE/TRD/trddyn.F90 | 6 +- src/OCE/TRD/trdglo.F90 | 8 +- src/OCE/TRD/trdken.F90 | 22 +- src/OCE/TRD/trdvor.F90 | 35 +- src/OCE/USR/usrdef_sbc.F90 | 12 +- src/OCE/ZDF/zdfosm.F90 | 4 +- src/OCE/ZDF/zdftke.F90 | 10 +- src/OCE/stp2d.F90 | 8 +- src/SAS/diawri.F90 | 12 +- src/SWE/stpmlf.F90 | 4 +- src/SWE/stprk3.F90 | 12 +- tests/BENCH/MY_SRC/usrdef_sbc.F90 | 2 +- tests/CANAL/EXPREF/file_def_nemo-oce.xml | 4 +- tests/DIA_GPU/EXPREF/file_def_nemo-oce.xml | 4 +- tests/DONUT/EXPREF/file_def_nemo-oce.xml | 4 +- tests/ICE_RHEO/EXPREF/file_def_nemo-oce.xml | 7 +- tests/ICE_RHEO/MY_SRC/usrdef_sbc.F90 | 30 +- tests/SWG/EXPREF/file_def_nemo-oce.xml | 4 +- tests/SWG/MY_SRC/usrdef_sbc.F90 | 13 +- tests/WAD/EXPREF/file_def_nemo-oce.xml | 4 +- 53 files changed, 551 insertions(+), 712 deletions(-) diff --git a/cfgs/AGRIF_DEMO/EXPREF/file_def_nemo-oce.xml b/cfgs/AGRIF_DEMO/EXPREF/file_def_nemo-oce.xml index 0d0368d72..6cb04ff08 100644 --- a/cfgs/AGRIF_DEMO/EXPREF/file_def_nemo-oce.xml +++ b/cfgs/AGRIF_DEMO/EXPREF/file_def_nemo-oce.xml @@ -34,6 +34,8 @@ <field field_ref="qt_oce" name="qt_oce" /> <field field_ref="saltflx" name="sfx" /> <field field_ref="taum" name="taum" /> + <field field_ref="utau" name="tauuo" /> + <field field_ref="vtau" name="tauvo" /> <field field_ref="wspd" name="windsp" /> <field field_ref="precip" name="precip" /> <!-- ice and snow --> @@ -44,7 +46,6 @@ <field field_ref="e3u" /> <field field_ref="ssu" name="uos" /> <field field_ref="uoce" name="uo" operation="instant" freq_op="5d" > @uoce_e3u / @e3u </field> - <field field_ref="utau" name="tauuo" /> <field field_ref="uocetr_eff" name="uocetr_eff" /> <!-- available with diaar5 --> <field field_ref="u_masstr" name="vozomatr" /> @@ -56,7 +57,6 @@ <field field_ref="e3v" /> <field field_ref="ssv" name="vos" /> <field field_ref="voce" name="vo" operation="instant" freq_op="5d" > @voce_e3v / @e3v </field> - <field field_ref="vtau" name="tauvo" /> <field field_ref="vocetr_eff" name="vocetr_eff" /> <!-- available with diaar5 --> <field field_ref="v_masstr" name="vomematr" /> diff --git a/cfgs/AMM12/EXPREF/file_def_nemo-oce.xml b/cfgs/AMM12/EXPREF/file_def_nemo-oce.xml index 776def802..d8c729141 100644 --- a/cfgs/AMM12/EXPREF/file_def_nemo-oce.xml +++ b/cfgs/AMM12/EXPREF/file_def_nemo-oce.xml @@ -100,6 +100,8 @@ <field field_ref="qsr" name="rsntds" /> <field field_ref="qt" name="tohfls" /> <field field_ref="taum" /> + <field field_ref="utau" name="tauuo" /> + <field field_ref="vtau" name="tauvo" /> <field field_ref="mldkz5" /> <field field_ref="mldr10_1" /> </file> @@ -107,13 +109,11 @@ <file id="file5" name_suffix="_grid_U" description="ocean U grid variables" > <field field_ref="uoce" name="uo" /> <field field_ref="ssu" name="uos" /> - <field field_ref="utau" name="tauuo" /> </file> <file id="file6" name_suffix="_grid_V" description="ocean V grid variables" > <field field_ref="voce" name="vo" /> <field field_ref="ssv" name="vos" /> - <field field_ref="vtau" name="tauvo" /> </file> <file id="file7" name_suffix="_grid_W" description="ocean W grid variables" > diff --git a/cfgs/GYRE_PISCES/EXPREF/file_def_nemo.xml b/cfgs/GYRE_PISCES/EXPREF/file_def_nemo.xml index e24e2271a..57144055d 100644 --- a/cfgs/GYRE_PISCES/EXPREF/file_def_nemo.xml +++ b/cfgs/GYRE_PISCES/EXPREF/file_def_nemo.xml @@ -32,16 +32,16 @@ <field field_ref="qt" name="sohefldo" /> <field field_ref="mldr10_1" name="somxl010" /> <field field_ref="mldkz5" name="somixhgt" /> + <field field_ref="utau" name="sozotaux" /> + <field field_ref="vtau" name="sometauy" /> </file> <file id="file2" name_suffix="_grid_U" description="ocean U grid variables" > <field field_ref="uoce" name="vozocrtx" /> - <field field_ref="utau" name="sozotaux" /> </file> <file id="file3" name_suffix="_grid_V" description="ocean V grid variables" > <field field_ref="voce" name="vomecrty" /> - <field field_ref="vtau" name="sometauy" /> </file> <file id="file4" name_suffix="_grid_W" description="ocean W grid variables" > diff --git a/cfgs/ORCA2_ICE_ABL/EXPREF/file_def_nemo-oce.xml b/cfgs/ORCA2_ICE_ABL/EXPREF/file_def_nemo-oce.xml index 6d6df82d2..9ba0c304f 100644 --- a/cfgs/ORCA2_ICE_ABL/EXPREF/file_def_nemo-oce.xml +++ b/cfgs/ORCA2_ICE_ABL/EXPREF/file_def_nemo-oce.xml @@ -30,6 +30,8 @@ <field field_ref="qt_oce" name="qt_oce" /> <field field_ref="saltflx" name="sfx" /> <field field_ref="taum" name="taum" /> + <field field_ref="utau" name="tauuo" /> + <field field_ref="vtau" name="tauvo" /> <field field_ref="wspd" name="windsp" /> <field field_ref="precip" name="precip" /> <!-- ice and snow --> @@ -40,14 +42,12 @@ <field field_ref="e3u" /> <field field_ref="ssu" name="uos" /> <field field_ref="uoce" name="uo" /> - <field field_ref="utau" name="tauuo" /> </file> <file id="file13" name_suffix="_grid_V" description="ocean V grid variables" > <field field_ref="e3v" /> <field field_ref="ssv" name="vos" /> <field field_ref="voce" name="vo" /> - <field field_ref="vtau" name="tauvo" /> </file> <file id="file14" name_suffix="_grid_ABL" description="ABL grid variables" > diff --git a/cfgs/ORCA2_ICE_PISCES/EXPREF/file_def_nemo-oce.xml b/cfgs/ORCA2_ICE_PISCES/EXPREF/file_def_nemo-oce.xml index c043f56fe..04590b55f 100644 --- a/cfgs/ORCA2_ICE_PISCES/EXPREF/file_def_nemo-oce.xml +++ b/cfgs/ORCA2_ICE_PISCES/EXPREF/file_def_nemo-oce.xml @@ -34,6 +34,8 @@ <field field_ref="qt_oce" name="qt_oce" /> <field field_ref="saltflx" name="sfx" /> <field field_ref="taum" name="taum" /> + <field field_ref="utau" name="tauuo" /> + <field field_ref="vtau" name="tauvo" /> <field field_ref="wspd" name="windsp" /> <field field_ref="precip" name="precip" /> <!-- ice and snow --> @@ -44,7 +46,6 @@ <field field_ref="e3u" /> <field field_ref="ssu" name="uos" /> <field field_ref="uoce" name="uo" operation="instant" freq_op="5d" > @uoce_e3u / @e3u </field> - <field field_ref="utau" name="tauuo" /> <field field_ref="uocetr_eff" name="uocetr_eff" /> <!-- available with diaar5 --> <field field_ref="u_masstr" name="vozomatr" /> @@ -56,7 +57,6 @@ <field field_ref="e3v" /> <field field_ref="ssv" name="vos" /> <field field_ref="voce" name="vo" operation="instant" freq_op="5d" > @voce_e3v / @e3v </field> - <field field_ref="vtau" name="tauvo" /> <field field_ref="vocetr_eff" name="vocetr_eff" /> <!-- available with diaar5 --> <field field_ref="v_masstr" name="vomematr" /> diff --git a/cfgs/ORCA2_SAS_ICE/EXPREF/file_def_nemo-oce.xml b/cfgs/ORCA2_SAS_ICE/EXPREF/file_def_nemo-oce.xml index 5c75b8920..2cdea4e23 100644 --- a/cfgs/ORCA2_SAS_ICE/EXPREF/file_def_nemo-oce.xml +++ b/cfgs/ORCA2_SAS_ICE/EXPREF/file_def_nemo-oce.xml @@ -34,6 +34,8 @@ <field field_ref="qt_oce" name="qt_oce" /> <field field_ref="saltflx" name="sfx" /> <field field_ref="taum" name="taum" /> + <field field_ref="utau" name="tauuo" /> + <field field_ref="vtau" name="tauvo" /> <field field_ref="wspd" name="windsp" /> <field field_ref="precip" name="precip" /> <!-- ice and snow --> @@ -44,7 +46,6 @@ <field field_ref="e3u" /> <field field_ref="ssu" name="uos" /> <field field_ref="uoce" name="uo" operation="instant" freq_op="5d" > @uoce_e3u / @e3u </field> - <field field_ref="utau" name="tauuo" /> <field field_ref="uocetr_eff" name="uocetr_eff" /> <!-- available with diaar5 --> <field field_ref="u_masstr" name="vozomatr" /> @@ -56,7 +57,6 @@ <field field_ref="e3v" /> <field field_ref="ssv" name="vos" /> <field field_ref="voce" name="vo" operation="instant" freq_op="5d" > @voce_e3v / @e3v </field> - <field field_ref="vtau" name="tauvo" /> <field field_ref="vocetr_eff" name="vocetr_eff" /> <!-- available with diaar5 --> <field field_ref="v_masstr" name="vomematr" /> diff --git a/cfgs/SHARED/field_def_nemo-oce.xml b/cfgs/SHARED/field_def_nemo-oce.xml index ea75b3ec7..e5bff5cc7 100644 --- a/cfgs/SHARED/field_def_nemo-oce.xml +++ b/cfgs/SHARED/field_def_nemo-oce.xml @@ -432,6 +432,8 @@ that are available in the tidal-forcing implementation (see <field id="erp" long_name="Surface Water Flux: Damping" standard_name="water_flux_out_of_sea_water_due_to_newtonian_relaxation" unit="kg/m2/s" /> <field id="taum" long_name="wind stress module" standard_name="magnitude_of_surface_downward_stress" unit="N/m2" /> <field id="wspd" long_name="wind speed module" standard_name="wind_speed" unit="m/s" /> + <field id="utau" long_name="Wind Stress along i-axis" standard_name="surface_downward_x_stress" unit="N/m2" /> + <field id="vtau" long_name="Wind Stress along j-axis" standard_name="surface_downward_y_stress" unit="N/m2" /> <!-- * variable relative to atmospheric pressure forcing : available with ln_apr_dyn --> <field id="ssh_ib" long_name="Inverse barometer sea surface height" standard_name="sea_surface_height_correction_due_to_air_pressure_at_low_frequency" unit="m" /> @@ -583,7 +585,6 @@ that are available in the tidal-forcing implementation (see <field id="e2u" long_name="U-cell width in meridional direction" standard_name="cell_width" unit="m" /> <field id="e3u" long_name="U-cell thickness" standard_name="cell_thickness" unit="m" grid_ref="grid_U_3D" /> <field id="e3u_0" long_name="Initial U-cell thickness" standard_name="ref_cell_thickness" unit="m" grid_ref="grid_U_3D"/> - <field id="utau" long_name="Wind Stress along i-axis" standard_name="surface_downward_x_stress" unit="N/m2" /> <field id="uoce" long_name="ocean current along i-axis" standard_name="sea_water_x_velocity" unit="m/s" grid_ref="grid_U_3D" /> <field id="uoce_e3u" long_name="ocean current along i-axis (thickness weighted)" unit="m/s" grid_ref="grid_U_3D" > uoce * e3u </field> <field id="uoce_e3u_vsum" long_name="ocean current along i-axis * e3u summed on the vertical" field_ref="uoce_e3u" unit="m3/s" grid_ref="grid_U_vsum"/> @@ -651,7 +652,6 @@ that are available in the tidal-forcing implementation (see <field id="e3v" long_name="V-cell thickness" standard_name="cell_thickness" unit="m" grid_ref="grid_V_3D" /> <field id="e3v_0" long_name="Initial V-cell thickness" standard_name="ref_cell_thickness" unit="m" grid_ref="grid_V_3D" /> <field id="hv" long_name="water column height at V point" standard_name="water_column_height_V" unit="m" /> - <field id="vtau" long_name="Wind Stress along j-axis" standard_name="surface_downward_y_stress" unit="N/m2" /> <field id="voce" long_name="ocean current along j-axis" standard_name="sea_water_y_velocity" unit="m/s" grid_ref="grid_V_3D" /> <field id="voce_e3v" long_name="ocean current along j-axis (thickness weighted)" unit="m/s" grid_ref="grid_V_3D" > voce * e3v </field> <field id="ssv" long_name="ocean surface current along j-axis" unit="m/s" /> @@ -1189,6 +1189,8 @@ that are available in the tidal-forcing implementation (see <field field_ref="qsr" name="rsntds" long_name="surface_net_downward_shortwave_flux" /> <field field_ref="qt" name="tohfls" long_name="surface_net_downward_total_heat_flux" /> <field field_ref="taum" /> + <field field_ref="utau" name="tauuo" long_name="surface_downward_x_stress" /> + <field field_ref="vtau" name="tauvo" long_name="surface_downward_y_stress" /> <field field_ref="20d" /> <field field_ref="mldkz5" /> <field field_ref="mldr10_1" /> @@ -1203,12 +1205,10 @@ that are available in the tidal-forcing implementation (see <field_group id="groupU" > <field field_ref="uoce" name="uo" long_name="sea_water_x_velocity" /> - <field field_ref="utau" name="tauuo" long_name="surface_downward_x_stress" /> </field_group> <field_group id="groupV" > <field field_ref="voce" name="vo" long_name="sea_water_y_velocity" /> - <field field_ref="vtau" name="tauvo" long_name="surface_downward_y_stress" /> </field_group> <field_group id="groupW" > diff --git a/cfgs/SHARED/namelist_ref b/cfgs/SHARED/namelist_ref index ae9d42ee7..dc77ca9a5 100644 --- a/cfgs/SHARED/namelist_ref +++ b/cfgs/SHARED/namelist_ref @@ -369,7 +369,7 @@ !*** receive *** sn_rcv_w10m = 'none' , 'no' , '' , '' , '' sn_rcv_taumod = 'coupled' , 'no' , '' , '' , '' - sn_rcv_tau = 'oce only' , 'no' , 'cartesian' , 'eastward-northward' , 'U,V' + sn_rcv_tau = 'oce only' , 'no' , 'cartesian' , 'eastward-northward' , '' sn_rcv_dqnsdt = 'coupled' , 'no' , '' , '' , '' sn_rcv_qsr = 'oce and ice' , 'no' , '' , '' , '' sn_rcv_qns = 'oce and ice' , 'no' , '' , '' , '' @@ -380,21 +380,22 @@ sn_rcv_iceflx = 'none' , 'no' , '' , '' , '' sn_rcv_mslp = 'none' , 'no' , '' , '' , '' sn_rcv_ts_ice = 'none' , 'no' , '' , '' , '' + sn_rcv_qtrice = 'none' , 'no' , '' , '' , '' sn_rcv_isf = 'none' , 'no' , '' , '' , '' sn_rcv_icb = 'none' , 'no' , '' , '' , '' - sn_rcv_hsig = 'none' , 'no' , '' , '' , 'T' - sn_rcv_phioc = 'none' , 'no' , '' , '' , 'T' - sn_rcv_sdrfx = 'none' , 'no' , '' , '' , 'T' - sn_rcv_sdrfy = 'none' , 'no' , '' , '' , 'T' - sn_rcv_wper = 'none' , 'no' , '' , '' , 'T' - sn_rcv_wnum = 'none' , 'no' , '' , '' , 'T' - sn_rcv_wstrf = 'none' , 'no' , '' , '' , 'T' - sn_rcv_wdrag = 'none' , 'no' , '' , '' , 'T' - sn_rcv_charn = 'none' , 'no' , '' , '' , 'T' - sn_rcv_taw = 'none' , 'no' , '' , '' , 'U,V' - sn_rcv_bhd = 'none' , 'no' , '' , '' , 'T' - sn_rcv_tusd = 'none' , 'no' , '' , '' , 'T' - sn_rcv_tvsd = 'none' , 'no' , '' , '' , 'T' + sn_rcv_hsig = 'none' , 'no' , '' , '' , '' + sn_rcv_phioc = 'none' , 'no' , '' , '' , '' + sn_rcv_sdrfx = 'none' , 'no' , '' , '' , '' + sn_rcv_sdrfy = 'none' , 'no' , '' , '' , '' + sn_rcv_wper = 'none' , 'no' , '' , '' , '' + sn_rcv_wnum = 'none' , 'no' , '' , '' , '' + sn_rcv_wstrf = 'none' , 'no' , '' , '' , '' + sn_rcv_wdrag = 'none' , 'no' , '' , '' , '' + sn_rcv_charn = 'none' , 'no' , '' , '' , '' + sn_rcv_taw = 'none' , 'no' , '' , '' , '' + sn_rcv_bhd = 'none' , 'no' , '' , '' , '' + sn_rcv_tusd = 'none' , 'no' , '' , '' , '' + sn_rcv_tvsd = 'none' , 'no' , '' , '' , '' / !----------------------------------------------------------------------- &namsbc_sas ! Stand-Alone Surface module: ocean data (SAS_SRC only) diff --git a/cfgs/SPITZ12/EXPREF/file_def_nemo-oce.xml b/cfgs/SPITZ12/EXPREF/file_def_nemo-oce.xml index aa6deff58..ee0982291 100644 --- a/cfgs/SPITZ12/EXPREF/file_def_nemo-oce.xml +++ b/cfgs/SPITZ12/EXPREF/file_def_nemo-oce.xml @@ -34,6 +34,8 @@ <field field_ref="qt_oce" name="qt_oce" /> <field field_ref="saltflx" name="sfx" /> <field field_ref="taum" name="taum" /> + <field field_ref="utau" name="tauuo" /> + <field field_ref="vtau" name="tauvo" /> <field field_ref="wspd" name="windsp" /> <field field_ref="precip" name="precip" /> <!-- ice and snow --> @@ -44,7 +46,6 @@ <field field_ref="e3u" /> <field field_ref="ssu" name="uos" /> <field field_ref="uoce" name="uo" operation="instant" freq_op="5d" > @uoce_e3u / @e3u </field> - <field field_ref="utau" name="tauuo" /> <field field_ref="uocetr_eff" name="uocetr_eff" /> <!-- available with diaar5 --> <field field_ref="u_masstr" name="vozomatr" /> @@ -56,7 +57,6 @@ <field field_ref="e3v" /> <field field_ref="ssv" name="vos" /> <field field_ref="voce" name="vo" operation="instant" freq_op="5d" > @voce_e3v / @e3v </field> - <field field_ref="vtau" name="tauvo" /> <field field_ref="vocetr_eff" name="vocetr_eff" /> <!-- available with diaar5 --> <field field_ref="v_masstr" name="vomematr" /> diff --git a/cfgs/WED025/EXPREF/file_def_nemo-oce.xml b/cfgs/WED025/EXPREF/file_def_nemo-oce.xml index 8d9d89e5f..d08f13e1f 100644 --- a/cfgs/WED025/EXPREF/file_def_nemo-oce.xml +++ b/cfgs/WED025/EXPREF/file_def_nemo-oce.xml @@ -31,6 +31,8 @@ <field field_ref="qt_oce" name="qt_oce" /> <field field_ref="saltflx" name="sfx" /> <field field_ref="taum" name="taum" /> + <field field_ref="utau" name="tauuo" /> + <field field_ref="vtau" name="tauvo" /> <field field_ref="wspd" name="windsp" /> <field field_ref="precip" name="precip" /> <field field_ref="snowpre" /> @@ -67,14 +69,12 @@ <field field_ref="e3u" /> <field field_ref="ssu" name="uos" /> <field field_ref="uoce" name="uo" operation="instant" freq_op="5d" > @uoce_e3u / @e3u </field> - <field field_ref="utau" name="tauuo" /> </file> <file id="file13" name_suffix="_grid_V" description="ocean V grid variables" > <field field_ref="e3v" /> <field field_ref="ssv" name="vos" /> <field field_ref="voce" name="vo" operation="instant" freq_op="5d" > @voce_e3v / @e3v </field> - <field field_ref="vtau" name="tauvo" /> </file> <file id="file14" name_suffix="_grid_W" description="ocean W grid variables" > diff --git a/sette/BATCH_TEMPLATE/batch-X64_JEANZAY b/sette/BATCH_TEMPLATE/batch-X64_JEANZAY index e7af8cea4..cafaff218 100644 --- a/sette/BATCH_TEMPLATE/batch-X64_JEANZAY +++ b/sette/BATCH_TEMPLATE/batch-X64_JEANZAY @@ -2,6 +2,7 @@ #SBATCH -A GROUP_IDRIS@cpu #SBATCH --job-name=SETTE_JOB # nom du job #SBATCH --partition=cpu_p1 # Nom de la partition d'exécution +#SBATCH --qos=qos_cpu-dev # quality of service #SBATCH --ntasks=NPROCS # Nombre total de processus MPI #SBATCH --ntasks-per-node=40 # Nombre de processus MPI par noeud # /!\ Attention, la ligne suivante est trompeuse mais dans le vocabulaire diff --git a/src/ABL/ablmod.F90 b/src/ABL/ablmod.F90 index b5b66b8a6..ffaddbbcd 100644 --- a/src/ABL/ablmod.F90 +++ b/src/ABL/ablmod.F90 @@ -98,8 +98,8 @@ CONTAINS REAL(wp) , INTENT(in ), DIMENSION(:,: ) :: pevp_ice ! Ce x Du over ice (T-point) REAL(wp) , INTENT(in ), DIMENSION(:,: ) :: pwndm_ice ! ||uwnd - uice|| REAL(wp) , INTENT(in ), DIMENSION(:,: ) :: pfrac_oce ! ocean fraction - REAL(wp) , INTENT( out), DIMENSION(:,: ) :: ptaui_ice ! ice-surface taux stress (U-point) - REAL(wp) , INTENT( out), DIMENSION(:,: ) :: ptauj_ice ! ice-surface tauy stress (V-point) + REAL(wp) , INTENT( out), DIMENSION(:,: ) :: ptaui_ice ! ice-surface taux stress (T-point) + REAL(wp) , INTENT( out), DIMENSION(:,: ) :: ptauj_ice ! ice-surface tauy stress (T-point) #endif ! REAL(wp), DIMENSION(1:jpi,1:jpj ) :: zwnd_i, zwnd_j @@ -617,25 +617,17 @@ CONTAINS zwnd_i(ji,jj) = zztmp * zwnd_i(ji,jj) zwnd_j(ji,jj) = zztmp * zwnd_j(ji,jj) END_2D - ! ... utau, vtau at U- and V_points, resp. - ! Note the use of 0.5*(2-umask) in order to unmask the stress along coastlines - ! Note the use of MAX(tmask(i,j),tmask(i+1,j) is to mask tau over ice shelves - DO_2D( 0, 0, 0, 0 ) - zcff = 0.5_wp * ( 2._wp - msk_abl(ji,jj)*msk_abl(ji+1,jj) ) - zztmp = MAX(msk_abl(ji,jj),msk_abl(ji+1,jj)) - ptaui(ji,jj) = zcff * zztmp * ( zwnd_i(ji,jj) + zwnd_i(ji+1,jj ) ) - zcff = 0.5_wp * ( 2._wp - msk_abl(ji,jj)*msk_abl(ji,jj+1) ) - zztmp = MAX(msk_abl(ji,jj),msk_abl(ji,jj+1)) - ptauj(ji,jj) = zcff * zztmp * ( zwnd_j(ji,jj) + zwnd_j(ji ,jj+1) ) + ! ... utau, vtau at T-points + DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) + ptaui(ji,jj) = zwnd_i(ji,jj) * msk_abl(ji,jj) !!clem tau: check + ptauj(ji,jj) = zwnd_j(ji,jj) * msk_abl(ji,jj) !!clem tau: check END_2D - ! - CALL lbc_lnk( 'ablmod', ptaui(:,:), 'U', -1.0_wp, ptauj(:,:), 'V', -1.0_wp ) CALL iom_put( "taum_oce", ptaum ) IF(sn_cfctl%l_prtctl) THEN - CALL prt_ctl( tab2d_1=ptaui , clinfo1=' abl_stp: utau : ', mask1=umask, & - & tab2d_2=ptauj , clinfo2=' vtau : ', mask2=vmask ) + CALL prt_ctl( tab2d_1=ptaui , clinfo1=' abl_stp: utau : ', mask1=tmask, & + & tab2d_2=ptauj , clinfo2=' vtau : ', mask2=tmask ) CALL prt_ctl( tab2d_1=pwndm , clinfo1=' abl_stp: wndm : ' ) ENDIF @@ -657,23 +649,14 @@ CONTAINS ! ------------------------------------------------------------ ! ! Wind stress relative to the moving ice ( U10m - U_ice ) ! ! ------------------------------------------------------------ ! - DO_2D( 0, 0, 0, 0 ) - - zztmp1 = 0.5_wp * ( u_abl(ji+1,jj ,2,nt_a) + u_abl(ji,jj,2,nt_a) ) - zztmp2 = 0.5_wp * ( v_abl(ji ,jj+1,2,nt_a) + v_abl(ji,jj,2,nt_a) ) - - ptaui_ice(ji,jj) = 0.5_wp * ( rhoa(ji+1,jj) * pCd_du_ice(ji+1,jj) & - & + rhoa(ji ,jj) * pCd_du_ice(ji ,jj) ) & - & * ( zztmp1 - pssu_ice(ji,jj) * rn_vfac ) - ptauj_ice(ji,jj) = 0.5_wp * ( rhoa(ji,jj+1) * pCd_du_ice(ji,jj+1) & - & + rhoa(ji,jj ) * pCd_du_ice(ji,jj ) ) & - & * ( zztmp2 - pssv_ice(ji,jj) * rn_vfac ) + DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) + ptaui_ice(ji,jj) = rhoa(ji,jj) * pCd_du_ice(ji,jj) * ( u_abl(ji,jj,2,nt_a) - pssu_ice(ji,jj) * rn_vfac ) + ptauj_ice(ji,jj) = rhoa(ji,jj) * pCd_du_ice(ji,jj) * ( v_abl(ji,jj,2,nt_a) - pssv_ice(ji,jj) * rn_vfac ) END_2D - CALL lbc_lnk( 'ablmod', ptaui_ice, 'U', -1.0_wp, ptauj_ice,'V', -1.0_wp ) ! IF(sn_cfctl%l_prtctl) THEN - CALL prt_ctl( tab2d_1=ptaui_ice , clinfo1=' abl_stp: utau_ice : ', mask1=umask, & - & tab2d_2=ptauj_ice , clinfo2=' vtau_ice : ', mask2=vmask ) + CALL prt_ctl( tab2d_1=ptaui_ice , clinfo1=' abl_stp: utau_ice : ', mask1=tmask, & + & tab2d_2=ptauj_ice , clinfo2=' vtau_ice : ', mask2=tmask ) END IF #endif ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< diff --git a/src/ABL/sbcabl.F90 b/src/ABL/sbcabl.F90 index f087614e3..c2c1d13c1 100644 --- a/src/ABL/sbcabl.F90 +++ b/src/ABL/sbcabl.F90 @@ -308,8 +308,8 @@ CONTAINS !! - Perform 1 time-step of the ABL model !! - Finalize flux computation in blk_oce_2 !! - !! ** Outputs : - utau : i-component of the stress at U-point (N/m2) - !! - vtau : j-component of the stress at V-point (N/m2) + !! ** Outputs : - utau : i-component of the stress at T-point (N/m2) + !! - vtau : j-component of the stress at T-point (N/m2) !! - taum : Wind stress module at T-point (N/m2) !! - wndm : Wind speed module at T-point (m/s) !! - qsr : Solar heat flux over the ocean (W/m2) diff --git a/src/ICE/icectl.F90 b/src/ICE/icectl.F90 index ad20e6734..9593be4f8 100644 --- a/src/ICE/icectl.F90 +++ b/src/ICE/icectl.F90 @@ -733,10 +733,10 @@ CONTAINS CALL prt_ctl_info(' ') CALL prt_ctl_info(' - Stresses : ') CALL prt_ctl_info(' ~~~~~~~~~~ ') - CALL prt_ctl(tab2d_1=utau , clinfo1= ' utau : ', mask1 = umask, & - & tab2d_2=vtau , clinfo2= ' vtau : ', mask2 = vmask) - CALL prt_ctl(tab2d_1=utau_ice , clinfo1= ' utau_ice : ', mask1 = umask, & - & tab2d_2=vtau_ice , clinfo2= ' vtau_ice : ', mask2 = vmask) + CALL prt_ctl(tab2d_1=utau , clinfo1= ' utau : ', mask1 = tmask, & + & tab2d_2=vtau , clinfo2= ' vtau : ', mask2 = tmask) + CALL prt_ctl(tab2d_1=utau_ice , clinfo1= ' utau_ice : ', mask1 = tmask, & + & tab2d_2=vtau_ice , clinfo2= ' vtau_ice : ', mask2 = tmask) END SUBROUTINE ice_prt3D diff --git a/src/ICE/icedyn_rhg_eap.F90 b/src/ICE/icedyn_rhg_eap.F90 index 2dcfc5759..797ceba71 100644 --- a/src/ICE/icedyn_rhg_eap.F90 +++ b/src/ICE/icedyn_rhg_eap.F90 @@ -311,8 +311,12 @@ CONTAINS zmV_t(ji,jj) = zmassV * z1_dtevp ! Drag ice-atm. - ztaux_ai(ji,jj) = zaU(ji,jj) * utau_ice(ji,jj) - ztauy_ai(ji,jj) = zaV(ji,jj) * vtau_ice(ji,jj) + ! Note the use of 0.5*(2-umask) in order to unmask the stress along coastlines + ! and the use of MAX(tmask(i,j),tmask(i+1,j) is to mask tau over ice shelves + ztaux_ai(ji,jj) = zaU(ji,jj) * 0.5_wp * ( utau_ice(ji,jj) + utau_ice(ji+1,jj) ) * & + & ( 2. - umask(ji,jj,1) ) * MAX( tmask(ji,jj,1), tmask(ji+1,jj,1) ) + ztauy_ai(ji,jj) = zaV(ji,jj) * 0.5_wp * ( vtau_ice(ji,jj) + vtau_ice(ji,jj+1) ) * & + & ( 2. - vmask(ji,jj,1) ) * MAX( tmask(ji,jj,1), tmask(ji,jj+1,1) ) ! Surface pressure gradient (- m*g*GRAD(ssh)) at U-V points zspgU(ji,jj) = - zmassU * grav * ( zsshdyn(ji+1,jj) - zsshdyn(ji,jj) ) * r1_e1u(ji,jj) diff --git a/src/ICE/icedyn_rhg_evp.F90 b/src/ICE/icedyn_rhg_evp.F90 index 52f15aae8..a260797a2 100644 --- a/src/ICE/icedyn_rhg_evp.F90 +++ b/src/ICE/icedyn_rhg_evp.F90 @@ -292,8 +292,12 @@ CONTAINS zmV_t(ji,jj) = zmassV * z1_dtevp ! Drag ice-atm. - ztaux_ai(ji,jj) = zaU(ji,jj) * utau_ice(ji,jj) - ztauy_ai(ji,jj) = zaV(ji,jj) * vtau_ice(ji,jj) + ! Note the use of 0.5*(2-umask) in order to unmask the stress along coastlines + ! and the use of MAX(tmask(i,j),tmask(i+1,j) is to mask tau over ice shelves + ztaux_ai(ji,jj) = zaU(ji,jj) * 0.5_wp * ( utau_ice(ji,jj) + utau_ice(ji+1,jj) ) * & + & ( 2. - umask(ji,jj,1) ) * MAX( tmask(ji,jj,1), tmask(ji+1,jj,1) ) + ztauy_ai(ji,jj) = zaV(ji,jj) * 0.5_wp * ( vtau_ice(ji,jj) + vtau_ice(ji,jj+1) ) * & + & ( 2. - vmask(ji,jj,1) ) * MAX( tmask(ji,jj,1), tmask(ji,jj+1,1) ) ! Surface pressure gradient (- m*g*GRAD(ssh)) at U-V points zspgU(ji,jj) = - zmassU * grav * ( zsshdyn(ji+1,jj) - zsshdyn(ji,jj) ) * r1_e1u(ji,jj) diff --git a/src/ICE/icedyn_rhg_vp.F90 b/src/ICE/icedyn_rhg_vp.F90 index c087f97a8..300336ccb 100644 --- a/src/ICE/icedyn_rhg_vp.F90 +++ b/src/ICE/icedyn_rhg_vp.F90 @@ -332,7 +332,10 @@ CONTAINS v_oceU(ji,jj) = 0.25_wp * ( (v_oce(ji,jj) + v_oce(ji,jj-1)) + (v_oce(ji+1,jj) + v_oce(ji+1,jj-1)) ) * umask(ji,jj,1) ! Wind stress - ztaux_ai(ji,jj) = za_iU(ji,jj) * utau_ice(ji,jj) + ! Note the use of 0.5*(2-umask) in order to unmask the stress along coastlines + ! and the use of MAX(tmask(i,j),tmask(i+1,j) is to mask tau over ice shelves + ztaux_ai(ji,jj) = za_iU(ji,jj) * 0.5_wp * ( utau_ice(ji,jj) + utau_ice(ji+1,jj) ) * & + & ( 2. - umask(ji,jj,1) ) * MAX( tmask(ji,jj,1), tmask(ji+1,jj,1) ) ! Force due to sea surface tilt(- m*g*GRAD(ssh)) zspgU(ji,jj) = - zmassU * grav * ( zsshdyn(ji+1,jj) - zsshdyn(ji,jj) ) * r1_e1u(ji,jj) @@ -369,7 +372,10 @@ CONTAINS u_oceV(ji,jj) = 0.25_wp * ( (u_oce(ji,jj) + u_oce(ji-1,jj)) + (u_oce(ji,jj+1) + u_oce(ji-1,jj+1)) ) * vmask(ji,jj,1) ! Wind stress - ztauy_ai(ji,jj) = za_iV(ji,jj) * vtau_ice(ji,jj) + ! Note the use of 0.5*(2-umask) in order to unmask the stress along coastlines + ! and the use of MAX(tmask(i,j),tmask(i+1,j) is to mask tau over ice shelves + ztauy_ai(ji,jj) = za_iV(ji,jj) * 0.5_wp * ( vtau_ice(ji,jj) + vtau_ice(ji,jj+1) ) * & + & ( 2. - vmask(ji,jj,1) ) * MAX( tmask(ji,jj,1), tmask(ji,jj+1,1) ) ! Force due to sea surface tilt(- m*g*GRAD(ssh)) zspgV(ji,jj) = - zmassV * grav * ( zsshdyn(ji,jj+1) - zsshdyn(ji,jj) ) * r1_e2v(ji,jj) diff --git a/src/ICE/icesbc.F90 b/src/ICE/icesbc.F90 index 819154052..553e44997 100644 --- a/src/ICE/icesbc.F90 +++ b/src/ICE/icesbc.F90 @@ -272,15 +272,13 @@ CONTAINS zvel (ji,jj) = 0.5_wp * SQRT( ( u_ice(ji-1,jj ) + u_ice(ji,jj) ) * ( u_ice(ji-1,jj ) + u_ice(ji,jj) ) + & & ( v_ice(ji ,jj-1) + v_ice(ji,jj) ) * ( v_ice(ji ,jj-1) + v_ice(ji,jj) ) ) END_2D + CALL lbc_lnk( 'icesbc', zfric, 'T', 1.0_wp, zvel, 'T', 1.0_wp ) ELSE ! if no ice dynamics => transfer directly the atmospheric stress to the ocean - DO_2D( 0, 0, 0, 0 ) - zfric(ji,jj) = r1_rho0 * SQRT( 0.5_wp * & - & ( utau(ji,jj) * utau(ji,jj) + utau(ji-1,jj) * utau(ji-1,jj) & - & + vtau(ji,jj) * vtau(ji,jj) + vtau(ji,jj-1) * vtau(ji,jj-1) ) ) * tmask(ji,jj,1) - zvel(ji,jj) = 0._wp + DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) + zfric(ji,jj) = r1_rho0 * SQRT( utau(ji,jj)*utau(ji,jj) + vtau(ji,jj)*vtau(ji,jj) ) * tmask(ji,jj,1) + zvel (ji,jj) = 0._wp END_2D ENDIF - CALL lbc_lnk( 'icesbc', zfric, 'T', 1.0_wp, zvel, 'T', 1.0_wp ) ! !--------------------------------------------------------------------! ! Partial computation of forcing for the thermodynamic sea ice model diff --git a/src/ICE/icethd_do.F90 b/src/ICE/icethd_do.F90 index 3cd1bcd7b..ace1000aa 100644 --- a/src/ICE/icethd_do.F90 +++ b/src/ICE/icethd_do.F90 @@ -372,8 +372,8 @@ CONTAINS DO_2D( 0, 0, 0, 0 ) IF ( qlead(ji,jj) < 0._wp ) THEN ! cooling ! -- Wind stress -- ! - ztaux = ( utau_ice(ji-1,jj ) * umask(ji-1,jj ,1) + utau_ice(ji,jj) * umask(ji,jj,1) ) * 0.5_wp - ztauy = ( vtau_ice(ji ,jj-1) * vmask(ji ,jj-1,1) + vtau_ice(ji,jj) * vmask(ji,jj,1) ) * 0.5_wp + ztaux = utau_ice(ji,jj) * tmask(ji,jj,1) + ztauy = vtau_ice(ji,jj) * tmask(ji,jj,1) ! Square root of wind stress ztenagm = SQRT( SQRT( ztaux * ztaux + ztauy * ztauy ) ) diff --git a/src/ICE/iceupdate.F90 b/src/ICE/iceupdate.F90 index a973312f0..f6d047adf 100644 --- a/src/ICE/iceupdate.F90 +++ b/src/ICE/iceupdate.F90 @@ -314,7 +314,7 @@ CONTAINS !! !! ** Action : * at each ice time step (every nn_fsbc time step): !! - compute the modulus of ice-ocean relative velocity - !! (*rho*Cd) at T-point (C-grid) or I-point (B-grid) + !! (*rho*Cd) at T-point (C-grid) !! tmod_io = rhoco * | U_ice-U_oce | !! - update the modulus of stress at ocean surface !! taum = (1-a) * taum + a * tmod_io * | U_ice-U_oce | @@ -325,19 +325,19 @@ CONTAINS !! !! NB: - ice-ocean rotation angle no more allowed !! - here we make an approximation: taum is only computed every ice time step - !! This avoids mutiple average to pass from T -> U,V grids and next from U,V grids - !! to T grid. taum is used in TKE and GLS, which should not be too sensitive to this approximaton... + !! This avoids mutiple average to pass from U,V grids to T grids + !! taum is used in TKE and GLS, which should not be too sensitive to this approximaton... !! - !! ** Outputs : - utau, vtau : surface ocean i- and j-stress (u- & v-pts) updated with ice-ocean fluxes + !! ** Outputs : - utau, vtau : surface ocean i- and j-stress (T-pts) updated with ice-ocean fluxes !! - taum : modulus of the surface ocean stress (T-point) updated with ice-ocean fluxes !!--------------------------------------------------------------------- INTEGER , INTENT(in) :: kt ! ocean time-step index REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pu_oce, pv_oce ! surface ocean currents ! INTEGER :: ji, jj ! dummy loop indices - REAL(wp) :: zat_u, zutau_ice, zu_t, zmodt ! local scalar - REAL(wp) :: zat_v, zvtau_ice, zv_t, zrhoco ! - - - REAL(wp) :: zflagi ! - - + REAL(wp) :: zutau_ice, zu_t, zmodt ! local scalar + REAL(wp) :: zvtau_ice, zv_t, zrhoco ! - - + REAL(wp) :: zflagi ! - - !!--------------------------------------------------------------------- IF( ln_timing ) CALL timing_start('iceupdate') @@ -352,8 +352,8 @@ CONTAINS IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN !== Ice time-step only ==! (i.e. surface module time-step) DO_2D( 0, 0, 0, 0 ) !* update the modulus of stress at ocean surface (T-point) ! ! 2*(U_ice-U_oce) at T-point - zu_t = u_ice(ji,jj) + u_ice(ji-1,jj) - u_oce(ji,jj) - u_oce(ji-1,jj) - zv_t = v_ice(ji,jj) + v_ice(ji,jj-1) - v_oce(ji,jj) - v_oce(ji,jj-1) + zu_t = u_ice(ji,jj) + u_ice(ji-1,jj) - u_oce(ji,jj) - u_oce(ji-1,jj) ! u_oce = ssu_m + zv_t = v_ice(ji,jj) + v_ice(ji,jj-1) - v_oce(ji,jj) - v_oce(ji,jj-1) ! v_oce = ssv_m ! ! |U_ice-U_oce|^2 zmodt = 0.25_wp * ( zu_t * zu_t + zv_t * zv_t ) ! ! update the ocean stress modulus @@ -377,19 +377,14 @@ CONTAINS ENDIF ! DO_2D( 0, 0, 0, 0 ) !* update the stress WITHOUT an ice-ocean rotation angle - ! ice area at u and v-points - zat_u = ( at_i(ji,jj) * tmask(ji,jj,1) + at_i (ji+1,jj ) * tmask(ji+1,jj ,1) ) & - & / MAX( 1.0_wp , tmask(ji,jj,1) + tmask(ji+1,jj ,1) ) - zat_v = ( at_i(ji,jj) * tmask(ji,jj,1) + at_i (ji ,jj+1 ) * tmask(ji ,jj+1,1) ) & - & / MAX( 1.0_wp , tmask(ji,jj,1) + tmask(ji ,jj+1,1) ) ! ! linearized quadratic drag formulation - zutau_ice = 0.5_wp * ( tmod_io(ji,jj) + tmod_io(ji+1,jj) ) * ( u_ice(ji,jj) - pu_oce(ji,jj) ) - zvtau_ice = 0.5_wp * ( tmod_io(ji,jj) + tmod_io(ji,jj+1) ) * ( v_ice(ji,jj) - pv_oce(ji,jj) ) + zutau_ice = 0.5_wp * tmod_io(ji,jj) * ( u_ice(ji,jj) + u_ice(ji-1,jj) - pu_oce(ji,jj) - pu_oce(ji-1,jj) ) + zvtau_ice = 0.5_wp * tmod_io(ji,jj) * ( v_ice(ji,jj) + v_ice(ji,jj-1) - pv_oce(ji,jj) - pv_oce(ji,jj-1) ) ! ! stresses at the ocean surface - utau(ji,jj) = ( 1._wp - zat_u ) * utau_oce(ji,jj) + zat_u * zutau_ice - vtau(ji,jj) = ( 1._wp - zat_v ) * vtau_oce(ji,jj) + zat_v * zvtau_ice + utau(ji,jj) = ( 1._wp - at_i(ji,jj) ) * utau_oce(ji,jj) + at_i(ji,jj) * zutau_ice + vtau(ji,jj) = ( 1._wp - at_i(ji,jj) ) * vtau_oce(ji,jj) + at_i(ji,jj) * zvtau_ice END_2D - CALL lbc_lnk( 'iceupdate', utau, 'U', -1.0_wp, vtau, 'V', -1.0_wp ) ! lateral boundary condition + CALL lbc_lnk( 'iceupdate', utau, 'T', -1.0_wp, vtau, 'T', -1.0_wp ) ! lateral boundary condition ! IF( ln_timing ) CALL timing_stop('iceupdate') ! diff --git a/src/OCE/CRS/crsfld.F90 b/src/OCE/CRS/crsfld.F90 index b8115048d..efa11cf50 100644 --- a/src/OCE/CRS/crsfld.F90 +++ b/src/OCE/CRS/crsfld.F90 @@ -211,8 +211,8 @@ CONTAINS ! sbc fields CALL crs_dom_ope( ssh(:,:,Kmm) , 'VOL', 'T', tmask, sshn_crs , p_e12=e1e2t, p_e3=ze3t , psgn=1.0_wp ) - CALL crs_dom_ope( utau , 'SUM', 'U', umask, utau_crs , p_e12=e2u , p_surf_crs=e2u_crs , psgn=1.0_wp ) - CALL crs_dom_ope( vtau , 'SUM', 'V', vmask, vtau_crs , p_e12=e1v , p_surf_crs=e1v_crs , psgn=1.0_wp ) + CALL crs_dom_ope( utau , 'SUM', 'T', tmask, utau_crs , p_e12=e2t , p_surf_crs=e2t_crs , psgn=1.0_wp ) !clem tau: check psgn ?? + CALL crs_dom_ope( vtau , 'SUM', 'T', tmask, vtau_crs , p_e12=e1t , p_surf_crs=e1t_crs , psgn=1.0_wp ) ! CALL crs_dom_ope( wndm , 'SUM', 'T', tmask, wndm_crs , p_e12=e1e2t, p_surf_crs=e1e2t_crs, psgn=1.0_wp ) CALL crs_dom_ope( rnf , 'MAX', 'T', tmask, rnf_crs , psgn=1.0_wp ) CALL crs_dom_ope( qsr , 'SUM', 'T', tmask, qsr_crs , p_e12=e1e2t, p_surf_crs=e1e2t_crs, psgn=1.0_wp ) diff --git a/src/OCE/DIA/diawri.F90 b/src/OCE/DIA/diawri.F90 index db5da93a7..0deeb3370 100644 --- a/src/OCE/DIA/diawri.F90 +++ b/src/OCE/DIA/diawri.F90 @@ -868,7 +868,11 @@ CONTAINS CALL histdef( nid_T, "sohtc300", "Heat content 300 m" , "J/m2" , & ! htc3 & jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) #endif - + CALL histdef( nid_T, "sozotaux", "Wind Stress along i-axis" , "N/m2" , & ! utau + & jpi, jpj, nh_T, 1 , 1, 1 , - 99, 32, clop, zsto, zout ) + CALL histdef( nid_T, "sometauy", "Wind Stress along j-axis" , "N/m2" , & ! vtau + & jpi, jpj, nh_T, 1 , 1, 1 , - 99, 32, clop, zsto, zout ) + ! CALL histend( nid_T, snc4chunks=snc4set ) ! !!! nid_U : 3D @@ -878,10 +882,7 @@ CONTAINS CALL histdef( nid_U, "sdzocrtx", "Stokes Drift Zonal Current" , "m/s" , & ! usd & jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout ) ENDIF - ! !!! nid_U : 2D - CALL histdef( nid_U, "sozotaux", "Wind Stress along i-axis" , "N/m2" , & ! utau - & jpi, jpj, nh_U, 1 , 1, 1 , - 99, 32, clop, zsto, zout ) - + ! CALL histend( nid_U, snc4chunks=snc4set ) ! !!! nid_V : 3D @@ -891,10 +892,7 @@ CONTAINS CALL histdef( nid_V, "sdmecrty", "Stokes Drift Meridional Current" , "m/s" , & ! vsd & jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout ) ENDIF - ! !!! nid_V : 2D - CALL histdef( nid_V, "sometauy", "Wind Stress along j-axis" , "N/m2" , & ! vtau - & jpi, jpj, nh_V, 1 , 1, 1 , - 99, 32, clop, zsto, zout ) - + ! CALL histend( nid_V, snc4chunks=snc4set ) ! !!! nid_W : 3D @@ -1066,12 +1064,12 @@ CONTAINS CALL histwrite( nid_T, "so28chgt", it, hd28 , ndim_hT, ndex_hT ) ! depth of the 28 isotherm CALL histwrite( nid_T, "sohtc300", it, htc3 , ndim_hT, ndex_hT ) ! first 300m heaat content #endif + CALL histwrite( nid_T, "sozotaux", it, utau , ndim_hT, ndex_hT ) ! i-wind stress + CALL histwrite( nid_T, "sometauy", it, vtau , ndim_hT, ndex_hT ) ! j-wind stress CALL histwrite( nid_U, "vozocrtx", it, uu(:,:,:,Kmm) , ndim_U , ndex_U ) ! i-current - CALL histwrite( nid_U, "sozotaux", it, utau , ndim_hU, ndex_hU ) ! i-wind stress CALL histwrite( nid_V, "vomecrty", it, vv(:,:,:,Kmm) , ndim_V , ndex_V ) ! j-current - CALL histwrite( nid_V, "sometauy", it, vtau , ndim_hV, ndex_hV ) ! j-wind stress IF( ln_zad_Aimp ) THEN DO_3D( 0, 0, 0, 0, 1, jpk ) diff --git a/src/OCE/DYN/dynatf.F90 b/src/OCE/DYN/dynatf.F90 index 997692694..1a5df9967 100644 --- a/src/OCE/DYN/dynatf.F90 +++ b/src/OCE/DYN/dynatf.F90 @@ -331,7 +331,7 @@ CONTAINS ALLOCATE(zutau(jpi,jpj)) DO_2D( 0, 0, 0, 0 ) jk = miku(ji,jj) - zutau(ji,jj) = utau(ji,jj) + 0.5_wp * rho0 * ( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * puu(ji,jj,jk,Kaa) + zutau(ji,jj) = utau(ji,jj) + 0.5_wp * rho0 * rCdU_top(ji,jj) * ( puu(ji-1,jj,jk,Kaa) + puu(ji,jj,jk,Kaa) ) END_2D CALL iom_put( "utau", zutau(:,:) ) DEALLOCATE(zutau) @@ -345,7 +345,7 @@ CONTAINS ALLOCATE(zvtau(jpi,jpj)) DO_2D( 0, 0, 0, 0 ) jk = mikv(ji,jj) - zvtau(ji,jj) = vtau(ji,jj) + 0.5_wp * rho0 * ( rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) * pvv(ji,jj,jk,Kaa) + zvtau(ji,jj) = vtau(ji,jj) + 0.5_wp * rho0 * rCdU_top(ji,jj) * ( pvv(ji,jj-1,jk,Kaa) + pvv(ji,jj,jk,Kaa) ) END_2D CALL iom_put( "vtau", zvtau(:,:) ) DEALLOCATE(zvtau) diff --git a/src/OCE/DYN/dynatf_qco.F90 b/src/OCE/DYN/dynatf_qco.F90 index f9ed9f464..3dd46032e 100644 --- a/src/OCE/DYN/dynatf_qco.F90 +++ b/src/OCE/DYN/dynatf_qco.F90 @@ -248,7 +248,7 @@ CONTAINS ALLOCATE(zutau(jpi,jpj)) DO_2D( 0, 0, 0, 0 ) jk = miku(ji,jj) - zutau(ji,jj) = utau(ji,jj) + 0.5_wp * rho0 * ( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * puu(ji,jj,jk,Kaa) + zutau(ji,jj) = utau(ji,jj) + 0.5_wp * rho0 * rCdU_top(ji,jj) * ( puu(ji-1,jj,jk,Kaa) + puu(ji,jj,jk,Kaa) ) END_2D CALL iom_put( "utau", zutau(:,:) ) DEALLOCATE(zutau) @@ -262,7 +262,7 @@ CONTAINS ALLOCATE(zvtau(jpi,jpj)) DO_2D( 0, 0, 0, 0 ) jk = mikv(ji,jj) - zvtau(ji,jj) = vtau(ji,jj) + 0.5_wp * rho0 * ( rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) * pvv(ji,jj,jk,Kaa) + zvtau(ji,jj) = vtau(ji,jj) + 0.5_wp * rho0 * rCdU_top(ji,jj) * ( pvv(ji,jj-1,jk,Kaa) + pvv(ji,jj,jk,Kaa) ) END_2D CALL iom_put( "vtau", zvtau(:,:) ) DEALLOCATE(zvtau) diff --git a/src/OCE/DYN/dynspg_ts.F90 b/src/OCE/DYN/dynspg_ts.F90 index 83b2d99ca..16c2fcc1f 100644 --- a/src/OCE/DYN/dynspg_ts.F90 +++ b/src/OCE/DYN/dynspg_ts.F90 @@ -334,14 +334,14 @@ CONTAINS ! ! ------------------ ! IF( ln_bt_fw ) THEN DO_2D( 0, 0, 0, 0 ) - zu_frc(ji,jj) = zu_frc(ji,jj) + r1_rho0 * utau(ji,jj) * r1_hu(ji,jj,Kmm) - zv_frc(ji,jj) = zv_frc(ji,jj) + r1_rho0 * vtau(ji,jj) * r1_hv(ji,jj,Kmm) + zu_frc(ji,jj) = zu_frc(ji,jj) + r1_rho0 * utauU(ji,jj) * r1_hu(ji,jj,Kmm) + zv_frc(ji,jj) = zv_frc(ji,jj) + r1_rho0 * vtauV(ji,jj) * r1_hv(ji,jj,Kmm) END_2D ELSE zztmp = r1_rho0 * r1_2 DO_2D( 0, 0, 0, 0 ) - zu_frc(ji,jj) = zu_frc(ji,jj) + zztmp * ( utau_b(ji,jj) + utau(ji,jj) ) * r1_hu(ji,jj,Kmm) - zv_frc(ji,jj) = zv_frc(ji,jj) + zztmp * ( vtau_b(ji,jj) + vtau(ji,jj) ) * r1_hv(ji,jj,Kmm) + zu_frc(ji,jj) = zu_frc(ji,jj) + zztmp * ( utau_b(ji,jj) + utauU(ji,jj) ) * r1_hu(ji,jj,Kmm) + zv_frc(ji,jj) = zv_frc(ji,jj) + zztmp * ( vtau_b(ji,jj) + vtauV(ji,jj) ) * r1_hv(ji,jj,Kmm) END_2D ENDIF ! diff --git a/src/OCE/DYN/dynzdf.F90 b/src/OCE/DYN/dynzdf.F90 index 1e0260219..98c65e460 100644 --- a/src/OCE/DYN/dynzdf.F90 +++ b/src/OCE/DYN/dynzdf.F90 @@ -267,10 +267,10 @@ CONTAINS DO_2D( 0, 0, 0, 0 ) !== second recurrence: SOLk = RHSk - Lk / Dk-1 Lk-1 ==! #if defined key_RK3 ! ! RK3: use only utau (not utau_b) - puu(ji,jj,1,Kaa) = puu(ji,jj,1,Kaa) + rDt * utau(ji,jj) & + puu(ji,jj,1,Kaa) = puu(ji,jj,1,Kaa) + rDt * utauU(ji,jj) & & / ( e3u(ji,jj,1,Kaa) * rho0 ) * umask(ji,jj,1) #else - puu(ji,jj,1,Kaa) = puu(ji,jj,1,Kaa) + zDt_2 * ( utau_b(ji,jj) + utau(ji,jj) ) & + puu(ji,jj,1,Kaa) = puu(ji,jj,1,Kaa) + zDt_2 * ( utau_b(ji,jj) + utauU(ji,jj) ) & & / ( e3u(ji,jj,1,Kaa) * rho0 ) * umask(ji,jj,1) #endif END_2D @@ -397,10 +397,10 @@ CONTAINS DO_2D( 0, 0, 0, 0 ) !== second recurrence: SOLk = RHSk - Lk / Dk-1 Lk-1 ==! #if defined key_RK3 ! ! RK3: use only vtau (not vtau_b) - pvv(ji,jj,1,Kaa) = pvv(ji,jj,1,Kaa) + rDt * vtau(ji,jj) & + pvv(ji,jj,1,Kaa) = pvv(ji,jj,1,Kaa) + rDt * vtauV(ji,jj) & & / ( e3v(ji,jj,1,Kaa) * rho0 ) * vmask(ji,jj,1) #else - pvv(ji,jj,1,Kaa) = pvv(ji,jj,1,Kaa) + zDt_2*( vtau_b(ji,jj) + vtau(ji,jj) ) & + pvv(ji,jj,1,Kaa) = pvv(ji,jj,1,Kaa) + zDt_2 * ( vtau_b(ji,jj) + vtauV(ji,jj) ) & & / ( e3v(ji,jj,1,Kaa) * rho0 ) * vmask(ji,jj,1) #endif END_2D diff --git a/src/OCE/SBC/cpl_oasis3.F90 b/src/OCE/SBC/cpl_oasis3.F90 index 0147ffdac..3c7e15ba4 100644 --- a/src/OCE/SBC/cpl_oasis3.F90 +++ b/src/OCE/SBC/cpl_oasis3.F90 @@ -419,6 +419,7 @@ CONTAINS IF( .NOT. ll_1st ) THEN CALL lbc_lnk( 'cpl_oasis3', pdata(:,:,jc), srcv(kid)%clgrid, srcv(kid)%nsgn ) ENDIF + !!clem: mettre T instead of clgrid ENDDO ! diff --git a/src/OCE/SBC/sbc_ice.F90 b/src/OCE/SBC/sbc_ice.F90 index acc6adf1f..cac087f79 100644 --- a/src/OCE/SBC/sbc_ice.F90 +++ b/src/OCE/SBC/sbc_ice.F90 @@ -50,8 +50,8 @@ MODULE sbc_ice REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: qcn_ice !: heat conduction flux in the layer below surface [W/m2] REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: qtr_ice_top !: solar flux transmitted below the ice surface [W/m2] - REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: utau_ice !: atmos-ice u-stress. VP: I-pt ; EVP: U,V-pts [N/m2] - REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: vtau_ice !: atmos-ice v-stress. VP: I-pt ; EVP: U,V-pts [N/m2] + REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: utau_ice !: atmos-ice u-stress. T-pts [N/m2] + REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: vtau_ice !: atmos-ice v-stress. T-pts [N/m2] REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: emp_ice !: sublimation - precip over sea ice [kg/m2/s] REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: topmelt !: category topmelt diff --git a/src/OCE/SBC/sbc_oce.F90 b/src/OCE/SBC/sbc_oce.F90 index 3098b5865..f4cd6c97e 100644 --- a/src/OCE/SBC/sbc_oce.F90 +++ b/src/OCE/SBC/sbc_oce.F90 @@ -103,34 +103,36 @@ MODULE sbc_oce INTEGER , PUBLIC :: ncpl_qsr_freq = 0 !: qsr coupling frequency per days from atmosphere (used by top) ! !! !! now ! before !! - REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: utau , utau_b !: sea surface i-stress (ocean referential) [N/m2] - REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: vtau , vtau_b !: sea surface j-stress (ocean referential) [N/m2] - REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: utau_icb, vtau_icb !: sea surface (i,j)-stress used by icebergs [N/m2] - REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: taum !: module of sea surface stress (at T-point) [N/m2] + REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: utau !: sea surface i-stress (ocean referential) T-pt [N/m2] + REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: vtau !: sea surface j-stress (ocean referential) T-pt [N/m2] + REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: utauU , utau_b !: sea surface i-stress (ocean referential) U-pt [N/m2] + REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: vtauV , vtau_b !: sea surface j-stress (ocean referential) V-pt [N/m2] + REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: utau_icb, vtau_icb !: sea surface (i,j)-stress used by icebergs [N/m2] + REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: taum !: module of sea surface stress (at T-point) [N/m2] !! wndm is used compute surface gases exchanges in ice-free ocean or leads - REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: wndm !: wind speed module at T-point (=|U10m-Uoce|) [m/s] - REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: rhoa !: air density at "rn_zu" m above the sea [kg/m3] - REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: qsr !: sea heat flux: solar [W/m2] - REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: qns , qns_b !: sea heat flux: non solar [W/m2] - REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: qsr_tot !: total solar heat flux (over sea and ice) [W/m2] - REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: qns_tot !: total non solar heat flux (over sea and ice) [W/m2] - REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: emp , emp_b !: freshwater budget: volume flux [Kg/m2/s] - REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: sfx , sfx_b !: salt flux [PSS.kg/m2/s] - REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: emp_tot !: total E-P over ocean and ice [Kg/m2/s] - REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: fmmflx !: freshwater budget: freezing/melting [Kg/m2/s] - REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: rnf , rnf_b !: river runoff [Kg/m2/s] - REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: fwficb , fwficb_b !: iceberg melting [Kg/m2/s] + REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: wndm !: wind speed module at T-point (=|U10m-Uoce|) [m/s] + REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: rhoa !: air density at "rn_zu" m above the sea [kg/m3] + REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: qsr !: sea heat flux: solar [W/m2] + REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: qns , qns_b !: sea heat flux: non solar [W/m2] + REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: qsr_tot !: total solar heat flux (over sea and ice) [W/m2] + REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: qns_tot !: total non solar heat flux (over sea and ice) [W/m2] + REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: emp , emp_b !: freshwater budget: volume flux [Kg/m2/s] + REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: sfx , sfx_b !: salt flux [PSS.kg/m2/s] + REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: emp_tot !: total E-P over ocean and ice [Kg/m2/s] + REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: fmmflx !: freshwater budget: freezing/melting [Kg/m2/s] + REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: rnf , rnf_b !: river runoff [Kg/m2/s] + REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: fwficb , fwficb_b !: iceberg melting [Kg/m2/s] !! REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: sbc_tsc, sbc_tsc_b !: sbc content trend [K.m/s] jpi,jpj,jpts REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: qsr_hc , qsr_hc_b !: heat content trend due to qsr flux [K.m/s] jpi,jpj,jpk !! !! - REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: tprecip !: total precipitation [Kg/m2/s] - REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: sprecip !: solid precipitation [Kg/m2/s] - REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: fr_i !: ice fraction = 1 - lead fraction (between 0 to 1) - REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: atm_co2 !: atmospheric pCO2 [ppm] - REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: xcplmask !: coupling mask for ln_mixcpl (warning: allocated in sbccpl) - REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: cloud_fra !: cloud cover (fraction of cloud in a gridcell) [-] + REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: tprecip !: total precipitation [Kg/m2/s] + REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: sprecip !: solid precipitation [Kg/m2/s] + REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: fr_i !: ice fraction = 1 - lead fraction (between 0 to 1) + REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: atm_co2 !: atmospheric pCO2 [ppm] + REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: xcplmask !: coupling mask for ln_mixcpl (warning: allocated in sbccpl) + REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: cloud_fra !: cloud cover (fraction of cloud in a gridcell) [-] !!--------------------------------------------------------------------- !! ABL Vertical Domain size @@ -177,8 +179,8 @@ CONTAINS !!--------------------------------------------------------------------- ierr(:) = 0 ! - ALLOCATE( utau(jpi,jpj) , utau_b(jpi,jpj) , taum(jpi,jpj) , & - & vtau(jpi,jpj) , vtau_b(jpi,jpj) , wndm(jpi,jpj) , rhoa(jpi,jpj) , STAT=ierr(1) ) + ALLOCATE( utau(jpi,jpj) , utau_b(jpi,jpj) , utauU(jpi,jpj) , taum(jpi,jpj) , & + & vtau(jpi,jpj) , vtau_b(jpi,jpj) , vtauV(jpi,jpj) , wndm(jpi,jpj) , rhoa(jpi,jpj) , STAT=ierr(1) ) ! ALLOCATE( qns_tot(jpi,jpj) , qns (jpi,jpj) , qns_b(jpi,jpj), & & qsr_tot(jpi,jpj) , qsr (jpi,jpj) , & @@ -205,9 +207,10 @@ CONTAINS END FUNCTION sbc_oce_alloc + !!clem => this subroutine is never used in nemo SUBROUTINE sbc_tau2wnd !!--------------------------------------------------------------------- - !! *** ROUTINE sbc_tau2wnd *** + !! *** ROUTINE *** !! !! ** Purpose : Estimation of wind speed as a function of wind stress !! @@ -217,17 +220,14 @@ CONTAINS USE lbclnk ! ocean lateral boundary conditions (or mpp link) REAL(wp) :: zrhoa = 1.22 ! Air density kg/m3 REAL(wp) :: zcdrag = 1.5e-3 ! drag coefficient - REAL(wp) :: ztx, zty, ztau, zcoef ! temporary variables + REAL(wp) :: ztau, zcoef ! temporary variables INTEGER :: ji, jj ! dummy indices !!--------------------------------------------------------------------- zcoef = 0.5 / ( zrhoa * zcdrag ) - DO_2D( 0, 0, 0, 0 ) - ztx = utau(ji-1,jj ) + utau(ji,jj) - zty = vtau(ji ,jj-1) + vtau(ji,jj) - ztau = SQRT( ztx * ztx + zty * zty ) + DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) + ztau = SQRT( utau(ji,jj)*utau(ji,jj) + vtau(ji,jj)*vtau(ji,jj) ) wndm(ji,jj) = SQRT ( ztau * zcoef ) * tmask(ji,jj,1) END_2D - CALL lbc_lnk( 'sbc_oce', wndm(:,:) , 'T', 1.0_wp ) ! END SUBROUTINE sbc_tau2wnd diff --git a/src/OCE/SBC/sbcblk.F90 b/src/OCE/SBC/sbcblk.F90 index aa3668c33..397155a58 100644 --- a/src/OCE/SBC/sbcblk.F90 +++ b/src/OCE/SBC/sbcblk.F90 @@ -493,7 +493,7 @@ CONTAINS !! the stress is assumed to be in the (i,j) mesh referential !! !! ** Action : defined at each time-step at the air-sea interface - !! - utau, vtau i- and j-component of the wind stress + !! - utau, vtau i- and j-component of the wind stress at T-point !! - taum wind stress module at T-point !! - wndm wind speed module at T-point over free ocean or leads in presence of sea-ice !! - qns, qsr non-solar and solar heat fluxes @@ -611,10 +611,10 @@ CONTAINS END SUBROUTINE sbc_blk - SUBROUTINE blk_oce_1( kt, pwndi, pwndj, ptair, pqair, & ! inp - & pslp , pst , pu , pv, & ! inp - & puatm, pvatm, pdqsr , pdqlw , & ! inp - & ptsk , pssq , pcd_du, psen, plat, pevp ) ! out + SUBROUTINE blk_oce_1( kt, pwndi, pwndj, ptair, pqair, & ! <<= in + & pslp , pst , pu , pv, & ! <<= in + & puatm, pvatm, pdqsr , pdqlw , & ! <<= in + & ptsk , pssq , pcd_du, psen, plat, pevp ) ! =>> out !!--------------------------------------------------------------------- !! *** ROUTINE blk_oce_1 *** !! @@ -657,7 +657,6 @@ CONTAINS #if defined key_cyclone REAL(wp), DIMENSION(jpi,jpj) :: zwnd_i, zwnd_j ! wind speed components at T-point #endif - REAL(wp), DIMENSION(jpi,jpj) :: ztau_i, ztau_j ! wind stress components at T-point REAL(wp), DIMENSION(jpi,jpj) :: zU_zu ! bulk wind speed at height zu [m/s] REAL(wp), DIMENSION(jpi,jpj) :: zcd_oce ! momentum transfert coefficient over ocean REAL(wp), DIMENSION(jpi,jpj) :: zch_oce ! sensible heat transfert coefficient over ocean @@ -695,22 +694,20 @@ CONTAINS #else ! ... scalar wind module at T-point (not masked) DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) - wndm(ji,jj) = SQRT( pwndi(ji,jj) * pwndi(ji,jj) + pwndj(ji,jj) * pwndj(ji,jj) ) + wndm(ji,jj) = SQRT( pwndi(ji,jj) * pwndi(ji,jj) + pwndj(ji,jj) * pwndj(ji,jj) ) END_2D #endif ! ----------------------------------------------------------------------------- ! ! I Solar FLUX ! ! ----------------------------------------------------------------------------- ! - ! ocean albedo assumed to be constant + modify now Qsr to include the diurnal cycle ! Short Wave - zztmp = 1. - albo + ! ocean albedo assumed to be constant + modify now Qsr to include the diurnal cycle IF( ln_dm2dc ) THEN - qsr(:,:) = zztmp * sbc_dcy( pdqsr(:,:) ) * tmask(:,:,1) + qsr(:,:) = ( 1._wp - albo ) * sbc_dcy( pdqsr(:,:) ) * tmask(:,:,1) ELSE - qsr(:,:) = zztmp * pdqsr(:,:) * tmask(:,:,1) + qsr(:,:) = ( 1._wp - albo ) * pdqsr(:,:) * tmask(:,:,1) ENDIF - ! ----------------------------------------------------------------------------- ! ! II Turbulent FLUXES ! ! ----------------------------------------------------------------------------- ! @@ -718,69 +715,62 @@ CONTAINS ! specific humidity at SST pssq(:,:) = rdct_qsat_salt * q_sat( ptsk(:,:), pslp(:,:) ) + ! Backup "bulk SST" and associated spec. hum. IF( ln_skin_cs .OR. ln_skin_wl ) THEN - !! Backup "bulk SST" and associated spec. hum. zztmp1(:,:) = zsspt(:,:) - zztmp2(:,:) = pssq(:,:) + zztmp2(:,:) = pssq (:,:) ENDIF - !! Time to call the user-selected bulk parameterization for - !! == transfer coefficients ==! Cd, Ch, Ce at T-point, and more... - SELECT CASE( nblk ) - + ! transfer coefficients (Cd, Ch, Ce at T-point, and more) + SELECT CASE( nblk ) ! user-selected bulk parameterization + ! CASE( np_NCAR ) CALL turb_ncar ( rn_zqt, rn_zu, zsspt, ptair, pssq, pqair, wndm, & - & zcd_oce, zch_oce, zce_oce, theta_zu, q_zu, zU_zu , & + & zcd_oce, zch_oce, zce_oce, theta_zu, q_zu, zU_zu , & & nb_iter=nn_iter_algo ) - ! CASE( np_COARE_3p0 ) CALL turb_coare3p0( kt, rn_zqt, rn_zu, zsspt, ptair, pssq, pqair, wndm, & - & ln_skin_cs, ln_skin_wl, & - & zcd_oce, zch_oce, zce_oce, theta_zu, q_zu, zU_zu, & - & nb_iter=nn_iter_algo, & + & ln_skin_cs, ln_skin_wl, & + & zcd_oce, zch_oce, zce_oce, theta_zu, q_zu, zU_zu, & + & nb_iter=nn_iter_algo, & & Qsw=qsr(:,:), rad_lw=pdqlw(:,:), slp=pslp(:,:) ) - ! CASE( np_COARE_3p6 ) CALL turb_coare3p6( kt, rn_zqt, rn_zu, zsspt, ptair, pssq, pqair, wndm, & - & ln_skin_cs, ln_skin_wl, & - & zcd_oce, zch_oce, zce_oce, theta_zu, q_zu, zU_zu, & - & nb_iter=nn_iter_algo, & + & ln_skin_cs, ln_skin_wl, & + & zcd_oce, zch_oce, zce_oce, theta_zu, q_zu, zU_zu, & + & nb_iter=nn_iter_algo, & & Qsw=qsr(:,:), rad_lw=pdqlw(:,:), slp=pslp(:,:) ) - ! CASE( np_ECMWF ) CALL turb_ecmwf ( kt, rn_zqt, rn_zu, zsspt, ptair, pssq, pqair, wndm, & - & ln_skin_cs, ln_skin_wl, & - & zcd_oce, zch_oce, zce_oce, theta_zu, q_zu, zU_zu, & - & nb_iter=nn_iter_algo, & + & ln_skin_cs, ln_skin_wl, & + & zcd_oce, zch_oce, zce_oce, theta_zu, q_zu, zU_zu, & + & nb_iter=nn_iter_algo, & & Qsw=qsr(:,:), rad_lw=pdqlw(:,:), slp=pslp(:,:) ) - ! CASE( np_ANDREAS ) CALL turb_andreas ( rn_zqt, rn_zu, zsspt, ptair, pssq, pqair, wndm, & - & zcd_oce, zch_oce, zce_oce, theta_zu, q_zu, zU_zu , & + & zcd_oce, zch_oce, zce_oce, theta_zu, q_zu, zU_zu, & & nb_iter=nn_iter_algo ) - ! CASE DEFAULT CALL ctl_stop( 'STOP', 'sbc_oce: non-existing bulk parameterizaton selected' ) - ! END SELECT - IF( iom_use('Cd_oce') ) CALL iom_put("Cd_oce", zcd_oce * tmask(:,:,1)) - IF( iom_use('Ce_oce') ) CALL iom_put("Ce_oce", zce_oce * tmask(:,:,1)) - IF( iom_use('Ch_oce') ) CALL iom_put("Ch_oce", zch_oce * tmask(:,:,1)) + ! outputs + IF( iom_use('Cd_oce') ) CALL iom_put( "Cd_oce", zcd_oce * tmask(:,:,1) ) + IF( iom_use('Ce_oce') ) CALL iom_put( "Ce_oce", zce_oce * tmask(:,:,1) ) + IF( iom_use('Ch_oce') ) CALL iom_put( "Ch_oce", zch_oce * tmask(:,:,1) ) !! LB: mainly here for debugging purpose: - IF( iom_use('theta_zt') ) CALL iom_put("theta_zt", (ptair-rt0) * tmask(:,:,1)) ! potential temperature at z=zt - IF( iom_use('q_zt') ) CALL iom_put("q_zt", pqair * tmask(:,:,1)) ! specific humidity " - IF( iom_use('theta_zu') ) CALL iom_put("theta_zu", (theta_zu -rt0) * tmask(:,:,1)) ! potential temperature at z=zu - IF( iom_use('q_zu') ) CALL iom_put("q_zu", q_zu * tmask(:,:,1)) ! specific humidity " - IF( iom_use('ssq') ) CALL iom_put("ssq", pssq * tmask(:,:,1)) ! saturation specific humidity at z=0 - IF( iom_use('wspd_blk') ) CALL iom_put("wspd_blk", zU_zu * tmask(:,:,1)) ! bulk wind speed at z=zu - + IF( iom_use('theta_zt') ) CALL iom_put( "theta_zt", (ptair-rt0) * tmask(:,:,1) ) ! potential temperature at z=zt + IF( iom_use('q_zt') ) CALL iom_put( "q_zt", pqair * tmask(:,:,1) ) ! specific humidity " + IF( iom_use('theta_zu') ) CALL iom_put( "theta_zu", (theta_zu -rt0) * tmask(:,:,1) ) ! potential temperature at z=zu + IF( iom_use('q_zu') ) CALL iom_put( "q_zu", q_zu * tmask(:,:,1) ) ! specific humidity " + IF( iom_use('ssq') ) CALL iom_put( "ssq", pssq * tmask(:,:,1) ) ! saturation specific humidity at z=0 + IF( iom_use('wspd_blk') ) CALL iom_put( "wspd_blk", zU_zu * tmask(:,:,1) ) ! bulk wind speed at z=zu + + ! In the presence of sea-ice we do not use the cool-skin/warm-layer update of zsspt, pssq & ptsk from turb_*() IF( ln_skin_cs .OR. ln_skin_wl ) THEN - !! In the presence of sea-ice we forget about the cool-skin/warm-layer update of zsspt, pssq & ptsk: WHERE ( fr_i(:,:) > 0.001_wp ) - ! sea-ice present, we forget about the update, using what we backed up before call to turb_*() zsspt(:,:) = zztmp1(:,:) - pssq(:,:) = zztmp2(:,:) + pssq (:,:) = zztmp2(:,:) END WHERE ! apply potential temperature increment to abolute SST ptsk(:,:) = ptsk(:,:) + ( zsspt(:,:) - zztmp1(:,:) ) @@ -809,10 +799,10 @@ CONTAINS END_2D CALL BULK_FORMULA( rn_zu, zsspt(:,:), pssq(:,:), theta_zu(:,:), q_zu(:,:), & - & zcd_oce(:,:), zch_oce(:,:), zce_oce(:,:), & - & wndm(:,:), zU_zu(:,:), pslp(:,:), rhoa(:,:), & - & taum(:,:), psen(:,:), plat(:,:), & - & pEvap=pevp(:,:), pfact_evap=rn_efac ) + & zcd_oce(:,:), zch_oce(:,:), zce_oce(:,:), & + & wndm(:,:), zU_zu(:,:), pslp(:,:), rhoa(:,:), & + & taum(:,:), psen(:,:), plat(:,:), & + & pEvap=pevp(:,:), pfact_evap=rn_efac ) psen(:,:) = psen(:,:) * tmask(:,:,1) plat(:,:) = plat(:,:) * tmask(:,:,1) @@ -821,57 +811,42 @@ CONTAINS DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) IF( wndm(ji,jj) > 0._wp ) THEN - zztmp = taum(ji,jj) / wndm(ji,jj) + zztmp = taum(ji,jj) / wndm(ji,jj) #if defined key_cyclone - ztau_i(ji,jj) = zztmp * zwnd_i(ji,jj) - ztau_j(ji,jj) = zztmp * zwnd_j(ji,jj) + utau(ji,jj) = zztmp * zwnd_i(ji,jj) + vtau(ji,jj) = zztmp * zwnd_j(ji,jj) #else - ztau_i(ji,jj) = zztmp * pwndi(ji,jj) - ztau_j(ji,jj) = zztmp * pwndj(ji,jj) + utau(ji,jj) = zztmp * pwndi(ji,jj) + vtau(ji,jj) = zztmp * pwndj(ji,jj) #endif ELSE - ztau_i(ji,jj) = 0._wp - ztau_j(ji,jj) = 0._wp + utau(ji,jj) = 0._wp + vtau(ji,jj) = 0._wp ENDIF END_2D IF( ln_crt_fbk ) THEN ! aply eq. 10 and 11 of Renault et al. 2020 (doi: 10.1029/2019MS001715) zstmax = MIN( rn_stau_a * 3._wp + rn_stau_b, 0._wp ) ! set the max value of Stau corresponding to a wind of 3 m/s (<0) - DO_2D( 0, 1, 0, 1 ) ! end at jpj and jpi, as ztau_j(ji,jj+1) ztau_i(ji+1,jj) used in the next loop - zstau = MIN( rn_stau_a * wndm(ji,jj) + rn_stau_b, zstmax ) ! stau (<0) must be smaller than zstmax - ztau_i(ji,jj) = ztau_i(ji,jj) + zstau * ( 0.5_wp * ( pu(ji-1,jj ) + pu(ji,jj) ) - puatm(ji,jj) ) - ztau_j(ji,jj) = ztau_j(ji,jj) + zstau * ( 0.5_wp * ( pv(ji ,jj-1) + pv(ji,jj) ) - pvatm(ji,jj) ) - taum(ji,jj) = SQRT( ztau_i(ji,jj) * ztau_i(ji,jj) + ztau_j(ji,jj) * ztau_j(ji,jj) ) + DO_2D( 0, 0, 0, 0 ) + zstau = MIN( rn_stau_a * wndm(ji,jj) + rn_stau_b, zstmax ) * tmask(ji,jj,1) ! stau (<0) must be smaller than zstmax + utau(ji,jj) = utau(ji,jj) + zstau * ( 0.5_wp * ( pu(ji-1,jj ) + pu(ji,jj) ) - puatm(ji,jj) ) + vtau(ji,jj) = vtau(ji,jj) + zstau * ( 0.5_wp * ( pv(ji ,jj-1) + pv(ji,jj) ) - pvatm(ji,jj) ) + taum(ji,jj) = SQRT( utau(ji,jj) * utau(ji,jj) + vtau(ji,jj) * vtau(ji,jj) ) END_2D + CALL lbc_lnk( 'sbcblk', utau, 'T', -1._wp, vtau, 'T', -1._wp, taum, 'T', 1._wp ) ENDIF - ! ... utau, vtau at U- and V_points, resp. - ! Note the use of 0.5*(2-umask) in order to unmask the stress along coastlines - ! Note that coastal wind stress is not used in the code... so this extra care has no effect - DO_2D( 0, 0, 0, 0 ) ! start loop at 2, in case ln_crt_fbk = T - utau(ji,jj) = 0.5 * ( 2. - umask(ji,jj,1) ) * ( ztau_i(ji,jj) + ztau_i(ji+1,jj ) ) & - & * MAX(tmask(ji,jj,1),tmask(ji+1,jj,1)) - vtau(ji,jj) = 0.5 * ( 2. - vmask(ji,jj,1) ) * ( ztau_j(ji,jj) + ztau_j(ji ,jj+1) ) & - & * MAX(tmask(ji,jj,1),tmask(ji,jj+1,1)) - END_2D - - IF( ln_crt_fbk ) THEN - CALL lbc_lnk( 'sbcblk', utau, 'U', -1._wp, vtau, 'V', -1._wp, taum, 'T', 1._wp ) - ELSE - CALL lbc_lnk( 'sbcblk', utau, 'U', -1._wp, vtau, 'V', -1._wp ) - ENDIF - - ! Saving open-ocean wind-stress (module and components) on T-points: - CALL iom_put( "taum_oce", taum(:,:)*tmask(:,:,1) ) ! output wind stress module - !#LB: These 2 lines below mostly here for 'STATION_ASF' test-case, otherwize "utau" (U-grid) and vtau" (V-grid) does the job in: [DYN/dynatf.F90]) - CALL iom_put( "utau_oce", ztau_i(:,:)*tmask(:,:,1) ) ! utau at T-points! - CALL iom_put( "vtau_oce", ztau_j(:,:)*tmask(:,:,1) ) ! vtau at T-points! + ! Saving open-ocean wind-stress (module and components) + CALL iom_put( "taum_oce", taum(:,:) ) ! wind stress module + ! ! LB: These 2 lines below mostly here for 'STATION_ASF' test-case + CALL iom_put( "utau_oce", utau(:,:) ) ! utau + CALL iom_put( "vtau_oce", vtau(:,:) ) ! vtau IF(sn_cfctl%l_prtctl) THEN CALL prt_ctl( tab2d_1=pssq , clinfo1=' blk_oce_1: pssq : ', mask1=tmask ) CALL prt_ctl( tab2d_1=wndm , clinfo1=' blk_oce_1: wndm : ', mask1=tmask ) - CALL prt_ctl( tab2d_1=utau , clinfo1=' blk_oce_1: utau : ', mask1=umask, & - & tab2d_2=vtau , clinfo2=' vtau : ', mask2=vmask ) + CALL prt_ctl( tab2d_1=utau , clinfo1=' blk_oce_1: utau : ', mask1=tmask, & + & tab2d_2=vtau , clinfo2=' vtau : ', mask2=tmask ) CALL prt_ctl( tab2d_1=zcd_oce, clinfo1=' blk_oce_1: Cd : ', mask1=tmask ) ENDIF ! @@ -896,8 +871,8 @@ CONTAINS !! at the ocean surface at each time step knowing Cd, Ch, Ce and !! atmospheric variables (from ABL or external data) !! - !! ** Outputs : - utau : i-component of the stress at U-point (N/m2) - !! - vtau : j-component of the stress at V-point (N/m2) + !! ** Outputs : - utau : i-component of the stress at T-point (N/m2) + !! - vtau : j-component of the stress at T-point (N/m2) !! - taum : Wind stress module at T-point (N/m2) !! - wndm : Wind speed module at T-point (m/s) !! - qsr : Solar heat flux over the ocean (W/m2) @@ -1015,11 +990,16 @@ CONTAINS REAL(wp) , INTENT( out), DIMENSION(:,: ), OPTIONAL :: pcd_dui ! if ln_abl ! INTEGER :: ji, jj ! dummy loop indices - REAL(wp) :: zootm_su ! sea-ice surface mean temperature - REAL(wp) :: zztmp1, zztmp2 ! temporary scalars - REAL(wp), DIMENSION(jpi,jpj) :: ztmp, zsipt ! temporary array + REAL(wp) :: zztmp ! temporary scalars + REAL(wp), DIMENSION(jpi,jpj) :: ztmp, zsipt ! temporary array + REAL(wp), DIMENSION(jpi,jpj) :: zmsk00 ! O% concentration ice mask !!--------------------------------------------------------------------- ! + ! treshold for outputs + DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) + zmsk00(ji,jj) = MAX( 0._wp , SIGN( 1._wp , fr_i(ji,jj) - 1.e-6_wp ) ) ! 1 if ice, 0 if no ice + END_2D + ! ------------------------------------------------------------ ! ! Wind module relative to the moving ice ( U10m - U_ice ) ! ! ------------------------------------------------------------ ! @@ -1032,9 +1012,9 @@ CONTAINS zsipt(:,:) = theta_exner( ptsui(:,:), pslp(:,:) ) ! sea-ice <-> atmosphere bulk transfer coefficients - SELECT CASE( nblk_ice ) - - CASE( np_ice_cst ) + SELECT CASE( nblk_ice ) ! user-selected bulk parameterization + ! + CASE( np_ice_cst ) ! Constant bulk transfer coefficients over sea-ice: Cd_ice(:,:) = rn_Cd_i Ch_ice(:,:) = rn_Ch_i @@ -1042,62 +1022,45 @@ CONTAINS ! no height adjustment, keeping zt values: theta_zu_i(:,:) = ptair(:,:) q_zu_i(:,:) = pqair(:,:) - - CASE( np_ice_an05 ) ! calculate new drag from Lupkes(2015) equations + ! + CASE( np_ice_an05 ) ! from Andreas(2005) equations ztmp(:,:) = q_sat( ptsui(:,:), pslp(:,:), l_ice=.TRUE. ) ! temporary array for SSQ CALL turb_ice_an05( rn_zqt, rn_zu, zsipt, ptair, ztmp, pqair, wndm_ice, & & Cd_ice, Ch_ice, Ce_ice, theta_zu_i, q_zu_i ) - !! - CASE( np_ice_lu12 ) + ! + CASE( np_ice_lu12 ) ! from Lupkes(2012) equations ztmp(:,:) = q_sat( ptsui(:,:), pslp(:,:), l_ice=.TRUE. ) ! temporary array for SSQ CALL turb_ice_lu12( rn_zqt, rn_zu, zsipt, ptair, ztmp, pqair, wndm_ice, fr_i, & & Cd_ice, Ch_ice, Ce_ice, theta_zu_i, q_zu_i ) - !! - CASE( np_ice_lg15 ) ! calculate new drag from Lupkes(2015) equations + ! + CASE( np_ice_lg15 ) ! from Lupkes and Gryanik (2015) equations ztmp(:,:) = q_sat( ptsui(:,:), pslp(:,:), l_ice=.TRUE. ) ! temporary array for SSQ CALL turb_ice_lg15( rn_zqt, rn_zu, zsipt, ptair, ztmp, pqair, wndm_ice, fr_i, & & Cd_ice, Ch_ice, Ce_ice, theta_zu_i, q_zu_i ) - !! + ! END SELECT - IF( iom_use('Cd_ice').OR.iom_use('Ce_ice').OR.iom_use('Ch_ice').OR.iom_use('taum_ice').OR.iom_use('utau_ice').OR.iom_use('vtau_ice') ) & - & ztmp(:,:) = ( 1._wp - MAX(0._wp, SIGN( 1._wp, 1.E-6_wp - fr_i )) )*tmask(:,:,1) ! mask for presence of ice ! - - IF( iom_use('Cd_ice') ) CALL iom_put("Cd_ice", Cd_ice*ztmp) - IF( iom_use('Ce_ice') ) CALL iom_put("Ce_ice", Ce_ice*ztmp) - IF( iom_use('Ch_ice') ) CALL iom_put("Ch_ice", Ch_ice*ztmp) - IF( ln_blk ) THEN ! ---------------------------------------------------- ! ! Wind stress relative to nonmoving ice ( U10m ) ! ! ---------------------------------------------------- ! ! supress moving ice in wind stress computation as we don't know how to do it properly... - DO_2D( 0, 1, 0, 1 ) ! at T point - zztmp1 = rhoa(ji,jj) * Cd_ice(ji,jj) * wndm_ice(ji,jj) - putaui(ji,jj) = zztmp1 * pwndi(ji,jj) - pvtaui(ji,jj) = zztmp1 * pwndj(ji,jj) + DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) + zztmp = rhoa(ji,jj) * Cd_ice(ji,jj) * wndm_ice(ji,jj) + putaui(ji,jj) = zztmp * pwndi(ji,jj) + pvtaui(ji,jj) = zztmp * pwndj(ji,jj) END_2D - !#LB: saving the module, and x-y components, of the ai wind-stress at T-points: NOT weighted by the ice concentration !!! - IF(iom_use('taum_ice')) CALL iom_put('taum_ice', SQRT( putaui*putaui + pvtaui*pvtaui )*ztmp ) - !#LB: These 2 lines below mostly here for 'STATION_ASF' test-case, otherwize "utau_oi" (U-grid) and vtau_oi" (V-grid) does the job in: [ICE/icedyn_rhg_evp.F90]) - IF(iom_use('utau_ice')) CALL iom_put("utau_ice", putaui*ztmp) ! utau at T-points! - IF(iom_use('vtau_ice')) CALL iom_put("vtau_ice", pvtaui*ztmp) ! vtau at T-points! - - ! - DO_2D( 0, 0, 0, 0 ) ! U & V-points (same as ocean). - !#LB: QUESTION?? so SI3 expects wind stress vector to be provided at U & V points? Not at T-points ? - ! take care of the land-sea mask to avoid "pollution" of coastal stress. p[uv]taui used in frazil and rheology - zztmp1 = 0.5_wp * ( 2. - umask(ji,jj,1) ) * MAX( tmask(ji,jj,1),tmask(ji+1,jj ,1) ) - zztmp2 = 0.5_wp * ( 2. - vmask(ji,jj,1) ) * MAX( tmask(ji,jj,1),tmask(ji ,jj+1,1) ) - putaui(ji,jj) = zztmp1 * ( putaui(ji,jj) + putaui(ji+1,jj ) ) - pvtaui(ji,jj) = zztmp2 * ( pvtaui(ji,jj) + pvtaui(ji ,jj+1) ) - END_2D - CALL lbc_lnk( 'sbcblk', putaui, 'U', -1._wp, pvtaui, 'V', -1._wp ) + ! outputs + ! LB: not weighted by the ice concentration + IF( iom_use('taum_ice') ) CALL iom_put( 'taum_ice', SQRT( putaui*putaui + pvtaui*pvtaui ) * zmsk00 ) + ! LB: These 2 lines below mostly here for 'STATION_ASF' test-case + IF( iom_use('utau_ice') ) CALL iom_put( "utau_ice", putaui * zmsk00 ) + IF( iom_use('vtau_ice') ) CALL iom_put( "vtau_ice", pvtaui * zmsk00 ) ! - IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab2d_1=putaui , clinfo1=' blk_ice: putaui : ', mask1=umask & - & , tab2d_2=pvtaui , clinfo2=' pvtaui : ', mask2=vmask ) + IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab2d_1=putaui , clinfo1=' blk_ice: putaui : ', mask1=tmask & + & , tab2d_2=pvtaui , clinfo2=' pvtaui : ', mask2=tmask ) ELSE ! ln_abl DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) @@ -1105,10 +1068,15 @@ CONTAINS pseni (ji,jj) = wndm_ice(ji,jj) * Ch_ice(ji,jj) pevpi (ji,jj) = wndm_ice(ji,jj) * Ce_ice(ji,jj) END_2D - pssqi(:,:) = q_sat( ptsui(:,:), pslp(:,:), l_ice=.TRUE. ) ; ! more accurate way to obtain ssq ! + pssqi(:,:) = q_sat( ptsui(:,:), pslp(:,:), l_ice=.TRUE. ) ! more accurate way to obtain ssq ENDIF ! ln_blk / ln_abl ! + ! outputs + IF( iom_use('Cd_ice') ) CALL iom_put( "Cd_ice", Cd_ice * zmsk00 ) + IF( iom_use('Ce_ice') ) CALL iom_put( "Ce_ice", Ce_ice * zmsk00 ) + IF( iom_use('Ch_ice') ) CALL iom_put( "Ch_ice", Ch_ice * zmsk00 ) + ! IF(sn_cfctl%l_prtctl) CALL prt_ctl(tab2d_1=wndm_ice , clinfo1=' blk_ice: wndm_ice : ', mask1=tmask ) ! END SUBROUTINE blk_ice_1 diff --git a/src/OCE/SBC/sbccpl.F90 b/src/OCE/SBC/sbccpl.F90 index eddf626a5..c3feee0dc 100644 --- a/src/OCE/SBC/sbccpl.F90 +++ b/src/OCE/SBC/sbccpl.F90 @@ -68,15 +68,15 @@ MODULE sbccpl INTEGER, PARAMETER :: jpr_otx1 = 1 ! 3 atmosphere-ocean stress components on grid 1 INTEGER, PARAMETER :: jpr_oty1 = 2 ! INTEGER, PARAMETER :: jpr_otz1 = 3 ! - INTEGER, PARAMETER :: jpr_otx2 = 4 ! 3 atmosphere-ocean stress components on grid 2 - INTEGER, PARAMETER :: jpr_oty2 = 5 ! - INTEGER, PARAMETER :: jpr_otz2 = 6 ! +!!$ INTEGER, PARAMETER :: jpr_otx2 = 4 ! 3 atmosphere-ocean stress components on grid 2 +!!$ INTEGER, PARAMETER :: jpr_oty2 = 5 ! +!!$ INTEGER, PARAMETER :: jpr_otz2 = 6 ! INTEGER, PARAMETER :: jpr_itx1 = 7 ! 3 atmosphere-ice stress components on grid 1 INTEGER, PARAMETER :: jpr_ity1 = 8 ! INTEGER, PARAMETER :: jpr_itz1 = 9 ! - INTEGER, PARAMETER :: jpr_itx2 = 10 ! 3 atmosphere-ice stress components on grid 2 - INTEGER, PARAMETER :: jpr_ity2 = 11 ! - INTEGER, PARAMETER :: jpr_itz2 = 12 ! +!!$ INTEGER, PARAMETER :: jpr_itx2 = 10 ! 3 atmosphere-ice stress components on grid 2 +!!$ INTEGER, PARAMETER :: jpr_ity2 = 11 ! +!!$ INTEGER, PARAMETER :: jpr_itz2 = 12 ! INTEGER, PARAMETER :: jpr_qsroce = 13 ! Qsr above the ocean INTEGER, PARAMETER :: jpr_qsrice = 14 ! Qsr above the ice INTEGER, PARAMETER :: jpr_qsrmix = 15 @@ -128,9 +128,9 @@ MODULE sbccpl INTEGER, PARAMETER :: jpr_isf = 60 INTEGER, PARAMETER :: jpr_icb = 61 INTEGER, PARAMETER :: jpr_ts_ice = 62 ! Sea ice surface temp - !!INTEGER, PARAMETER :: jpr_qtrice = 63 ! Transmitted solar thru sea-ice + INTEGER, PARAMETER :: jpr_qtrice = 63 ! Transmitted solar thru sea-ice - INTEGER, PARAMETER :: jprcv = 62 ! total number of fields received + INTEGER, PARAMETER :: jprcv = 63 ! total number of fields received INTEGER, PARAMETER :: jps_fice = 1 ! ice fraction sent to the atmosphere INTEGER, PARAMETER :: jps_toce = 2 ! ocean temperature @@ -194,7 +194,7 @@ MODULE sbccpl & sn_snd_thick1, sn_snd_cond, sn_snd_mpnd , sn_snd_sstfrz, sn_snd_ttilyr ! ! Received from the atmosphere TYPE(FLD_C) :: sn_rcv_w10m, sn_rcv_taumod, sn_rcv_tau, sn_rcv_dqnsdt, sn_rcv_qsr, & - & sn_rcv_qns , sn_rcv_emp , sn_rcv_rnf, sn_rcv_ts_ice + & sn_rcv_qns , sn_rcv_emp , sn_rcv_rnf, sn_rcv_ts_ice, sn_rcv_qtrice TYPE(FLD_C) :: sn_rcv_cal, sn_rcv_iceflx, sn_rcv_co2, sn_rcv_mslp, sn_rcv_icb, sn_rcv_isf ! ! Send to waves TYPE(FLD_C) :: sn_snd_ifrac, sn_snd_crtw, sn_snd_wlev @@ -202,7 +202,6 @@ MODULE sbccpl TYPE(FLD_C) :: sn_rcv_hsig, sn_rcv_phioc, sn_rcv_sdrfx, sn_rcv_sdrfy, sn_rcv_wper, sn_rcv_wnum, & & sn_rcv_wstrf, sn_rcv_wdrag, sn_rcv_charn, sn_rcv_taw, sn_rcv_bhd, sn_rcv_tusd, sn_rcv_tvsd ! ! Other namelist parameters -!! TYPE(FLD_C) :: sn_rcv_qtrice INTEGER :: nn_cplmodel ! Maximum number of models to/from which NEMO is potentialy sending/receiving data LOGICAL :: ln_usecplmask ! use a coupling mask file to merge data received from several models ! -> file cplmask.nc with the float variable called cplmask (jpi,jpj,nn_cplmodel) @@ -281,7 +280,7 @@ CONTAINS & sn_rcv_sdrfx , sn_rcv_sdrfy , sn_rcv_wper , sn_rcv_wnum , sn_rcv_wstrf , & & sn_rcv_charn , sn_rcv_taw , sn_rcv_bhd , sn_rcv_tusd , sn_rcv_tvsd, & & sn_rcv_wdrag , sn_rcv_qns , sn_rcv_emp , sn_rcv_rnf , sn_rcv_cal , & - & sn_rcv_iceflx, sn_rcv_co2 , sn_rcv_icb , sn_rcv_isf , sn_rcv_ts_ice, & !!, sn_rcv_qtrice + & sn_rcv_iceflx, sn_rcv_co2 , sn_rcv_icb , sn_rcv_isf , sn_rcv_ts_ice, sn_rcv_qtrice, & & sn_rcv_mslp !!--------------------------------------------------------------------- @@ -313,7 +312,6 @@ CONTAINS WRITE(numout,*)' surface stress = ', TRIM(sn_rcv_tau%cldes ), ' (', TRIM(sn_rcv_tau%clcat ), ')' WRITE(numout,*)' - referential = ', sn_rcv_tau%clvref WRITE(numout,*)' - orientation = ', sn_rcv_tau%clvor - WRITE(numout,*)' - mesh = ', sn_rcv_tau%clvgrd WRITE(numout,*)' non-solar heat flux sensitivity = ', TRIM(sn_rcv_dqnsdt%cldes), ' (', TRIM(sn_rcv_dqnsdt%clcat), ')' WRITE(numout,*)' solar heat flux = ', TRIM(sn_rcv_qsr%cldes ), ' (', TRIM(sn_rcv_qsr%clcat ), ')' WRITE(numout,*)' non-solar heat flux = ', TRIM(sn_rcv_qns%cldes ), ' (', TRIM(sn_rcv_qns%clcat ), ')' @@ -323,7 +321,7 @@ CONTAINS WRITE(numout,*)' iceberg = ', TRIM(sn_rcv_icb%cldes ), ' (', TRIM(sn_rcv_icb%clcat ), ')' WRITE(numout,*)' ice shelf = ', TRIM(sn_rcv_isf%cldes ), ' (', TRIM(sn_rcv_isf%clcat ), ')' WRITE(numout,*)' sea ice heat fluxes = ', TRIM(sn_rcv_iceflx%cldes), ' (', TRIM(sn_rcv_iceflx%clcat), ')' -!! WRITE(numout,*)' transmitted solar thru sea-ice = ', TRIM(sn_rcv_qtrice%cldes), ' (', TRIM(sn_rcv_qtrice%clcat), ')' + WRITE(numout,*)' transmitted solar thru sea-ice = ', TRIM(sn_rcv_qtrice%cldes), ' (', TRIM(sn_rcv_qtrice%clcat), ')' WRITE(numout,*)' atm co2 = ', TRIM(sn_rcv_co2%cldes ), ' (', TRIM(sn_rcv_co2%clcat ), ')' WRITE(numout,*)' Sea ice surface skin temperature= ', TRIM(sn_rcv_ts_ice%cldes), ' (', TRIM(sn_rcv_ts_ice%clcat), ')' WRITE(numout,*)' surface waves:' @@ -358,7 +356,7 @@ CONTAINS WRITE(numout,*)' - mesh = ', sn_snd_crtw%clvgrd ENDIF IF( lwp .AND. ln_wave) THEN ! control print - WRITE(numout,*)' surface waves:' + WRITE(numout,*)' surface waves:' WRITE(numout,*)' Significant wave heigth = ', TRIM(sn_rcv_hsig%cldes ), ' (', TRIM(sn_rcv_hsig%clcat ), ')' WRITE(numout,*)' Wave to oce energy flux = ', TRIM(sn_rcv_phioc%cldes ), ' (', TRIM(sn_rcv_phioc%clcat ), ')' WRITE(numout,*)' Surface Stokes drift grid u = ', TRIM(sn_rcv_sdrfx%cldes ), ' (', TRIM(sn_rcv_sdrfx%clcat ), ')' @@ -368,8 +366,8 @@ CONTAINS WRITE(numout,*)' Stress frac adsorbed by waves = ', TRIM(sn_rcv_wstrf%cldes ), ' (', TRIM(sn_rcv_wstrf%clcat ), ')' WRITE(numout,*)' Neutral surf drag coefficient = ', TRIM(sn_rcv_wdrag%cldes ), ' (', TRIM(sn_rcv_wdrag%clcat ), ')' WRITE(numout,*)' Charnock coefficient = ', TRIM(sn_rcv_charn%cldes ), ' (', TRIM(sn_rcv_charn%clcat ), ')' - WRITE(numout,*)' Transport associated to Stokes drift grid u = ', TRIM(sn_rcv_tusd%cldes ), ' (', TRIM(sn_rcv_tusd%clcat ), ')' - WRITE(numout,*)' Transport associated to Stokes drift grid v = ', TRIM(sn_rcv_tvsd%cldes ), ' (', TRIM(sn_rcv_tvsd%clcat ), ')' + WRITE(numout,*)' Transport associated to Stokes drift u = ', TRIM(sn_rcv_tusd%cldes ), ' (', TRIM(sn_rcv_tusd%clcat ), ')' + WRITE(numout,*)' Transport associated to Stokes drift v = ', TRIM(sn_rcv_tvsd%cldes ), ' (', TRIM(sn_rcv_tvsd%clcat ), ')' WRITE(numout,*)' Bernouilli pressure head = ', TRIM(sn_rcv_bhd%cldes ), ' (', TRIM(sn_rcv_bhd%clcat ), ')' WRITE(numout,*)'Wave to ocean momentum flux and Net wave-supported stress = ', TRIM(sn_rcv_taw%cldes ), ' (', TRIM(sn_rcv_taw%clcat ), ')' WRITE(numout,*)' Surface current to waves = ', TRIM(sn_snd_crtw%cldes ), ' (', TRIM(sn_snd_crtw%clcat ), ')' @@ -399,87 +397,26 @@ CONTAINS srcv(jpr_otx1)%clname = 'O_OTaux1' ! 1st ocean component on grid ONE (T or U) srcv(jpr_oty1)%clname = 'O_OTauy1' ! 2nd - - - - srcv(jpr_otz1)%clname = 'O_OTauz1' ! 3rd - - - - - srcv(jpr_otx2)%clname = 'O_OTaux2' ! 1st ocean component on grid TWO (V) - srcv(jpr_oty2)%clname = 'O_OTauy2' ! 2nd - - - - - srcv(jpr_otz2)%clname = 'O_OTauz2' ! 3rd - - - - ! srcv(jpr_itx1)%clname = 'O_ITaux1' ! 1st ice component on grid ONE (T, F, I or U) srcv(jpr_ity1)%clname = 'O_ITauy1' ! 2nd - - - - srcv(jpr_itz1)%clname = 'O_ITauz1' ! 3rd - - - - - srcv(jpr_itx2)%clname = 'O_ITaux2' ! 1st ice component on grid TWO (V) - srcv(jpr_ity2)%clname = 'O_ITauy2' ! 2nd - - - - - srcv(jpr_itz2)%clname = 'O_ITauz2' ! 3rd - - - - ! ! Vectors: change of sign at north fold ONLY if on the local grid IF( TRIM( sn_rcv_tau%cldes ) == 'oce only' .OR. TRIM( sn_rcv_tau%cldes ) == 'oce and ice' & .OR. TRIM( sn_rcv_tau%cldes ) == 'mixed oce-ice' ) THEN ! avoid working with the atmospheric fields if they are not coupled ! - IF( TRIM( sn_rcv_tau%clvor ) == 'local grid' ) srcv(jpr_otx1:jpr_itz2)%nsgn = -1. - - ! ! Set grid and action - SELECT CASE( TRIM( sn_rcv_tau%clvgrd ) ) ! 'T', 'U,V', 'U,V,I', 'U,V,F', 'T,I', 'T,F', or 'T,U,V' - CASE( 'T' ) - srcv(jpr_otx1:jpr_itz2)%clgrid = 'T' ! oce and ice components given at T-point - srcv(jpr_otx1:jpr_otz1)%laction = .TRUE. ! receive oce components on grid 1 - srcv(jpr_itx1:jpr_itz1)%laction = .TRUE. ! receive ice components on grid 1 - CASE( 'U,V' ) - srcv(jpr_otx1:jpr_otz1)%clgrid = 'U' ! oce components given at U-point - srcv(jpr_otx2:jpr_otz2)%clgrid = 'V' ! and V-point - srcv(jpr_itx1:jpr_itz1)%clgrid = 'U' ! ice components given at U-point - srcv(jpr_itx2:jpr_itz2)%clgrid = 'V' ! and V-point - srcv(jpr_otx1:jpr_itz2)%laction = .TRUE. ! receive oce and ice components on both grid 1 & 2 - CASE( 'U,V,T' ) - srcv(jpr_otx1:jpr_otz1)%clgrid = 'U' ! oce components given at U-point - srcv(jpr_otx2:jpr_otz2)%clgrid = 'V' ! and V-point - srcv(jpr_itx1:jpr_itz1)%clgrid = 'T' ! ice components given at T-point - srcv(jpr_otx1:jpr_otz2)%laction = .TRUE. ! receive oce components on grid 1 & 2 - srcv(jpr_itx1:jpr_itz1)%laction = .TRUE. ! receive ice components on grid 1 only - CASE( 'U,V,I' ) - srcv(jpr_otx1:jpr_otz1)%clgrid = 'U' ! oce components given at U-point - srcv(jpr_otx2:jpr_otz2)%clgrid = 'V' ! and V-point - srcv(jpr_itx1:jpr_itz1)%clgrid = 'I' ! ice components given at I-point - srcv(jpr_otx1:jpr_otz2)%laction = .TRUE. ! receive oce components on grid 1 & 2 - srcv(jpr_itx1:jpr_itz1)%laction = .TRUE. ! receive ice components on grid 1 only - CASE( 'U,V,F' ) - srcv(jpr_otx1:jpr_otz1)%clgrid = 'U' ! oce components given at U-point - srcv(jpr_otx2:jpr_otz2)%clgrid = 'V' ! and V-point - srcv(jpr_itx1:jpr_itz1)%clgrid = 'F' ! ice components given at F-point - srcv(jpr_otx1:jpr_otz2)%laction = .TRUE. ! receive oce components on grid 1 & 2 - srcv(jpr_itx1:jpr_itz1)%laction = .TRUE. ! receive ice components on grid 1 only - CASE( 'T,I' ) - srcv(jpr_otx1:jpr_itz2)%clgrid = 'T' ! oce and ice components given at T-point - srcv(jpr_itx1:jpr_itz1)%clgrid = 'I' ! ice components given at I-point - srcv(jpr_otx1:jpr_otz1)%laction = .TRUE. ! receive oce components on grid 1 - srcv(jpr_itx1:jpr_itz1)%laction = .TRUE. ! receive ice components on grid 1 - CASE( 'T,F' ) - srcv(jpr_otx1:jpr_itz2)%clgrid = 'T' ! oce and ice components given at T-point - srcv(jpr_itx1:jpr_itz1)%clgrid = 'F' ! ice components given at F-point - srcv(jpr_otx1:jpr_otz1)%laction = .TRUE. ! receive oce components on grid 1 - srcv(jpr_itx1:jpr_itz1)%laction = .TRUE. ! receive ice components on grid 1 - CASE( 'T,U,V' ) - srcv(jpr_otx1:jpr_otz1)%clgrid = 'T' ! oce components given at T-point - srcv(jpr_itx1:jpr_itz1)%clgrid = 'U' ! ice components given at U-point - srcv(jpr_itx2:jpr_itz2)%clgrid = 'V' ! and V-point - srcv(jpr_otx1:jpr_otz1)%laction = .TRUE. ! receive oce components on grid 1 only - srcv(jpr_itx1:jpr_itz2)%laction = .TRUE. ! receive ice components on grid 1 & 2 - CASE default - CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_rcv_tau%clvgrd' ) - END SELECT + IF( TRIM( sn_rcv_tau%clvor ) == 'local grid' ) srcv(jpr_otx1:jpr_itz1)%nsgn = -1. + + ! ! Set grid and action + srcv(jpr_otx1:jpr_otz1)%laction = .TRUE. ! receive oce components on grid 1 + srcv(jpr_itx1:jpr_itz1)%laction = .TRUE. ! receive ice components on grid 1 ! IF( TRIM( sn_rcv_tau%clvref ) == 'spherical' ) & ! spherical: 3rd component not received - & srcv( (/jpr_otz1, jpr_otz2, jpr_itz1, jpr_itz2/) )%laction = .FALSE. - ! - IF( TRIM( sn_rcv_tau%clvor ) == 'local grid' ) THEN ! already on local grid -> no need of the second grid - srcv(jpr_otx2:jpr_otz2)%laction = .FALSE. - srcv(jpr_itx2:jpr_itz2)%laction = .FALSE. - srcv(jpr_oty1)%clgrid = srcv(jpr_oty2)%clgrid ! not needed but cleaner... - srcv(jpr_ity1)%clgrid = srcv(jpr_ity2)%clgrid ! not needed but cleaner... - ENDIF + & srcv( (/jpr_otz1, jpr_itz1/) )%laction = .FALSE. ! - IF( TRIM( sn_rcv_tau%cldes ) /= 'oce and ice' ) THEN ! 'oce and ice' case ocean stress on ocean mesh used - srcv(jpr_itx1:jpr_itz2)%laction = .FALSE. ! ice components not received - srcv(jpr_itx1)%clgrid = 'U' ! ocean stress used after its transformation - srcv(jpr_ity1)%clgrid = 'V' ! i.e. it is always at U- & V-points for i- & j-comp. resp. + IF( TRIM( sn_rcv_tau%cldes ) /= 'oce and ice' ) THEN ! 'oce and ice' case ocean stress on ocean mesh used + srcv(jpr_itx1:jpr_itz1)%laction = .FALSE. ! ice components not received ENDIF ENDIF @@ -612,18 +549,18 @@ CONTAINS ENDIF srcv(jpr_topm:jpr_botm)%laction = .TRUE. ENDIF -!! ! ! --------------------------- ! -!! ! ! transmitted solar thru ice ! -!! ! ! --------------------------- ! -!! srcv(jpr_qtrice)%clname = 'OQtr' -!! IF( TRIM(sn_rcv_qtrice%cldes) == 'coupled' ) THEN -!! IF ( TRIM( sn_rcv_qtrice%clcat ) == 'yes' ) THEN -!! srcv(jpr_qtrice)%nct = nn_cats_cpl -!! ELSE -!! CALL ctl_stop( 'sbc_cpl_init: sn_rcv_qtrice%clcat should always be set to yes currently' ) -!! ENDIF -!! srcv(jpr_qtrice)%laction = .TRUE. -!! ENDIF + ! ! --------------------------- ! + ! ! transmitted solar thru ice ! + ! ! --------------------------- ! + srcv(jpr_qtrice)%clname = 'OQtr' + IF( TRIM(sn_rcv_qtrice%cldes) == 'coupled' ) THEN + IF ( TRIM( sn_rcv_qtrice%clcat ) == 'yes' ) THEN + srcv(jpr_qtrice)%nct = nn_cats_cpl + ELSE + CALL ctl_stop( 'sbc_cpl_init: sn_rcv_qtrice%clcat should always be set to yes currently' ) + ENDIF + srcv(jpr_qtrice)%laction = .TRUE. + ENDIF ! ! ------------------------- ! ! ! ice skin temperature ! ! ! ------------------------- ! @@ -725,11 +662,10 @@ CONTAINS srcv(:)%clgrid = 'T' ! force default definition in case of opa <-> sas coupling srcv(:)%nsgn = 1. ! force default definition in case of opa <-> sas coupling srcv( (/jpr_qsroce, jpr_qnsoce, jpr_oemp, jpr_sflx, jpr_fice, jpr_otx1, jpr_oty1, jpr_taum/) )%laction = .TRUE. - srcv(jpr_otx1)%clgrid = 'U' ! oce components given at U-point - srcv(jpr_oty1)%clgrid = 'V' ! and V-point + srcv(jpr_otx1)%clgrid = 'T' ! oce components given at T-point + srcv(jpr_oty1)%clgrid = 'T' ! Vectors: change of sign at north fold ONLY if on the local grid srcv( (/jpr_otx1,jpr_oty1/) )%nsgn = -1. - sn_rcv_tau%clvgrd = 'U,V' sn_rcv_tau%clvor = 'local grid' sn_rcv_tau%clvref = 'spherical' sn_rcv_emp%cldes = 'oce only' @@ -1162,7 +1098,7 @@ CONTAINS !! ** Method : receive all fields from the atmosphere and transform !! them into ocean surface boundary condition fields !! - !! ** Action : update utau, vtau ocean stress at U,V grid + !! ** Action : update utau, vtau ocean stress at T-point !! taum wind stress module at T-point !! wndm wind speed module at T-point over free ocean or leads in presence of sea-ice !! qns non solar heat fluxes including emp heat content (ocean only case) @@ -1221,39 +1157,20 @@ CONTAINS IF( TRIM( sn_rcv_tau%clvref ) == 'cartesian' ) THEN ! 2 components on the sphere ! ! (cartesian to spherical -> 3 to 2 components) ! - CALL geo2oce( frcv(jpr_otx1)%z3(:,:,1), frcv(jpr_oty1)%z3(:,:,1), frcv(jpr_otz1)%z3(:,:,1), & - & srcv(jpr_otx1)%clgrid, ztx, zty ) + CALL geo2oce( frcv(jpr_otx1)%z3(:,:,1), frcv(jpr_oty1)%z3(:,:,1), frcv(jpr_otz1)%z3(:,:,1), 'T', ztx, zty ) frcv(jpr_otx1)%z3(:,:,1) = ztx(:,:) ! overwrite 1st comp. on the 1st grid frcv(jpr_oty1)%z3(:,:,1) = zty(:,:) ! overwrite 2nd comp. on the 1st grid ! - IF( srcv(jpr_otx2)%laction ) THEN - CALL geo2oce( frcv(jpr_otx2)%z3(:,:,1), frcv(jpr_oty2)%z3(:,:,1), frcv(jpr_otz2)%z3(:,:,1), & - & srcv(jpr_otx2)%clgrid, ztx, zty ) - frcv(jpr_otx2)%z3(:,:,1) = ztx(:,:) ! overwrite 1st comp. on the 2nd grid - frcv(jpr_oty2)%z3(:,:,1) = zty(:,:) ! overwrite 2nd comp. on the 2nd grid - ENDIF - ! ENDIF ! IF( TRIM( sn_rcv_tau%clvor ) == 'eastward-northward' ) THEN ! 2 components oriented along the local grid ! ! (geographical to local grid -> rotate the components) - CALL rot_rep( frcv(jpr_otx1)%z3(:,:,1), frcv(jpr_oty1)%z3(:,:,1), srcv(jpr_otx1)%clgrid, 'en->i', ztx ) - IF( srcv(jpr_otx2)%laction ) THEN - CALL rot_rep( frcv(jpr_otx2)%z3(:,:,1), frcv(jpr_oty2)%z3(:,:,1), srcv(jpr_otx2)%clgrid, 'en->j', zty ) - ELSE - CALL rot_rep( frcv(jpr_otx1)%z3(:,:,1), frcv(jpr_oty1)%z3(:,:,1), srcv(jpr_otx1)%clgrid, 'en->j', zty ) - ENDIF + CALL rot_rep( frcv(jpr_otx1)%z3(:,:,1), frcv(jpr_oty1)%z3(:,:,1), 'T', 'en->i', ztx ) + CALL rot_rep( frcv(jpr_otx1)%z3(:,:,1), frcv(jpr_oty1)%z3(:,:,1), 'T', 'en->j', zty ) frcv(jpr_otx1)%z3(:,:,1) = ztx(:,:) ! overwrite 1st component on the 1st grid frcv(jpr_oty1)%z3(:,:,1) = zty(:,:) ! overwrite 2nd component on the 2nd grid ENDIF ! - IF( srcv(jpr_otx1)%clgrid == 'T' ) THEN - DO_2D( 0, 0, 0, 0 ) ! T ==> (U,V) - frcv(jpr_otx1)%z3(ji,jj,1) = 0.5 * ( frcv(jpr_otx1)%z3(ji+1,jj ,1) + frcv(jpr_otx1)%z3(ji,jj,1) ) - frcv(jpr_oty1)%z3(ji,jj,1) = 0.5 * ( frcv(jpr_oty1)%z3(ji ,jj+1,1) + frcv(jpr_oty1)%z3(ji,jj,1) ) - END_2D - CALL lbc_lnk( 'sbccpl', frcv(jpr_otx1)%z3(:,:,1), 'U', -1.0_wp, frcv(jpr_oty1)%z3(:,:,1), 'V', -1.0_wp ) - ENDIF llnewtx = .TRUE. ELSE llnewtx = .FALSE. @@ -1272,12 +1189,9 @@ CONTAINS IF( .NOT. srcv(jpr_taum)%laction ) THEN ! compute wind stress module from its components if not received ! => need to be done only when otx1 was changed IF( llnewtx ) THEN - DO_2D( 0, 0, 0, 0 ) - zzx = frcv(jpr_otx1)%z3(ji-1,jj ,1) + frcv(jpr_otx1)%z3(ji,jj,1) - zzy = frcv(jpr_oty1)%z3(ji ,jj-1,1) + frcv(jpr_oty1)%z3(ji,jj,1) - frcv(jpr_taum)%z3(ji,jj,1) = 0.5 * SQRT( zzx * zzx + zzy * zzy ) - END_2D - CALL lbc_lnk( 'sbccpl', frcv(jpr_taum)%z3(:,:,1), 'T', 1.0_wp ) + zzx = frcv(jpr_otx1)%z3(ji,jj,1) + zzy = frcv(jpr_oty1)%z3(ji,jj,1) + frcv(jpr_taum)%z3(ji,jj,1) = 0.5 * SQRT( zzx * zzx + zzy * zzy ) llnewtau = .TRUE. ELSE llnewtau = .FALSE. @@ -1594,7 +1508,7 @@ CONTAINS !! ** Action : return ptau_i, ptau_j, the stress over the ice !!---------------------------------------------------------------------- REAL(wp), INTENT(inout), DIMENSION(:,:) :: p_taui ! i- & j-components of atmos-ice stress [N/m2] - REAL(wp), INTENT(inout), DIMENSION(:,:) :: p_tauj ! at I-point (B-grid) or U & V-point (C-grid) + REAL(wp), INTENT(inout), DIMENSION(:,:) :: p_tauj ! at T-point !! INTEGER :: ji, jj ! dummy loop indices INTEGER :: itx ! index of taux over ice @@ -1616,28 +1530,16 @@ CONTAINS ! IF( TRIM( sn_rcv_tau%clvref ) == 'cartesian' ) THEN ! 2 components on the sphere ! ! (cartesian to spherical -> 3 to 2 components) - CALL geo2oce( frcv(jpr_itx1)%z3(:,:,1), frcv(jpr_ity1)%z3(:,:,1), frcv(jpr_itz1)%z3(:,:,1), & - & srcv(jpr_itx1)%clgrid, ztx, zty ) + CALL geo2oce( frcv(jpr_itx1)%z3(:,:,1), frcv(jpr_ity1)%z3(:,:,1), frcv(jpr_itz1)%z3(:,:,1), 'T', ztx, zty ) frcv(jpr_itx1)%z3(:,:,1) = ztx(:,:) ! overwrite 1st comp. on the 1st grid frcv(jpr_ity1)%z3(:,:,1) = zty(:,:) ! overwrite 2nd comp. on the 1st grid ! - IF( srcv(jpr_itx2)%laction ) THEN - CALL geo2oce( frcv(jpr_itx2)%z3(:,:,1), frcv(jpr_ity2)%z3(:,:,1), frcv(jpr_itz2)%z3(:,:,1), & - & srcv(jpr_itx2)%clgrid, ztx, zty ) - frcv(jpr_itx2)%z3(:,:,1) = ztx(:,:) ! overwrite 1st comp. on the 2nd grid - frcv(jpr_ity2)%z3(:,:,1) = zty(:,:) ! overwrite 2nd comp. on the 2nd grid - ENDIF - ! ENDIF ! IF( TRIM( sn_rcv_tau%clvor ) == 'eastward-northward' ) THEN ! 2 components oriented along the local grid ! ! (geographical to local grid -> rotate the components) - CALL rot_rep( frcv(jpr_itx1)%z3(:,:,1), frcv(jpr_ity1)%z3(:,:,1), srcv(jpr_itx1)%clgrid, 'en->i', ztx ) - IF( srcv(jpr_itx2)%laction ) THEN - CALL rot_rep( frcv(jpr_itx2)%z3(:,:,1), frcv(jpr_ity2)%z3(:,:,1), srcv(jpr_itx2)%clgrid, 'en->j', zty ) - ELSE - CALL rot_rep( frcv(jpr_itx1)%z3(:,:,1), frcv(jpr_ity1)%z3(:,:,1), srcv(jpr_itx1)%clgrid, 'en->j', zty ) - ENDIF + CALL rot_rep( frcv(jpr_itx1)%z3(:,:,1), frcv(jpr_ity1)%z3(:,:,1), 'T', 'en->i', ztx ) + CALL rot_rep( frcv(jpr_itx1)%z3(:,:,1), frcv(jpr_ity1)%z3(:,:,1), 'T', 'en->j', zty ) frcv(jpr_itx1)%z3(:,:,1) = ztx(:,:) ! overwrite 1st component on the 1st grid frcv(jpr_ity1)%z3(:,:,1) = zty(:,:) ! overwrite 2nd component on the 1st grid ENDIF @@ -1651,29 +1553,8 @@ CONTAINS ! ! ======================= ! ! ! put on ice grid ! ! ! ======================= ! - ! - ! j+1 j -----V---F - ! ice stress on ice velocity point ! | - ! (C-grid ==>(U,V)) j | T U - ! | | - ! j j-1 -I-------| - ! (for I) | | - ! i-1 i i - ! i i+1 (for I) - SELECT CASE ( srcv(jpr_itx1)%clgrid ) - CASE( 'U' ) - p_taui(:,:) = frcv(jpr_itx1)%z3(:,:,1) ! (U,V) ==> (U,V) - p_tauj(:,:) = frcv(jpr_ity1)%z3(:,:,1) - CASE( 'T' ) - DO_2D( 0, 0, 0, 0 ) ! T ==> (U,V) - ! take care of the land-sea mask to avoid "pollution" of coastal stress. p[uv]taui used in frazil and rheology - zztmp1 = 0.5_wp * ( 2. - umask(ji,jj,1) ) * MAX( tmask(ji,jj,1),tmask(ji+1,jj ,1) ) - zztmp2 = 0.5_wp * ( 2. - vmask(ji,jj,1) ) * MAX( tmask(ji,jj,1),tmask(ji ,jj+1,1) ) - p_taui(ji,jj) = zztmp1 * ( frcv(jpr_itx1)%z3(ji+1,jj ,1) + frcv(jpr_itx1)%z3(ji,jj,1) ) - p_tauj(ji,jj) = zztmp2 * ( frcv(jpr_ity1)%z3(ji ,jj+1,1) + frcv(jpr_ity1)%z3(ji,jj,1) ) - END_2D - CALL lbc_lnk( 'sbccpl', p_taui, 'U', -1., p_tauj, 'V', -1. ) - END SELECT + p_taui(:,:) = frcv(jpr_itx1)%z3(:,:,1) + p_tauj(:,:) = frcv(jpr_ity1)%z3(:,:,1) ENDIF ! @@ -1794,8 +1675,8 @@ CONTAINS END SELECT ! --- evaporation over ice (kg/m2/s) --- ! - IF (ln_scale_ice_flux) THEN ! typically met-office requirements - IF (sn_rcv_emp%clcat == 'yes') THEN + IF( ln_scale_ice_flux ) THEN ! typically met-office requirements + IF( sn_rcv_emp%clcat == 'yes' ) THEN WHERE( a_i(:,:,:) > 1.e-10 ) ; zevap_ice(:,:,:) = frcv(jpr_ievp)%z3(:,:,:) * a_i_last_couple(:,:,:) / a_i(:,:,:) ELSEWHERE ; zevap_ice(:,:,:) = 0._wp END WHERE @@ -1812,7 +1693,7 @@ CONTAINS ENDDO ENDIF ELSE - IF (sn_rcv_emp%clcat == 'yes') THEN + IF( sn_rcv_emp%clcat == 'yes' ) THEN zevap_ice(:,:,1:jpl) = frcv(jpr_ievp)%z3(:,:,1:jpl) WHERE( picefr(:,:) > 1.e-10 ) ; zevap_ice_total(:,:) = SUM( zevap_ice(:,:,:) * a_i(:,:,:), dim=3 ) / picefr(:,:) ELSEWHERE ; zevap_ice_total(:,:) = 0._wp @@ -1826,7 +1707,7 @@ CONTAINS ENDIF ENDIF - IF ( TRIM( sn_rcv_emp%cldes ) == 'conservative' ) THEN + IF( TRIM( sn_rcv_emp%cldes ) == 'conservative' ) THEN ! For conservative case zemp_ice has not been defined yet. Do it now. zemp_ice(:,:) = zevap_ice_total(:,:) * picefr(:,:) - frcv(jpr_snow)%z3(:,:,1) * picefr(:,:) ENDIF @@ -1926,13 +1807,13 @@ CONTAINS & - zevap_ice_total(:,:) * picefr(:,:) ) * tmask(:,:,1) ) ! ice-free oce evap (cell average) ! note: runoff output is done in sbcrnf (which includes icebergs too) and iceshelf output is done in sbcisf !! IF( srcv(jpr_rnf)%laction ) CALL iom_put( 'runoffs' , rnf(:,:) * tmask(:,:,1) ) ! runoff -!! IF( srcv(jpr_isf)%laction ) CALL iom_put( 'iceshelf_cea', fwfisf(:,:) * tmask(:,:,1) ) ! iceshelf +!! IF( srcv(jpr_isf)%laction ) CALL iom_put( 'iceshelf_cea', fwfisf(:,:) * tmask(:,:,1) ) ! iceshelf ! ! ! ========================= ! SELECT CASE( TRIM( sn_rcv_iceflx%cldes ) ) ! ice topmelt and botmelt ! ! ! ========================= ! - CASE ('coupled') - IF (ln_scale_ice_flux) THEN + CASE( 'coupled' ) + IF( ln_scale_ice_flux ) THEN WHERE( a_i(:,:,:) > 1.e-10_wp ) qml_ice(:,:,:) = frcv(jpr_topm)%z3(:,:,:) * a_i_last_couple(:,:,:) / a_i(:,:,:) qcn_ice(:,:,:) = frcv(jpr_botm)%z3(:,:,:) * a_i_last_couple(:,:,:) / a_i(:,:,:) @@ -2213,28 +2094,29 @@ CONTAINS ! ELSEIF( ln_cndflx .AND. .NOT.ln_cndemulate ) THEN !== conduction flux as surface forcing ==! ! -!! SELECT CASE( TRIM( sn_rcv_qtrice%cldes ) ) -!! ! -!! ! ! ===> here we receive the qtr_ice_top array from the coupler -!! CASE ('coupled') -!! IF (ln_scale_ice_flux) THEN -!! WHERE( a_i(:,:,:) > 1.e-10_wp ) -!! zqtr_ice_top(:,:,:) = frcv(jpr_qtrice)%z3(:,:,:) * a_i_last_couple(:,:,:) / a_i(:,:,:) -!! ELSEWHERE -!! zqtr_ice_top(:,:,:) = 0.0_wp -!! ENDWHERE -!! ELSE -!! zqtr_ice_top(:,:,:) = frcv(jpr_qtrice)%z3(:,:,:) -!! ENDIF -!! -!! ! Add retrieved transmitted solar radiation onto the ice and total solar radiation -!! zqsr_ice(:,:,:) = zqsr_ice(:,:,:) + zqtr_ice_top(:,:,:) -!! zqsr_tot(:,:) = zqsr_tot(:,:) + SUM( zqtr_ice_top(:,:,:) * a_i(:,:,:), dim=3 ) -!! -!! ! if we are not getting this data from the coupler then assume zero (fully opaque ice) -!! CASE ('none') - zqtr_ice_top(:,:,:) = 0._wp -!! END SELECT + ! + SELECT CASE( TRIM( sn_rcv_qtrice%cldes ) ) + ! + ! ! ===> here we receive the qtr_ice_top array from the coupler + CASE ('coupled') + IF (ln_scale_ice_flux) THEN + WHERE( a_i(:,:,:) > 1.e-10_wp ) + zqtr_ice_top(:,:,:) = frcv(jpr_qtrice)%z3(:,:,:) * a_i_last_couple(:,:,:) / a_i(:,:,:) + ELSEWHERE + zqtr_ice_top(:,:,:) = 0.0_wp + ENDWHERE + ELSE + zqtr_ice_top(:,:,:) = frcv(jpr_qtrice)%z3(:,:,:) + ENDIF + + ! Add retrieved transmitted solar radiation onto the ice and total solar radiation + zqsr_ice(:,:,:) = zqsr_ice(:,:,:) + zqtr_ice_top(:,:,:) + zqsr_tot(:,:) = zqsr_tot(:,:) + SUM( zqtr_ice_top(:,:,:) * a_i(:,:,:), dim=3 ) + + ! if we are not getting this data from the coupler then assume zero (fully opaque ice) + CASE ('none') + zqtr_ice_top(:,:,:) = 0._wp + END SELECT ! ENDIF @@ -2403,10 +2285,10 @@ CONTAINS CASE( 'weighted ice' ) ; SELECT CASE( sn_snd_alb%clcat ) CASE( 'yes' ) - ztmp3(:,:,1:jpl) = alb_ice(:,:,1:jpl) * a_i(:,:,1:jpl) + ztmp3(:,:,1:jpl) = alb_ice(:,:,1:jpl) * a_i(:,:,1:jpl) CASE( 'no' ) WHERE( fr_i (:,:) > 0. ) - ztmp1(:,:) = SUM ( alb_ice(:,:,1:jpl) * a_i(:,:,1:jpl), dim=3 ) + ztmp1(:,:) = SUM ( alb_ice(:,:,1:jpl) * a_i(:,:,1:jpl), dim=3 ) ELSEWHERE ztmp1(:,:) = 0. END WHERE @@ -2572,17 +2454,18 @@ CONTAINS IF( ssnd(jps_ocx1)%laction ) THEN ! Surface current ! ! ! ------------------------- ! ! - ! j+1 j -----V---F - ! surface velocity always sent from T point ! | - ! j | T U - ! | | - ! j j-1 -I-------| - ! (for I) | | - ! i-1 i i - ! i i+1 (for I) + ! j -----V---F + ! surface velocity always sent from T point ! | + ! j | T U + ! | | + ! j-1 -I-------| + ! | | + ! i-1 i i + !!clem: make a new variable at T-point to replace uu and vv => uuT and vvT for instance IF( nn_components == jp_iam_oce ) THEN zotx1(:,:) = uu(:,:,1,Kmm) zoty1(:,:) = vv(:,:,1,Kmm) + !!clem : should be demi sum, no? Or uuT and vvT ELSE SELECT CASE( TRIM( sn_snd_crt%cldes ) ) CASE( 'oce only' ) ! C-grid ==> T @@ -2652,42 +2535,42 @@ CONTAINS ! ! Surface current to waves ! ! ! ------------------------- ! IF( ssnd(jps_ocxw)%laction .OR. ssnd(jps_ocyw)%laction ) THEN - ! - ! j+1 j -----V---F - ! surface velocity always sent from T point ! | - ! j | T U - ! | | - ! j j-1 -I-------| - ! (for I) | | - ! i-1 i i - ! i i+1 (for I) - SELECT CASE( TRIM( sn_snd_crtw%cldes ) ) - CASE( 'oce only' ) ! C-grid ==> T - DO_2D( 0, 0, 0, 0 ) - zotx1(ji,jj) = 0.5 * ( uu(ji,jj,1,Kmm) + uu(ji-1,jj ,1,Kmm) ) - zoty1(ji,jj) = 0.5 * ( vv(ji,jj,1,Kmm) + vv(ji , jj-1,1,Kmm) ) - END_2D - CASE( 'weighted oce and ice' ) ! Ocean and Ice on C-grid ==> T - DO_2D( 0, 0, 0, 0 ) - zotx1(ji,jj) = 0.5 * ( uu (ji,jj,1,Kmm) + uu (ji-1,jj ,1,Kmm) ) * zfr_l(ji,jj) - zoty1(ji,jj) = 0.5 * ( vv (ji,jj,1,Kmm) + vv (ji ,jj-1,1,Kmm) ) * zfr_l(ji,jj) - zitx1(ji,jj) = 0.5 * ( u_ice(ji,jj ) + u_ice(ji-1,jj ) ) * fr_i(ji,jj) - zity1(ji,jj) = 0.5 * ( v_ice(ji,jj ) + v_ice(ji ,jj-1 ) ) * fr_i(ji,jj) - END_2D - CALL lbc_lnk( 'sbccpl', zitx1, 'T', -1.0_wp, zity1, 'T', -1.0_wp ) - CASE( 'mixed oce-ice' ) ! Ocean and Ice on C-grid ==> T - DO_2D( 0, 0, 0, 0 ) - zotx1(ji,jj) = 0.5 * ( uu (ji,jj,1,Kmm) + uu (ji-1,jj ,1,Kmm) ) * zfr_l(ji,jj) & - & + 0.5 * ( u_ice(ji,jj ) + u_ice(ji-1,jj ) ) * fr_i(ji,jj) - zoty1(ji,jj) = 0.5 * ( vv (ji,jj,1,Kmm) + vv (ji ,jj-1,1,Kmm) ) * zfr_l(ji,jj) & - & + 0.5 * ( v_ice(ji,jj ) + v_ice(ji ,jj-1 ) ) * fr_i(ji,jj) - END_2D - END SELECT + ! + ! j -----V---F + ! surface velocity always sent from T point ! | + ! j | T U + ! | | + ! j-1 -I-------| + ! | | + ! i-1 i i + !!clem: make a new variable at T-point to replace uu and vv => uuT and vvT for instance + SELECT CASE( TRIM( sn_snd_crtw%cldes ) ) + CASE( 'oce only' ) ! C-grid ==> T + DO_2D( 0, 0, 0, 0 ) + zotx1(ji,jj) = 0.5 * ( uu(ji,jj,1,Kmm) + uu(ji-1,jj ,1,Kmm) ) + zoty1(ji,jj) = 0.5 * ( vv(ji,jj,1,Kmm) + vv(ji , jj-1,1,Kmm) ) + END_2D + CASE( 'weighted oce and ice' ) ! Ocean and Ice on C-grid ==> T + DO_2D( 0, 0, 0, 0 ) + zotx1(ji,jj) = 0.5 * ( uu (ji,jj,1,Kmm) + uu (ji-1,jj ,1,Kmm) ) * zfr_l(ji,jj) + zoty1(ji,jj) = 0.5 * ( vv (ji,jj,1,Kmm) + vv (ji ,jj-1,1,Kmm) ) * zfr_l(ji,jj) + zitx1(ji,jj) = 0.5 * ( u_ice(ji,jj ) + u_ice(ji-1,jj ) ) * fr_i(ji,jj) + zity1(ji,jj) = 0.5 * ( v_ice(ji,jj ) + v_ice(ji ,jj-1 ) ) * fr_i(ji,jj) + END_2D + CALL lbc_lnk( 'sbccpl', zitx1, 'T', -1.0_wp, zity1, 'T', -1.0_wp ) + CASE( 'mixed oce-ice' ) ! Ocean and Ice on C-grid ==> T + DO_2D( 0, 0, 0, 0 ) + zotx1(ji,jj) = 0.5 * ( uu (ji,jj,1,Kmm) + uu (ji-1,jj ,1,Kmm) ) * zfr_l(ji,jj) & + & + 0.5 * ( u_ice(ji,jj ) + u_ice(ji-1,jj ) ) * fr_i(ji,jj) + zoty1(ji,jj) = 0.5 * ( vv (ji,jj,1,Kmm) + vv (ji ,jj-1,1,Kmm) ) * zfr_l(ji,jj) & + & + 0.5 * ( v_ice(ji,jj ) + v_ice(ji ,jj-1 ) ) * fr_i(ji,jj) + END_2D + END SELECT CALL lbc_lnk( 'sbccpl', zotx1, ssnd(jps_ocxw)%clgrid, -1.0_wp, zoty1, ssnd(jps_ocyw)%clgrid, -1.0_wp ) ! ! IF( TRIM( sn_snd_crtw%clvor ) == 'eastward-northward' ) THEN ! Rotation of the components - ! ! Ocean component + ! ! Ocean component CALL rot_rep( zotx1, zoty1, ssnd(jps_ocxw)%clgrid, 'ij->e', ztmp1 ) ! 1st component CALL rot_rep( zotx1, zoty1, ssnd(jps_ocxw)%clgrid, 'ij->n', ztmp2 ) ! 2nd component zotx1(:,:) = ztmp1(:,:) ! overwrite the components @@ -2700,18 +2583,18 @@ CONTAINS ENDIF ENDIF ! -! ! spherical coordinates to cartesian -> 2 components to 3 components -! IF( TRIM( sn_snd_crtw%clvref ) == 'cartesian' ) THEN -! ztmp1(:,:) = zotx1(:,:) ! ocean currents -! ztmp2(:,:) = zoty1(:,:) -! CALL oce2geo ( ztmp1, ztmp2, 'T', zotx1, zoty1, zotz1 ) -! ! -! IF( ssnd(jps_ivx1)%laction ) THEN ! ice velocities -! ztmp1(:,:) = zitx1(:,:) -! ztmp1(:,:) = zity1(:,:) -! CALL oce2geo ( ztmp1, ztmp2, 'T', zitx1, zity1, zitz1 ) -! ENDIF -! ENDIF + ! ! spherical coordinates to cartesian -> 2 components to 3 components + ! IF( TRIM( sn_snd_crtw%clvref ) == 'cartesian' ) THEN + ! ztmp1(:,:) = zotx1(:,:) ! ocean currents + ! ztmp2(:,:) = zoty1(:,:) + ! CALL oce2geo ( ztmp1, ztmp2, 'T', zotx1, zoty1, zotz1 ) + ! ! + ! IF( ssnd(jps_ivx1)%laction ) THEN ! ice velocities + ! ztmp1(:,:) = zitx1(:,:) + ! ztmp1(:,:) = zity1(:,:) + ! CALL oce2geo ( ztmp1, ztmp2, 'T', zitx1, zity1, zitz1 ) + ! ENDIF + ! ENDIF ! IF( ssnd(jps_ocxw)%laction ) CALL cpl_snd( jps_ocxw, isec, RESHAPE ( zotx1, (/jpi,jpj,1/) ), info ) ! ocean x current 1st grid IF( ssnd(jps_ocyw)%laction ) CALL cpl_snd( jps_ocyw, isec, RESHAPE ( zoty1, (/jpi,jpj,1/) ), info ) ! ocean y current 1st grid diff --git a/src/OCE/SBC/sbcflx.F90 b/src/OCE/SBC/sbcflx.F90 index c9ad4fcb5..7607f8e24 100644 --- a/src/OCE/SBC/sbcflx.F90 +++ b/src/OCE/SBC/sbcflx.F90 @@ -119,8 +119,8 @@ CONTAINS END DO ! ! fill sf with slf_i and control print CALL fld_fill( sf, slf_i, cn_dir, 'sbc_flx', 'flux formulation for ocean surface boundary condition', 'namsbc_flx' ) - sf(jp_utau)%cltype = 'U' ; sf(jp_utau)%zsgn = -1._wp ! vector field at U point: overwrite default definition of cltype and zsgn - sf(jp_vtau)%cltype = 'V' ; sf(jp_vtau)%zsgn = -1._wp ! vector field at V point: overwrite default definition of cltype and zsgn + sf(jp_utau)%cltype = 'T' ; sf(jp_utau)%zsgn = -1._wp ! vector field at T point: overwrite default definition of cltype and zsgn + sf(jp_vtau)%cltype = 'T' ; sf(jp_vtau)%zsgn = -1._wp ! vector field at T point: overwrite default definition of cltype and zsgn ! ENDIF @@ -147,8 +147,8 @@ CONTAINS ENDIF #endif DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) ! set the ocean fluxes from read fields - utau(ji,jj) = sf(jp_utau)%fnow(ji,jj,1) * umask(ji,jj,1) - vtau(ji,jj) = sf(jp_vtau)%fnow(ji,jj,1) * vmask(ji,jj,1) + utau(ji,jj) = sf(jp_utau)%fnow(ji,jj,1) * tmask(ji,jj,1) + vtau(ji,jj) = sf(jp_vtau)%fnow(ji,jj,1) * tmask(ji,jj,1) qns (ji,jj) = ( sf(jp_qtot)%fnow(ji,jj,1) - sf(jp_qsr)%fnow(ji,jj,1) ) * tmask(ji,jj,1) emp (ji,jj) = sf(jp_emp )%fnow(ji,jj,1) * tmask(ji,jj,1) !!sfx (ji,jj) = sf(jp_sfx )%fnow(ji,jj,1) * tmask(ji,jj,1) @@ -170,19 +170,15 @@ CONTAINS ENDIF ! ENDIF - ! ! module of wind stress and wind speed at T-point - ! Note the use of 0.5*(2-umask) in order to unmask the stress along coastlines + ! + ! module of wind stress and wind speed at T-point zcoef = 1. / ( zrhoa * zcdrag ) - DO_2D( 0, 0, 0, 0 ) - ztx = ( utau(ji-1,jj ) + utau(ji,jj) ) * 0.5_wp * ( 2._wp - MIN( umask(ji-1,jj ,1), umask(ji,jj,1) ) ) - zty = ( vtau(ji ,jj-1) + vtau(ji,jj) ) * 0.5_wp * ( 2._wp - MIN( vmask(ji ,jj-1,1), vmask(ji,jj,1) ) ) - zmod = SQRT( ztx * ztx + zty * zty ) * tmask(ji,jj,1) + DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) + zmod = SQRT( utau(ji,jj) * utau(ji,jj) + vtau(ji,jj) * vtau(ji,jj) ) * tmask(ji,jj,1) taum(ji,jj) = zmod wndm(ji,jj) = SQRT( zmod * zcoef ) !!clem: not used? END_2D ! - CALL lbc_lnk( 'sbcflx', taum, 'T', 1._wp, wndm, 'T', 1._wp ) - ! END SUBROUTINE sbc_flx !!====================================================================== diff --git a/src/OCE/SBC/sbcmod.F90 b/src/OCE/SBC/sbcmod.F90 index 8881e2f4a..81f018a49 100644 --- a/src/OCE/SBC/sbcmod.F90 +++ b/src/OCE/SBC/sbcmod.F90 @@ -158,8 +158,8 @@ CONTAINS ! IF( .NOT.ln_usr ) THEN ! the model calendar needs some specificities (except in user defined case) IF( MOD( rday , rn_Dt ) /= 0. ) CALL ctl_stop( 'the time step must devide the number of second of in a day' ) - IF( MOD( rday , 2. ) /= 0. ) CALL ctl_stop( 'the number of second of in a day must be an even number' ) - IF( MOD( rn_Dt , 2. ) /= 0. ) CALL ctl_stop( 'the time step (in second) must be an even number' ) + IF( MOD( rday , 2. ) /= 0. ) CALL ctl_stop( 'the number of second of in a day must be an even number' ) + IF( MOD( rn_Dt, 2. ) /= 0. ) CALL ctl_stop( 'the time step (in second) must be an even number' ) ENDIF ! !** check option consistency ! @@ -395,16 +395,16 @@ CONTAINS ! ! ---------------------------------------- ! IF( kt /= nit000 ) THEN ! Swap of forcing fields ! ! ! ---------------------------------------- ! - utau_b(:,:) = utau(:,:) ! Swap the ocean forcing fields - vtau_b(:,:) = vtau(:,:) ! (except at nit000 where before fields - qns_b (:,:) = qns (:,:) ! are set at the end of the routine) - emp_b (:,:) = emp (:,:) - sfx_b (:,:) = sfx (:,:) + utau_b(:,:) = utauU(:,:) ! Swap the ocean forcing fields + vtau_b(:,:) = vtauV(:,:) ! (except at nit000 where before fields + qns_b (:,:) = qns (:,:) ! are set at the end of the routine) + emp_b (:,:) = emp (:,:) + sfx_b (:,:) = sfx (:,:) IF( ln_rnf ) THEN rnf_b (:,: ) = rnf (:,: ) rnf_tsc_b(:,:,:) = rnf_tsc(:,:,:) ENDIF - ! + ! ENDIF ! ! ---------------------------------------- ! ! ! forcing field computation ! @@ -420,60 +420,62 @@ CONTAINS ! SELECT CASE( nsbc ) ! Compute ocean surface boundary condition ! ! (i.e. utau,vtau, qns, qsr, emp, sfx) - CASE( jp_usr ) ; CALL usrdef_sbc_oce( kt, Kbb ) ! user defined formulation - CASE( jp_flx ) ; CALL sbc_flx ( kt ) ! flux formulation + CASE( jp_usr ) ; CALL usrdef_sbc_oce( kt, Kbb ) ! user defined formulation + CASE( jp_flx ) ; CALL sbc_flx ( kt ) ! flux formulation CASE( jp_blk ) IF( ll_sas ) CALL sbc_cpl_rcv ( kt, nn_fsbc, nn_ice, Kbb, Kmm ) ! OCE-SAS coupling: SAS receiving fields from OCE !!!!!!!!!!! ATTENTION:ln_wave is not only used for oasis coupling !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! IF( ln_wave ) THEN - IF ( lk_oasis ) CALL sbc_cpl_rcv ( kt, nn_fsbc, nn_ice, Kbb, Kmm ) ! OCE-wave coupling + IF ( lk_oasis ) CALL sbc_cpl_rcv ( kt, nn_fsbc, nn_ice, Kbb, Kmm ) ! OCE-wave coupling CALL sbc_wave ( kt, Kmm ) ENDIF - CALL sbc_blk ( kt ) ! bulk formulation for the ocean + CALL sbc_blk ( kt ) ! bulk formulation for the ocean ! CASE( jp_abl ) IF( ll_sas ) CALL sbc_cpl_rcv ( kt, nn_fsbc, nn_ice, Kbb, Kmm ) ! OCE-SAS coupling: SAS receiving fields from OCE - CALL sbc_abl ( kt ) ! ABL formulation for the ocean + CALL sbc_abl ( kt ) ! ABL formulation for the ocean ! CASE( jp_purecpl ) ; CALL sbc_cpl_rcv ( kt, nn_fsbc, nn_ice, Kbb, Kmm ) ! pure coupled formulation CASE( jp_none ) - IF( ll_opa ) CALL sbc_cpl_rcv ( kt, nn_fsbc, nn_ice, Kbb, Kmm ) ! OCE-SAS coupling: OCE receiving fields from SAS + IF( ll_opa ) CALL sbc_cpl_rcv ( kt, nn_fsbc, nn_ice, Kbb, Kmm ) ! OCE-SAS coupling: OCE receiving fields from SAS END SELECT ! - IF( ln_mixcpl ) CALL sbc_cpl_rcv ( kt, nn_fsbc, nn_ice, Kbb, Kmm ) ! forced-coupled mixed formulation after forcing + IF( ln_mixcpl ) CALL sbc_cpl_rcv ( kt, nn_fsbc, nn_ice, Kbb, Kmm ) ! forced-coupled mixed formulation after forcing ! IF( ln_wave .AND. ln_tauoc ) THEN ! Wave stress reduction - DO_2D( 0, 0, 0, 0) - utau(ji,jj) = utau(ji,jj) * ( tauoc_wave(ji,jj) + tauoc_wave(ji+1,jj) ) * 0.5_wp - vtau(ji,jj) = vtau(ji,jj) * ( tauoc_wave(ji,jj) + tauoc_wave(ji,jj+1) ) * 0.5_wp - END_2D ! - CALL lbc_lnk( 'sbcwave', utau, 'U', -1. ) - CALL lbc_lnk( 'sbcwave', vtau, 'V', -1. ) - ! - taum(:,:) = taum(:,:)*tauoc_wave(:,:) + DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) + utau(ji,jj) = utau(ji,jj) * tauoc_wave(ji,jj) + vtau(ji,jj) = vtau(ji,jj) * tauoc_wave(ji,jj) + taum(ji,jj) = taum(ji,jj) * tauoc_wave(ji,jj) + END_2D ! IF( kt == nit000 ) CALL ctl_warn( 'sbc: You are subtracting the wave stress to the ocean.', & & 'If not requested select ln_tauoc=.false.' ) ! ELSEIF( ln_wave .AND. ln_taw ) THEN ! Wave stress reduction - utau(:,:) = utau(:,:) - tawx(:,:) + twox(:,:) - vtau(:,:) = vtau(:,:) - tawy(:,:) + twoy(:,:) - CALL lbc_lnk( 'sbcwave', utau, 'U', -1. ) - CALL lbc_lnk( 'sbcwave', vtau, 'V', -1. ) - ! - DO_2D( 0, 0, 0, 0) - taum(ji,jj) = sqrt((.5*(utau(ji-1,jj)+utau(ji,jj)))**2 + (.5*(vtau(ji,jj-1)+vtau(ji,jj)))**2) + DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) + utau(ji,jj) = utau(ji,jj) - tawx(ji,jj) + twox(ji,jj) + vtau(ji,jj) = vtau(ji,jj) - tawy(ji,jj) + twoy(ji,jj) + taum(ji,jj) = SQRT( utau(ji,jj)*utau(ji,jj) + vtau(ji,jj)*vtau(ji,jj) ) END_2D ! IF( kt == nit000 ) CALL ctl_warn( 'sbc: You are subtracting the wave stress to the ocean.', & & 'If not requested select ln_taw=.false.' ) ! ENDIF - CALL lbc_lnk( 'sbcmod', taum(:,:), 'T', 1. ) + ! ! IF( ln_icebergs ) THEN ! save pure stresses (with no ice-ocean stress) for use by icebergs - utau_icb(:,:) = utau(:,:) ; vtau_icb(:,:) = vtau(:,:) + ! Note the use of 0.5*(2-umask) in order to unmask the stress along coastlines + ! and the use of MAX(tmask(i,j),tmask(i+1,j) is to mask tau over ice shelves + DO_2D( 0, 0, 0, 0 ) + utau_icb(ji,jj) = 0.5_wp * ( utau(ji,jj) + utau(ji+1,jj) ) * & + & ( 2. - umask(ji,jj,1) ) * MAX( tmask(ji,jj,1), tmask(ji+1,jj,1) ) + vtau_icb(ji,jj) = 0.5_wp * ( vtau(ji,jj) + vtau(ji,jj+1) ) * & + & ( 2. - vmask(ji,jj,1) ) * MAX( tmask(ji,jj,1), tmask(ji,jj+1,1) ) + END_2D + CALL lbc_lnk( 'sbcmod', utau_icb, 'U', -1.0_wp, vtau_icb, 'V', -1.0_wp ) ENDIF ! ! !== Misc. Options ==! @@ -507,9 +509,6 @@ CONTAINS ! Should not be run if ln_diurnal_only IF( l_sbc_clo ) CALL sbc_clo( kt ) -!!$!RBbug do not understand why see ticket 667 -!!$!clem: it looks like it is necessary for the north fold (in certain circumstances). Don't know why. -!!$ CALL lbc_lnk( 'sbcmod', emp, 'T', 1.0_wp ) IF( ll_wd ) THEN ! If near WAD point limit the flux for now zthscl = atanh(rn_wd_sbcfra) ! taper frac default is .999 zwdht(:,:) = ssh(:,:,Kmm) + ht_0(:,:) - rn_wdmin1 ! do this calc of water @@ -534,6 +533,16 @@ CONTAINS emp (:,:) = emp(:,:) * zwght(:,:) END WHERE ENDIF + ! --- calculate utau and vtau on U,V-points --- ! + ! Note the use of 0.5*(2-umask) in order to unmask the stress along coastlines + ! and the use of MAX(tmask(i,j),tmask(i+1,j) is to mask tau over ice shelves + DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) + utauU (ji,jj) = 0.5_wp * ( utau(ji,jj) + utau(ji+1,jj) ) * & + & ( 2. - umask(ji,jj,1) ) * MAX( tmask(ji,jj,1), tmask(ji+1,jj,1) ) + vtauV (ji,jj) = 0.5_wp * ( vtau(ji,jj) + vtau(ji,jj+1) ) * & + & ( 2. - vmask(ji,jj,1) ) * MAX( tmask(ji,jj,1), tmask(ji,jj+1,1) ) + END_2D + IF( nn_hls == 1 ) CALL lbc_lnk( 'sbcmod', utauU, 'U', -1.0_wp, vtauV, 'V', -1.0_wp ) ! IF( kt == nit000 ) THEN ! set the forcing field at nit000 - 1 ! ! ! ---------------------------------------- ! @@ -543,27 +552,28 @@ CONTAINS IF( ln_rstart .AND. .NOT.l_1st_euler ) THEN !* MLF: Restart: read in restart file #endif IF(lwp) WRITE(numout,*) ' nit000-1 surface forcing fields read in the restart file' - CALL iom_get( numror, jpdom_auto, 'utau_b', utau_b ) ! i-stress - CALL iom_get( numror, jpdom_auto, 'vtau_b', vtau_b ) ! j-stress - CALL iom_get( numror, jpdom_auto, 'qns_b', qns_b ) ! non solar heat flux - CALL iom_get( numror, jpdom_auto, 'emp_b', emp_b ) ! freshwater flux + CALL iom_get( numror, jpdom_auto, 'utau_b', utau_b, cd_type = 'U', psgn = -1._wp ) ! i-stress + CALL iom_get( numror, jpdom_auto, 'vtau_b', vtau_b, cd_type = 'V', psgn = -1._wp ) ! j-stress + CALL iom_get( numror, jpdom_auto, 'qns_b', qns_b, cd_type = 'T', psgn = 1._wp ) ! non solar heat flux + CALL iom_get( numror, jpdom_auto, 'emp_b', emp_b, cd_type = 'T', psgn = 1._wp ) ! freshwater flux ! NB: The 3D heat content due to qsr forcing (qsr_hc_b) is treated in traqsr ! To ensure restart capability with 3.3x/3.4 restart files !! to be removed in v3.6 IF( iom_varid( numror, 'sfx_b', ldstop = .FALSE. ) > 0 ) THEN - CALL iom_get( numror, jpdom_auto, 'sfx_b', sfx_b ) ! before salt flux (T-point) + CALL iom_get( numror, jpdom_auto, 'sfx_b', sfx_b, cd_type = 'T', psgn = 1._wp ) ! before salt flux (T-point) ELSE sfx_b (:,:) = sfx(:,:) ENDIF ELSE !* no restart: set from nit000 values IF(lwp) WRITE(numout,*) ' nit000-1 surface forcing fields set to nit000' - utau_b(:,:) = utau(:,:) - vtau_b(:,:) = vtau(:,:) + utau_b(:,:) = utauU(:,:) + vtau_b(:,:) = vtauV(:,:) qns_b (:,:) = qns (:,:) emp_b (:,:) = emp (:,:) sfx_b (:,:) = sfx (:,:) ENDIF ENDIF ! + ! #if defined key_RK3 ! ! ---------------------------------------- ! IF( lrst_oce .AND. lk_SWE ) THEN ! RK3: Write in the ocean restart file ! @@ -578,13 +588,13 @@ CONTAINS IF(lwp) WRITE(numout,*) 'sbc : ocean surface forcing fields written in ocean restart file ', & & 'at it= ', kt,' date= ', ndastp IF(lwp) WRITE(numout,*) '~~~~' - CALL iom_rstput( kt, nitrst, numrow, 'utau_b' , utau ) - CALL iom_rstput( kt, nitrst, numrow, 'vtau_b' , vtau ) - CALL iom_rstput( kt, nitrst, numrow, 'qns_b' , qns ) + CALL iom_rstput( kt, nitrst, numrow, 'utau_b' , utauU ) + CALL iom_rstput( kt, nitrst, numrow, 'vtau_b' , vtauV ) + CALL iom_rstput( kt, nitrst, numrow, 'qns_b' , qns ) ! The 3D heat content due to qsr forcing is treated in traqsr ! CALL iom_rstput( kt, nitrst, numrow, 'qsr_b' , qsr ) - CALL iom_rstput( kt, nitrst, numrow, 'emp_b' , emp ) - CALL iom_rstput( kt, nitrst, numrow, 'sfx_b' , sfx ) + CALL iom_rstput( kt, nitrst, numrow, 'emp_b' , emp ) + CALL iom_rstput( kt, nitrst, numrow, 'sfx_b' , sfx ) ENDIF ! ! ---------------------------------------- ! ! ! Outputs and control print ! @@ -628,8 +638,8 @@ CONTAINS CALL prt_ctl(tab3d_1=tmask , clinfo1=' tmask - : ', mask1=tmask, kdim=jpk ) CALL prt_ctl(tab3d_1=ts(:,:,:,jp_tem,Kmm), clinfo1=' sst - : ', mask1=tmask, kdim=1 ) CALL prt_ctl(tab3d_1=ts(:,:,:,jp_sal,Kmm), clinfo1=' sss - : ', mask1=tmask, kdim=1 ) - CALL prt_ctl(tab2d_1=utau , clinfo1=' utau - : ', mask1=umask, & - & tab2d_2=vtau , clinfo2=' vtau - : ', mask2=vmask ) + CALL prt_ctl(tab2d_1=utau , clinfo1=' utau - : ', mask1=tmask, & + & tab2d_2=vtau , clinfo2=' vtau - : ', mask2=tmask ) ENDIF IF( kt == nitend ) CALL sbc_final ! Close down surface module if necessary diff --git a/src/OCE/TRD/trddyn.F90 b/src/OCE/TRD/trddyn.F90 index 8816056e5..32e95a865 100644 --- a/src/OCE/TRD/trddyn.F90 +++ b/src/OCE/TRD/trddyn.F90 @@ -140,8 +140,10 @@ CONTAINS ! ! ! wind stress trends ALLOCATE( z2dx(jpi,jpj) , z2dy(jpi,jpj) ) - z2dx(:,:) = ( utau_b(:,:) + utau(:,:) ) / ( e3u(:,:,1,Kmm) * rho0 ) - z2dy(:,:) = ( vtau_b(:,:) + vtau(:,:) ) / ( e3v(:,:,1,Kmm) * rho0 ) + 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 ) + END_2D CALL iom_put( "utrd_tau", z2dx ) CALL iom_put( "vtrd_tau", z2dy ) DEALLOCATE( z2dx , z2dy ) diff --git a/src/OCE/TRD/trdglo.F90 b/src/OCE/TRD/trdglo.F90 index 0e4a2aea0..a19770158 100644 --- a/src/OCE/TRD/trdglo.F90 +++ b/src/OCE/TRD/trdglo.F90 @@ -127,10 +127,10 @@ 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 ) - zvt = ( utau_b(ji,jj) + utau(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) + vtau(ji,jj) ) * tmask_i(ji,jj+1) * tmask_i(ji,jj) * vmask(ji,jj,jk) & - & * z1_2rho0 * e1e2v(ji,jj) + 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) & + & * z1_2rho0 * e1e2v(ji,jj) umo(jpdyn_tau) = umo(jpdyn_tau) + zvt vmo(jpdyn_tau) = vmo(jpdyn_tau) + zvs hke(jpdyn_tau) = hke(jpdyn_tau) + uu(ji,jj,1,Kmm) * zvt + vv(ji,jj,1,Kmm) * zvs diff --git a/src/OCE/TRD/trdken.F90 b/src/OCE/TRD/trdken.F90 index d97b3eb80..1920e3e72 100644 --- a/src/OCE/TRD/trdken.F90 +++ b/src/OCE/TRD/trdken.F90 @@ -118,16 +118,18 @@ 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) ) - z2dx(:,:) = uu(:,:,1,Kmm) * ( utau_b(:,:) + utau(:,:) ) * e1e2u(:,:) * umask(:,:,1) - z2dy(:,:) = vv(:,:,1,Kmm) * ( vtau_b(:,:) + vtau(:,:) ) * e1e2v(:,:) * vmask(:,:,1) - zke2d(1:nn_hls,:) = 0._wp ; zke2d(:,1:nn_hls) = 0._wp - DO_2D( 0, nn_hls, 0, nn_hls ) - 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( 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 ) + ! 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.... diff --git a/src/OCE/TRD/trdvor.F90 b/src/OCE/TRD/trdvor.F90 index 974f5b1a9..8515d666d 100644 --- a/src/OCE/TRD/trdvor.F90 +++ b/src/OCE/TRD/trdvor.F90 @@ -105,9 +105,9 @@ CONTAINS CASE( jpdyn_zad ) ; CALL trd_vor_zint( putrd, pvtrd, jpvor_zad, Kmm ) ! Vertical Advection CASE( jpdyn_spg ) ; CALL trd_vor_zint( putrd, pvtrd, jpvor_spg, Kmm ) ! Surface Pressure Grad. CASE( jpdyn_zdf ) ! Vertical Diffusion - DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) ! wind stress trends - ztswu(ji,jj) = 0.5 * ( utau_b(ji,jj) + utau(ji,jj) ) / ( e3u(ji,jj,1,Kmm) * rho0 ) - ztswv(ji,jj) = 0.5 * ( vtau_b(ji,jj) + vtau(ji,jj) ) / ( e3v(ji,jj,1,Kmm) * rho0 ) + DO_2D( 0, 1, 0, 1 ) ! wind stress trends + ztswu(ji,jj) = 0.5 * ( utau_b(ji,jj) + utauU(ji,jj) ) / ( e3u(ji,jj,1,Kmm) * rho0 ) + ztswv(ji,jj) = 0.5 * ( vtau_b(ji,jj) + vtauV(ji,jj) ) / ( e3v(ji,jj,1,Kmm) * rho0 ) END_2D CALL trd_vor_zint( putrd, pvtrd, jpvor_zdf, Kmm ) ! zdf trend including surf./bot. stresses CALL trd_vor_zint( ztswu, ztswv, jpvor_swf, Kmm ) ! surface wind stress @@ -165,7 +165,7 @@ CONTAINS SELECT CASE( ktrd ) ! CASE( jpvor_bfr ) ! bottom friction - DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) + DO_2D( 0, 1, 0, 1 ) ikbu = mbkv(ji,jj) ikbv = mbkv(ji,jj) zudpvor(ji,jj) = putrdvor(ji,jj) * e3u(ji,jj,ikbu,Kmm) * e1u(ji,jj) * umask(ji,jj,ikbu) @@ -173,14 +173,18 @@ CONTAINS END_2D ! CASE( jpvor_swf ) ! wind stress - zudpvor(:,:) = putrdvor(:,:) * e3u(:,:,1,Kmm) * e1u(:,:) * umask(:,:,1) - zvdpvor(:,:) = pvtrdvor(:,:) * e3v(:,:,1,Kmm) * e2v(:,:) * vmask(:,:,1) + DO_2D( 0, 1, 0, 1 ) + zudpvor(ji,jj) = putrdvor(ji,jj) * e3u(ji,jj,1,Kmm) * e1u(ji,jj) * umask(ji,jj,1) + zvdpvor(ji,jj) = pvtrdvor(ji,jj) * e3v(ji,jj,1,Kmm) * e2v(ji,jj) * vmask(ji,jj,1) + END_2D ! END SELECT ! Average except for Beta.V - zudpvor(:,:) = zudpvor(:,:) * r1_hu(:,:,Kmm) - zvdpvor(:,:) = zvdpvor(:,:) * r1_hv(:,:,Kmm) + DO_2D( 0, 1, 0, 1 ) + zudpvor(ji,jj) = zudpvor(ji,jj) * r1_hu(ji,jj,Kmm) + zvdpvor(ji,jj) = zvdpvor(ji,jj) * r1_hv(ji,jj,Kmm) + END_2D ! Curl DO_2D( 0, 0, 0, 0 ) @@ -238,10 +242,11 @@ CONTAINS ! I vertical integration of 3D trends ! ===================================== ! putrdvor and pvtrdvor terms - DO jk = 1,jpk - zudpvor(:,:) = zudpvor(:,:) + putrdvor(:,:,jk) * e3u(:,:,jk,Kmm) * e1u(:,:) * umask(:,:,jk) - zvdpvor(:,:) = zvdpvor(:,:) + pvtrdvor(:,:,jk) * e3v(:,:,jk,Kmm) * e2v(:,:) * vmask(:,:,jk) - END DO + zudpvor(:,:) = 0._wp ; zvdpvor(:,:) = 0._wp + DO_3D( 0, 1, 0, 1, 1, jpk ) + zudpvor(ji,jj) = zudpvor(ji,jj) + putrdvor(ji,jj,jk) * e3u(ji,jj,jk,Kmm) * e1u(ji,jj) * umask(ji,jj,jk) + zvdpvor(ji,jj) = zvdpvor(ji,jj) + pvtrdvor(ji,jj,jk) * e3v(ji,jj,jk,Kmm) * e2v(ji,jj) * vmask(ji,jj,jk) + END_3D ! Planetary vorticity: 2nd computation (Beta.V term) store the vertical sum ! as Beta.V term need intergration, not average @@ -254,8 +259,10 @@ CONTAINS ENDIF ! ! Average - zudpvor(:,:) = zudpvor(:,:) * r1_hu(:,:,Kmm) - zvdpvor(:,:) = zvdpvor(:,:) * r1_hv(:,:,Kmm) + DO_2D( 0, 1, 0, 1 ) + zudpvor(ji,jj) = zudpvor(ji,jj) * r1_hu(ji,jj,Kmm) + zvdpvor(ji,jj) = zvdpvor(ji,jj) * r1_hv(ji,jj,Kmm) + END_2D ! ! Curl DO_2D( 0, 0, 0, 0 ) diff --git a/src/OCE/USR/usrdef_sbc.F90 b/src/OCE/USR/usrdef_sbc.F90 index ce5b785a9..2af3f30d2 100644 --- a/src/OCE/USR/usrdef_sbc.F90 +++ b/src/OCE/USR/usrdef_sbc.F90 @@ -166,19 +166,17 @@ CONTAINS ! seasonal oscillation intensity ztau_sais = 0.015 ztaun = ztau - ztau_sais * COS( (ztime - ztimemax) / (ztimemin - ztimemax) * rpi ) - DO_2D( 1, 1, 1, 1 ) + DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) ! domain from 15deg to 50deg and 1/2 period along 14deg ! so 5/4 of half period with seasonal cycle - utau(ji,jj) = - ztaun * SIN( rpi * (gphiu(ji,jj) - 15.) / (29.-15.) ) - vtau(ji,jj) = ztaun * SIN( rpi * (gphiv(ji,jj) - 15.) / (29.-15.) ) + utau(ji,jj) = - ztaun * SIN( rpi * (gphit(ji,jj) - 15.) / (29.-15.) ) + vtau(ji,jj) = ztaun * SIN( rpi * (gphit(ji,jj) - 15.) / (29.-15.) ) END_2D ! module of wind stress and wind speed at T-point zcoef = 1. / ( zrhoa * zcdrag ) - DO_2D( 0, 0, 0, 0 ) - ztx = utau(ji-1,jj ) + utau(ji,jj) - zty = vtau(ji ,jj-1) + vtau(ji,jj) - zmod = 0.5 * SQRT( ztx * ztx + zty * zty ) + DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) + zmod = SQRT( utau(ji,jj) * utau(ji,jj) + vtau(ji,jj) * vtau(ji,jj) ) taum(ji,jj) = zmod wndm(ji,jj) = SQRT( zmod * zcoef ) END_2D diff --git a/src/OCE/ZDF/zdfosm.F90 b/src/OCE/ZDF/zdfosm.F90 index df648a8c0..490a5d56d 100644 --- a/src/OCE/ZDF/zdfosm.F90 +++ b/src/OCE/ZDF/zdfosm.F90 @@ -486,8 +486,8 @@ CONTAINS & grav * zbeta * swsav(ji,jj) ! OBSBL END_2D DO_2D_OVR( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) - suw0(ji,jj) = -0.5_wp * (utau(ji-1,jj) + utau(ji,jj)) * r1_rho0 * tmask(ji,jj,1) ! Surface upward velocity fluxes - zvw0 = -0.5_wp * (vtau(ji,jj-1) + vtau(ji,jj)) * r1_rho0 * tmask(ji,jj,1) + suw0(ji,jj) = - utau(ji,jj) * r1_rho0 * tmask(ji,jj,1) ! Surface upward velocity fluxes + zvw0 = - vtau(ji,jj) * r1_rho0 * tmask(ji,jj,1) sustar(ji,jj) = MAX( SQRT( SQRT( suw0(ji,jj) * suw0(ji,jj) + zvw0 * zvw0 ) ), & ! Friction velocity (sustar), at & 1e-8_wp ) ! T-point : LMD94 eq. 2 scos_wind(ji,jj) = -1.0_wp * suw0(ji,jj) / ( sustar(ji,jj) * sustar(ji,jj) ) diff --git a/src/OCE/ZDF/zdftke.F90 b/src/OCE/ZDF/zdftke.F90 index 68dd9f6ed..4d3adbd3e 100644 --- a/src/OCE/ZDF/zdftke.F90 +++ b/src/OCE/ZDF/zdftke.F90 @@ -210,7 +210,7 @@ CONTAINS REAL(wp) :: zcdrag = 1.5e-3 ! drag coefficient REAL(wp) :: zbbrau, zbbirau, zri ! local scalars REAL(wp) :: zfact1, zfact2, zfact3 ! - - - REAL(wp) :: ztx2 , zty2 , zcof ! - - + REAL(wp) :: zcof ! - - REAL(wp) :: ztau , zdif ! - - REAL(wp) :: zus , zwlc , zind ! - - REAL(wp) :: zzd_up, zzd_lw ! - - @@ -302,8 +302,8 @@ CONTAINS ! Projection of Stokes drift in the wind stress direction ! DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) - ztaui = 0.5_wp * ( utau(ji,jj) + utau(ji-1,jj) ) - ztauj = 0.5_wp * ( vtau(ji,jj) + vtau(ji,jj-1) ) + ztaui = utau(ji,jj) + ztauj = vtau(ji,jj) z1_norm = 1._wp / MAX( SQRT(ztaui*ztaui+ztauj*ztauj), 1.e-12 ) * tmask(ji,jj,1) zWlc2(ji,jj) = 0.5_wp * z1_norm * ( MAX( ut0sd(ji,jj)*ztaui + vt0sd(ji,jj)*ztauj, 0._wp ) )**2 END_2D @@ -483,9 +483,7 @@ CONTAINS END_2D ELSEIF( nn_etau == 3 ) THEN !* penetration belox the mixed layer (HF variability) DO_3D_OVR( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 ) - ztx2 = utau(ji-1,jj ) + utau(ji,jj) - zty2 = vtau(ji ,jj-1) + vtau(ji,jj) - ztau = 0.5_wp * SQRT( ztx2 * ztx2 + zty2 * zty2 ) * tmask(ji,jj,1) ! module of the mean stress + ztau = SQRT( utau(ji,jj)*utau(ji,jj) + vtau(ji,jj)*vtau(ji,jj) ) * tmask(ji,jj,1) ! module of the mean stress zdif = taum(ji,jj) - ztau ! mean of modulus - modulus of the mean zdif = rhftau_scl * MAX( 0._wp, zdif + rhftau_add ) ! apply some modifications... en(ji,jj,jk) = en(ji,jj,jk) + zbbrau * zdif * EXP( -gdepw(ji,jj,jk,Kmm) / htau(ji,jj) ) & diff --git a/src/OCE/stp2d.F90 b/src/OCE/stp2d.F90 index 96170dc20..06c08b84b 100644 --- a/src/OCE/stp2d.F90 +++ b/src/OCE/stp2d.F90 @@ -151,14 +151,14 @@ CONTAINS ! !* wind forcing *! IF( ln_bt_fw ) THEN DO_2D( 0, 0, 0, 0 ) - Ue_rhs(ji,jj) = Ue_rhs(ji,jj) + r1_rho0 * utau(ji,jj) * r1_hu(ji,jj,Kbb) - Ve_rhs(ji,jj) = Ve_rhs(ji,jj) + r1_rho0 * vtau(ji,jj) * r1_hv(ji,jj,Kbb) + Ue_rhs(ji,jj) = Ue_rhs(ji,jj) + r1_rho0 * utauU(ji,jj) * r1_hu(ji,jj,Kbb) + Ve_rhs(ji,jj) = Ve_rhs(ji,jj) + r1_rho0 * vtauV(ji,jj) * r1_hv(ji,jj,Kbb) END_2D ELSE zztmp = r1_rho0 * r1_2 DO_2D( 0, 0, 0, 0 ) - Ue_rhs(ji,jj) = Ue_rhs(ji,jj) + zztmp * ( utau_b(ji,jj) + utau(ji,jj) ) * r1_hu(ji,jj,Kbb) - Ve_rhs(ji,jj) = Ve_rhs(ji,jj) + zztmp * ( vtau_b(ji,jj) + vtau(ji,jj) ) * r1_hv(ji,jj,Kbb) + Ue_rhs(ji,jj) = Ue_rhs(ji,jj) + zztmp * ( utau_b(ji,jj) + utauU(ji,jj) ) * r1_hu(ji,jj,Kbb) + Ve_rhs(ji,jj) = Ve_rhs(ji,jj) + zztmp * ( vtau_b(ji,jj) + vtauV(ji,jj) ) * r1_hv(ji,jj,Kbb) END_2D ENDIF ! diff --git a/src/SAS/diawri.F90 b/src/SAS/diawri.F90 index 68c13c0b4..8d0c2ece3 100644 --- a/src/SAS/diawri.F90 +++ b/src/SAS/diawri.F90 @@ -317,22 +317,22 @@ CONTAINS ! ENDIF ! + CALL histdef( nid_T, "sozotaux", "Wind Stress along i-axis" , "N/m2" , & ! utau + & jpi, jpj, nh_T, 1 , 1, 1 , - 99, 32, clop, zsto, zout ) + CALL histdef( nid_T, "sometauy", "Wind Stress along j-axis" , "N/m2" , & ! vtau + & jpi, jpj, nh_T, 1 , 1, 1 , - 99, 32, clop, zsto, zout ) CALL histend( nid_T, snc4chunks=snc4set ) ! !!! nid_U : 3D CALL histdef( nid_U, "ssu_m", "Velocity component in x-direction", "m/s" , & ! ssu & jpi, jpj, nh_U, 1 , 1, 1 , - 99, 32, clop, zsto, zout ) - CALL histdef( nid_U, "sozotaux", "Wind Stress along i-axis" , "N/m2" , & ! utau - & jpi, jpj, nh_U, 1 , 1, 1 , - 99, 32, clop, zsto, zout ) CALL histend( nid_U, snc4chunks=snc4set ) ! !!! nid_V : 3D CALL histdef( nid_V, "ssv_m", "Velocity component in y-direction", "m/s", & ! ssv_m & jpi, jpj, nh_V, 1 , 1, 1 , - 99, 32, clop, zsto, zout ) - CALL histdef( nid_V, "sometauy", "Wind Stress along j-axis" , "N/m2" , & ! vtau - & jpi, jpj, nh_V, 1 , 1, 1 , - 99, 32, clop, zsto, zout ) CALL histend( nid_V, snc4chunks=snc4set ) @@ -366,6 +366,8 @@ CONTAINS CALL histwrite( nid_T, "soshfldo", it, qsr , ndim_hT, ndex_hT ) ! solar heat flux CALL histwrite( nid_T, "soicecov", it, fr_i , ndim_hT, ndex_hT ) ! ice fraction CALL histwrite( nid_T, "sowindsp", it, wndm , ndim_hT, ndex_hT ) ! wind speed + CALL histwrite( nid_T, "sozotaux", it, utau , ndim_hT, ndex_hT ) ! i-wind stress + CALL histwrite( nid_T, "sometauy", it, vtau , ndim_hT, ndex_hT ) ! j-wind stress ! IF( ln_abl ) THEN ALLOCATE( zw3d_abl(jpi,jpj,jpka) ) @@ -393,11 +395,9 @@ CONTAINS ! Write fields on U grid CALL histwrite( nid_U, "ssu_m" , it, ssu_m , ndim_hU, ndex_hU ) ! i-current speed - CALL histwrite( nid_U, "sozotaux", it, utau , ndim_hU, ndex_hU ) ! i-wind stress ! Write fields on V grid CALL histwrite( nid_V, "ssv_m" , it, ssv_m , ndim_hV, ndex_hV ) ! j-current speed - CALL histwrite( nid_V, "sometauy", it, vtau , ndim_hV, ndex_hV ) ! j-wind stress ! 3. Close all files ! --------------------------------------- diff --git a/src/SWE/stpmlf.F90 b/src/SWE/stpmlf.F90 index 5494248c2..4c49e55aa 100644 --- a/src/SWE/stpmlf.F90 +++ b/src/SWE/stpmlf.F90 @@ -114,9 +114,9 @@ CONTAINS zrhs_u = - grav * ( ssh(ji+1,jj,Nnn) - ssh(ji,jj,Nnn) ) * r1_e1u(ji,jj) zrhs_v = - grav * ( ssh(ji,jj+1,Nnn) - ssh(ji,jj,Nnn) ) * r1_e2v(ji,jj) ! ! wind stress and layer friction - zrhs_u = zrhs_u + z1_2rho0 * ( utau_b(ji,jj) + utau(ji,jj) ) / e3u(ji,jj,jk,Nnn) & + zrhs_u = zrhs_u + z1_2rho0 * ( utau_b(ji,jj) + utauU(ji,jj) ) / e3u(ji,jj,jk,Nnn) & & - rn_rfr * uu(ji,jj,jk,Nbb) - zrhs_v = zrhs_v + z1_2rho0 * ( vtau_b(ji,jj) + vtau(ji,jj) ) / e3v(ji,jj,jk,Nnn) & + zrhs_v = zrhs_v + z1_2rho0 * ( vtau_b(ji,jj) + vtauV(ji,jj) ) / e3v(ji,jj,jk,Nnn) & & - rn_rfr * vv(ji,jj,jk,Nbb) ! ! ==> RHS uu(ji,jj,jk,Nrhs) = uu(ji,jj,jk,Nrhs) + zrhs_u diff --git a/src/SWE/stprk3.F90 b/src/SWE/stprk3.F90 index 0e1a1910f..2fded6aea 100644 --- a/src/SWE/stprk3.F90 +++ b/src/SWE/stprk3.F90 @@ -135,9 +135,9 @@ CONTAINS zrhs_v = - grav * ( ssh(ji,jj+1,Nbb) - ssh(ji,jj,Nbb) ) * r1_e2v(ji,jj) #if defined key_RK3all ! ! wind stress and layer friction - zrhs_u = zrhs_u + r1_rho0 * ( z5_6*utau_b(ji,jj) + (1._wp - z5_6)*utau(ji,jj) ) / e3u(ji,jj,jk,Nbb) & + zrhs_u = zrhs_u + r1_rho0 * ( z5_6*utau_b(ji,jj) + (1._wp - z5_6)*utauU(ji,jj) ) / e3u(ji,jj,jk,Nbb) & & - rn_rfr * uu(ji,jj,jk,Nbb) - zrhs_v = zrhs_v + r1_rho0 * ( z5_6*vtau_b(ji,jj) + (1._wp - z5_6)*vtau(ji,jj) ) / e3v(ji,jj,jk,Nbb) & + zrhs_v = zrhs_v + r1_rho0 * ( z5_6*vtau_b(ji,jj) + (1._wp - z5_6)*vtauV(ji,jj) ) / e3v(ji,jj,jk,Nbb) & & - rn_rfr * vv(ji,jj,jk,Nbb) #endif ! ! ==> RHS @@ -201,9 +201,9 @@ CONTAINS zrhs_v = - grav * ( ssh(ji,jj+1,Nnn) - ssh(ji,jj,Nnn) ) * r1_e2v(ji,jj) #if defined key_RK3all ! ! wind stress and layer friction - zrhs_u = zrhs_u + r1_rho0 * ( z3_4*utau_b(ji,jj) + (1._wp - z3_4)*utau(ji,jj) ) / e3u(ji,jj,jk,Nnn) & + zrhs_u = zrhs_u + r1_rho0 * ( z3_4*utau_b(ji,jj) + (1._wp - z3_4)*utauU(ji,jj) ) / e3u(ji,jj,jk,Nnn) & & - rn_rfr * uu(ji,jj,jk,Nbb) - zrhs_v = zrhs_v + r1_rho0 * ( z3_4*vtau_b(ji,jj) + (1._wp - z3_4)*vtau(ji,jj) ) / e3v(ji,jj,jk,Nnn) & + zrhs_v = zrhs_v + r1_rho0 * ( z3_4*vtau_b(ji,jj) + (1._wp - z3_4)*vtauV(ji,jj) ) / e3v(ji,jj,jk,Nnn) & & - rn_rfr * vv(ji,jj,jk,Nbb) #endif ! ! ==> RHS @@ -265,9 +265,9 @@ CONTAINS zrhs_u = - grav * ( ssh(ji+1,jj,Nnn) - ssh(ji,jj,Nnn) ) * r1_e1u(ji,jj) zrhs_v = - grav * ( ssh(ji,jj+1,Nnn) - ssh(ji,jj,Nnn) ) * r1_e2v(ji,jj) ! ! wind stress and layer friction - zrhs_u = zrhs_u + z1_2rho0 * ( utau_b(ji,jj) + utau(ji,jj) ) / e3u(ji,jj,jk,Nnn) & + zrhs_u = zrhs_u + z1_2rho0 * ( utau_b(ji,jj) + utauU(ji,jj) ) / e3u(ji,jj,jk,Nnn) & & - rn_rfr * uu(ji,jj,jk,Nbb) - zrhs_v = zrhs_v + z1_2rho0 * ( vtau_b(ji,jj) + vtau(ji,jj) ) / e3v(ji,jj,jk,Nnn) & + zrhs_v = zrhs_v + z1_2rho0 * ( vtau_b(ji,jj) + vtauV(ji,jj) ) / e3v(ji,jj,jk,Nnn) & & - rn_rfr * vv(ji,jj,jk,Nbb) ! ! ==> RHS uu(ji,jj,jk,Nrhs) = uu(ji,jj,jk,Nrhs) + zrhs_u diff --git a/tests/BENCH/MY_SRC/usrdef_sbc.F90 b/tests/BENCH/MY_SRC/usrdef_sbc.F90 index e9ac757b6..4bd129228 100644 --- a/tests/BENCH/MY_SRC/usrdef_sbc.F90 +++ b/tests/BENCH/MY_SRC/usrdef_sbc.F90 @@ -109,7 +109,7 @@ CONTAINS vtau_ice(ji,jj) = 0.1_wp + zztmp END_2D - CALL lbc_lnk( 'usrdef_sbc', utau_ice, 'U', -1., vtau_ice, 'V', -1. ) + CALL lbc_lnk( 'usrdef_sbc', utau_ice, 'T', -1., vtau_ice, 'T', -1. ) #endif ! END SUBROUTINE usrdef_sbc_ice_tau diff --git a/tests/CANAL/EXPREF/file_def_nemo-oce.xml b/tests/CANAL/EXPREF/file_def_nemo-oce.xml index 689b95dd6..45712ad6d 100644 --- a/tests/CANAL/EXPREF/file_def_nemo-oce.xml +++ b/tests/CANAL/EXPREF/file_def_nemo-oce.xml @@ -20,15 +20,15 @@ <field field_ref="ssrelpotvor" /> <field field_ref="saltc" /> <field field_ref="salt2c" /> + <field field_ref="utau" /> + <field field_ref="vtau" /> </file> <file id="file3" name_suffix="_grid_U" description="ocean U grid variables" > - <field field_ref="utau" /> <field field_ref="uoce" /> </file> <file id="file4" name_suffix="_grid_V" description="ocean V grid variables" > - <field field_ref="vtau" /> <field field_ref="voce" /> </file> diff --git a/tests/DIA_GPU/EXPREF/file_def_nemo-oce.xml b/tests/DIA_GPU/EXPREF/file_def_nemo-oce.xml index c043f56fe..04590b55f 100644 --- a/tests/DIA_GPU/EXPREF/file_def_nemo-oce.xml +++ b/tests/DIA_GPU/EXPREF/file_def_nemo-oce.xml @@ -34,6 +34,8 @@ <field field_ref="qt_oce" name="qt_oce" /> <field field_ref="saltflx" name="sfx" /> <field field_ref="taum" name="taum" /> + <field field_ref="utau" name="tauuo" /> + <field field_ref="vtau" name="tauvo" /> <field field_ref="wspd" name="windsp" /> <field field_ref="precip" name="precip" /> <!-- ice and snow --> @@ -44,7 +46,6 @@ <field field_ref="e3u" /> <field field_ref="ssu" name="uos" /> <field field_ref="uoce" name="uo" operation="instant" freq_op="5d" > @uoce_e3u / @e3u </field> - <field field_ref="utau" name="tauuo" /> <field field_ref="uocetr_eff" name="uocetr_eff" /> <!-- available with diaar5 --> <field field_ref="u_masstr" name="vozomatr" /> @@ -56,7 +57,6 @@ <field field_ref="e3v" /> <field field_ref="ssv" name="vos" /> <field field_ref="voce" name="vo" operation="instant" freq_op="5d" > @voce_e3v / @e3v </field> - <field field_ref="vtau" name="tauvo" /> <field field_ref="vocetr_eff" name="vocetr_eff" /> <!-- available with diaar5 --> <field field_ref="v_masstr" name="vomematr" /> diff --git a/tests/DONUT/EXPREF/file_def_nemo-oce.xml b/tests/DONUT/EXPREF/file_def_nemo-oce.xml index 40a640ffb..0f5dc0b4e 100644 --- a/tests/DONUT/EXPREF/file_def_nemo-oce.xml +++ b/tests/DONUT/EXPREF/file_def_nemo-oce.xml @@ -28,6 +28,8 @@ <field field_ref="qt_oce" name="qt_oce" /> <field field_ref="saltflx" name="sfx" /> <field field_ref="taum" name="taum" /> + <field field_ref="utau" name="tauuo" /> + <field field_ref="vtau" name="tauvo" /> <field field_ref="wspd" name="windsp" /> <field field_ref="precip" name="precip" /> </file> @@ -36,14 +38,12 @@ <field field_ref="e3u" /> <field field_ref="ssu" name="uos" /> <field field_ref="uoce" name="uo" operation="instant" freq_op="1d" > @uoce_e3u / @e3u </field> - <field field_ref="utau" name="tauuo" /> </file> <file id="file13" name_suffix="_grid_V" description="ocean V grid variables" > <field field_ref="e3v" /> <field field_ref="ssv" name="vos" /> <field field_ref="voce" name="vo" operation="instant" freq_op="1d" > @voce_e3v / @e3v </field> - <field field_ref="vtau" name="tauvo" /> </file> </file_group> diff --git a/tests/ICE_RHEO/EXPREF/file_def_nemo-oce.xml b/tests/ICE_RHEO/EXPREF/file_def_nemo-oce.xml index 0dbbd2144..b307458f5 100644 --- a/tests/ICE_RHEO/EXPREF/file_def_nemo-oce.xml +++ b/tests/ICE_RHEO/EXPREF/file_def_nemo-oce.xml @@ -13,14 +13,17 @@ <file_group id="1ts" output_freq="1ts" output_level="10" enabled=".TRUE."/> <!-- 1 time step files --> <file_group id="1h" output_freq="1h" output_level="10" enabled=".TRUE."> <!-- 1h files --> + <file id="file1" name_suffix="_grid_T" description="ocean T grid variables" > + <field field_ref="utau" name="sozotaux" /> + <field field_ref="vtau" name="sometauy" /> + </file> + <file id="file2" name_suffix="_grid_U" description="ocean U grid variables" > <field field_ref="uoce" name="vozocrtx" /> - <field field_ref="utau" name="sozotaux" /> </file> <file id="file3" name_suffix="_grid_V" description="ocean V grid variables" > <field field_ref="voce" name="vomecrty" /> - <field field_ref="vtau" name="sometauy" /> </file> </file_group> diff --git a/tests/ICE_RHEO/MY_SRC/usrdef_sbc.F90 b/tests/ICE_RHEO/MY_SRC/usrdef_sbc.F90 index 8d5fe4eb6..2a22cfc80 100644 --- a/tests/ICE_RHEO/MY_SRC/usrdef_sbc.F90 +++ b/tests/ICE_RHEO/MY_SRC/usrdef_sbc.F90 @@ -125,37 +125,21 @@ CONTAINS windu(ji,jj) = Umax/sqrt(d*1000)*(d-2*mig(ji)*res)/((d-2*mig(ji)*res)**2+(d-2*mjg(jj)*res)**2*Rwind**2)**(1/4)*min(kt*30./21600,1.) windv(ji,jj) = Umax/sqrt(d*1000)*(d-2*mjg(jj)*res)/((d-2*mig(ji)*res)**2+(d-2*mjg(jj)*res)**2*Rwind**2)**(1/4)*Rwind*min(kt*30./21600,1.) END_2D - CALL lbc_lnk( 'usrdef_sbc', windu, 'U', -1., windv, 'V', -1. ) - - wndm_ice(:,:) = 0._wp !!gm brutal.... ! ------------------------------------------------------------ ! - ! Wind module relative to the moving ice ( U10m - U_ice ) ! + ! Wind module and stress relative to the moving ice ( U10m - U_ice ) ! ! ------------------------------------------------------------ ! - ! C-grid ice dynamics : U & V-points (same as ocean) DO_2D( 0, 0, 0, 0 ) - zwndi_t = ( windu(ji,jj) - r_vfac * 0.5 * ( u_ice(ji-1,jj ) + u_ice(ji,jj) ) ) + zwndi_t = ( windu(ji,jj) - r_vfac * 0.5 * ( u_ice(ji-1,jj) + u_ice(ji,jj) ) ) zwndj_t = ( windv(ji,jj) - r_vfac * 0.5 * ( v_ice(ji,jj-1) + v_ice(ji,jj) ) ) + ! wndm_ice(ji,jj) = SQRT( zwndi_t * zwndi_t + zwndj_t * zwndj_t ) * tmask(ji,jj,1) + ! + utau_ice(ji,jj) = zrhoa * Cd_atm * wndm_ice(ji,jj) * zwndi_t + vtau_ice(ji,jj) = zrhoa * Cd_atm * wndm_ice(ji,jj) * zwndj_t END_2D - CALL lbc_lnk( 'usrdef_sbc', wndm_ice, 'T', 1. ) - - !!gm brutal.... - utau_ice (:,:) = 0._wp - vtau_ice (:,:) = 0._wp - !!gm end + CALL lbc_lnk( 'usrdef_sbc', wndm_ice, 'T', 1., utau_ice, 'T', -1., vtau_ice, 'T', -1. ) - ! ------------------------------------------------------------ ! - ! Wind stress relative to the moving ice ( U10m - U_ice ) ! - ! ------------------------------------------------------------ ! - ! C-grid ice dynamics : U & V-points (same as ocean) - DO_2D( 0, 0, 0, 0 ) - utau_ice(ji,jj) = 0.5 * zrhoa * Cd_atm * ( wndm_ice(ji+1,jj ) + wndm_ice(ji,jj) ) & - & * ( 0.5 * (windu(ji+1,jj) + windu(ji,jj) ) - r_vfac * u_ice(ji,jj) ) - vtau_ice(ji,jj) = 0.5 * zrhoa * Cd_atm * ( wndm_ice(ji,jj+1 ) + wndm_ice(ji,jj) ) & - & * ( 0.5 * (windv(ji,jj+1) + windv(ji,jj) ) - r_vfac * v_ice(ji,jj) ) - END_2D - CALL lbc_lnk( 'usrdef_sbc', utau_ice, 'U', -1., vtau_ice, 'V', -1. ) ! END SUBROUTINE usrdef_sbc_ice_tau diff --git a/tests/SWG/EXPREF/file_def_nemo-oce.xml b/tests/SWG/EXPREF/file_def_nemo-oce.xml index fb16bbbd5..3b1c4d50b 100644 --- a/tests/SWG/EXPREF/file_def_nemo-oce.xml +++ b/tests/SWG/EXPREF/file_def_nemo-oce.xml @@ -27,6 +27,8 @@ <field field_ref="sKE" name="sKE" operation = "instant" /> <field field_ref="ssh" name="sossheig" operation = "instant" /> <field field_ref="wetdep" name="hswe_wd" operation = "instant" /> + <field field_ref="utau" name="sozotaux" operation = "instant" /> + <field field_ref="vtau" name="sometauy" operation = "instant" /> </file> <file id="file2" name_suffix="_grid_U" description="ocean U grid variables" > @@ -34,7 +36,6 @@ <field field_ref="e3u_0" name="e3u_0" operation = "instant" /> <field field_ref="hu" name="hu" operation = "instant" /> <field field_ref="ssu" name="ssu" operation = "instant" /> - <field field_ref="utau" name="sozotaux" operation = "instant" /> </file> <file id="file3" name_suffix="_grid_V" description="ocean V grid variables" > @@ -42,7 +43,6 @@ <field field_ref="e3v_0" name="e3v_0" operation = "instant" /> <field field_ref="hv" name="hv" operation = "instant" /> <field field_ref="ssv" name="ssv" operation = "instant" /> - <field field_ref="vtau" name="sometauy" operation = "instant" /> </file> <file id="file5" name_suffix="_grid_F" description="ocean F grid variables" > diff --git a/tests/SWG/MY_SRC/usrdef_sbc.F90 b/tests/SWG/MY_SRC/usrdef_sbc.F90 index a30699547..66ea77ad2 100644 --- a/tests/SWG/MY_SRC/usrdef_sbc.F90 +++ b/tests/SWG/MY_SRC/usrdef_sbc.F90 @@ -61,7 +61,7 @@ CONTAINS REAL(wp) :: ztauu, ztauv ! wind intensity projeted REAL(wp) :: zrhoa = 1.22 ! Air density kg/m3 REAL(wp) :: zcdrag = 1.5e-3 ! drag coefficient - REAL(wp) :: ztx, zty, zmod, zcoef ! temporary variables + REAL(wp) :: zmod, zcoef ! temporary variables !!--------------------------------------------------------------------- ! ---------------------------- ! @@ -88,20 +88,17 @@ CONTAINS DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) ! length of the domain : 2000km x 2000km - utau(ji,jj) = - ztauu * COS( rpi * gphiu(ji,jj) / 2000000._wp) - vtau(ji,jj) = - ztauv * COS( rpi * gphiv(ji,jj) / 2000000._wp) + utau(ji,jj) = - ztauu * COS( rpi * gphit(ji,jj) / 2000000._wp) + vtau(ji,jj) = - ztauv * COS( rpi * gphit(ji,jj) / 2000000._wp) END_2D ! module of wind stress and wind speed at T-point zcoef = 1. / ( zrhoa * zcdrag ) - DO_2D( 0, 0, 0, 0 ) - ztx = utau(ji-1,jj ) + utau(ji,jj) - zty = vtau(ji ,jj-1) + vtau(ji,jj) - zmod = 0.5 * SQRT( ztx * ztx + zty * zty ) + DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) + zmod = SQRT( utau(ji,jj) * utau(ji,jj) + vtau(ji,jj) * vtau(ji,jj) ) taum(ji,jj) = zmod wndm(ji,jj) = SQRT( zmod * zcoef ) END_2D - CALL lbc_lnk( 'usrdef_sbc', taum(:,:), 'T', 1. , wndm(:,:), 'T', 1. ) ! END SUBROUTINE usrdef_sbc_oce diff --git a/tests/WAD/EXPREF/file_def_nemo-oce.xml b/tests/WAD/EXPREF/file_def_nemo-oce.xml index a9c10b7e0..90124a636 100644 --- a/tests/WAD/EXPREF/file_def_nemo-oce.xml +++ b/tests/WAD/EXPREF/file_def_nemo-oce.xml @@ -69,16 +69,16 @@ <field field_ref="qt" name="sohefldo" /> <field field_ref="mldr10_1" name="somxl010" /> <field field_ref="mldkz5" name="somixhgt" /> + <field field_ref="utau" name="sozotaux" /> + <field field_ref="vtau" name="sometauy" /> </file> <file id="file2" name_suffix="_grid_U" description="ocean U grid variables" > <field field_ref="uoce" name="vozocrtx" /> - <field field_ref="utau" name="sozotaux" /> </file> <file id="file3" name_suffix="_grid_V" description="ocean V grid variables" > <field field_ref="voce" name="vomecrty" /> - <field field_ref="vtau" name="sometauy" /> </file> <file id="file4" name_suffix="_grid_W" description="ocean W grid variables" > -- GitLab