From 917994a2a85a44afa79d0101cab6964ad7376365 Mon Sep 17 00:00:00 2001 From: Guillaume Samson <guillaume.samson@mercator-ocean.fr> Date: Thu, 29 Sep 2022 13:33:26 +0000 Subject: [PATCH] Resolve "Endless Summer 2022 - ABL halo cleanup" --- cfgs/ORCA2_ICE_PISCES/EXPREF/namelist_cfg | 11 + cfgs/SHARED/grid_def_nemo.xml | 6 +- cfgs/SHARED/namelist_ref | 12 +- makenemo | 4 +- sette/BATCH_TEMPLATE/batch-X64_AA_INTEL_OMPI | 2 +- sette/all_functions.sh | 2 +- sette/sette.sh | 5 +- src/ABL/abl.F90 | 26 +- src/ABL/ablmod.F90 | 542 +++++++++---------- src/ABL/ablrst.F90 | 6 +- src/ABL/sbcabl.F90 | 12 +- src/OCE/SBC/sbcblk.F90 | 1 + src/OCE/ZDF/zdfgls.F90 | 3 +- src/OCE/ZDF/zdftke.F90 | 3 +- 14 files changed, 315 insertions(+), 320 deletions(-) diff --git a/cfgs/ORCA2_ICE_PISCES/EXPREF/namelist_cfg b/cfgs/ORCA2_ICE_PISCES/EXPREF/namelist_cfg index 4b17c1029..c6e94339d 100644 --- a/cfgs/ORCA2_ICE_PISCES/EXPREF/namelist_cfg +++ b/cfgs/ORCA2_ICE_PISCES/EXPREF/namelist_cfg @@ -118,6 +118,17 @@ &namsbc_abl ! Atmospheric Boundary Layer formulation (ln_abl = T) !----------------------------------------------------------------------- cn_dom = 'dom_cfg_abl_L25Z10' + + nn_dyn_restore = 2 ! restoring option for dynamical ABL variables: = 0 no restoring + ! = 1 equatorial restoring + ! = 2 global restoring + ! = 1 equatorial restoring + ! = 2 global restoring + rn_ldyn_min = 24. ! dynamics nudging magnitude inside the ABL [hour] (~3 rn_Dt) + rn_ldyn_max = 6. ! dynamics nudging magnitude above the ABL [hour] (~1 rn_Dt) + rn_ltra_min = 24. ! tracers nudging magnitude inside the ABL [hour] (~3 rn_Dt) + rn_ltra_max = 6. ! tracers nudging magnitude above the ABL [hour] (~1 rn_Dt) + rn_vfac = 1. / !----------------------------------------------------------------------- &namtra_qsr ! penetrative solar radiation (ln_traqsr =T) diff --git a/cfgs/SHARED/grid_def_nemo.xml b/cfgs/SHARED/grid_def_nemo.xml index 374f67c65..0120bb3cc 100644 --- a/cfgs/SHARED/grid_def_nemo.xml +++ b/cfgs/SHARED/grid_def_nemo.xml @@ -318,14 +318,14 @@ <!-- ABL grid definition --> <grid id="grid_TA_2D"> - <domain domain_ref="grid_T" /> + <domain domain_ref="grid_T_inner" /> </grid> <grid id="grid_TA_3D"> - <domain domain_ref="grid_T" /> + <domain domain_ref="grid_T_inner" /> <axis id="ght_abl" /> </grid> <grid id="grid_WA_3D"> - <domain domain_ref="grid_T" /> + <domain domain_ref="grid_T_inner" /> <axis id="ghw_abl" /> </grid> <!-- --> diff --git a/cfgs/SHARED/namelist_ref b/cfgs/SHARED/namelist_ref index 28a67e4d5..00688501d 100644 --- a/cfgs/SHARED/namelist_ref +++ b/cfgs/SHARED/namelist_ref @@ -312,7 +312,7 @@ &namsbc_abl ! Atmospheric Boundary Layer formulation (ln_abl = T) !----------------------------------------------------------------------- cn_dir = './' ! root directory for the location of the ABL grid file - cn_dom = 'dom_cfg_abl.nc' + cn_dom = 'dom_cfg_abl' cn_ablrst_in = "restart_abl" ! suffix of abl restart name (input) cn_ablrst_out = "restart_abl" ! suffix of abl restart name (output) @@ -326,11 +326,11 @@ nn_dyn_restore = 0 ! restoring option for dynamical ABL variables: = 0 no restoring ! = 1 equatorial restoring ! = 2 global restoring - rn_vfac = 0. - rn_ldyn_min = 4.5 ! dynamics nudging magnitude inside the ABL [hour] (~3 rn_Dt) - rn_ldyn_max = 1.5 ! dynamics nudging magnitude above the ABL [hour] (~1 rn_Dt) - rn_ltra_min = 4.5 ! tracers nudging magnitude inside the ABL [hour] (~3 rn_Dt) - rn_ltra_max = 1.5 ! tracers nudging magnitude above the ABL [hour] (~1 rn_Dt) + rn_ldyn_min = 0. ! dynamics nudging magnitude inside the ABL [hour] (~3 rn_Dt) + rn_ldyn_max = 0. ! dynamics nudging magnitude above the ABL [hour] (~1 rn_Dt) + rn_ltra_min = 0. ! tracers nudging magnitude inside the ABL [hour] (~3 rn_Dt) + rn_ltra_max = 0. ! tracers nudging magnitude above the ABL [hour] (~1 rn_Dt) + rn_vfac = 0. nn_amxl = 0 ! mixing length: = 0 Deardorff 80 length-scale ! = 1 length-scale based on the distance to the PBL height ! = 2 Bougeault & Lacarrere 89 length-scale diff --git a/makenemo b/makenemo index 8f429bd2d..a7c8bd906 100755 --- a/makenemo +++ b/makenemo @@ -223,8 +223,8 @@ done if [ -n "$x_name" ] # is the configuration existing in cfgs or tests? or is it a new conf? then - incfg=$( find $CFGS_DIR -type d -name $x_name | wc -l ) # this configuration exists in CFGS_DIR - intst=$( find $TESTS_DIR -type d -name $x_name | wc -l ) # this configuration exists in TESTS_DIR + incfg=$( find $CFGS_DIR -type d -name $x_name | wc -l ) # this configuration exists in CFGS_DIR + intst=$( find $TESTS_DIR -type d -name $x_name | wc -l ) # this configuration exists in TESTS_DIR if [ $(( $incfg + $intst )) -eq 2 ] ; then echo -e "\033[0;31m\nERROR: the $x_name directory is found twice: in $CFGS_DIR and in $TESTS_DIR\033[0m\n" exit 4 diff --git a/sette/BATCH_TEMPLATE/batch-X64_AA_INTEL_OMPI b/sette/BATCH_TEMPLATE/batch-X64_AA_INTEL_OMPI index bdd075c6c..3eca6d7c6 100644 --- a/sette/BATCH_TEMPLATE/batch-X64_AA_INTEL_OMPI +++ b/sette/BATCH_TEMPLATE/batch-X64_AA_INTEL_OMPI @@ -98,7 +98,7 @@ EOF # # Run SPMD case # - time ./nemo + $MPIRUN --ntasks=TOTAL_NPROCS ./nemo fi # diff --git a/sette/all_functions.sh b/sette/all_functions.sh index 8bc6e975c..8737aa4be 100755 --- a/sette/all_functions.sh +++ b/sette/all_functions.sh @@ -172,7 +172,7 @@ clean_config() { # define validation dir set_valid_dir () { if [ ${DETACHED_HEAD} == "no" ] ; then - branchname=$( git branch --show-current ) + branchname=$( git -C ${MAIN_DIR} branch --show-current ) REVISION_NB=$( git -C ${MAIN_DIR} rev-list --abbrev-commit origin/$branchname | head -1l ) else REVISION_NB=${DETACHED_CMIT} diff --git a/sette/sette.sh b/sette/sette.sh index 7537af5e8..cc2846f44 100755 --- a/sette/sette.sh +++ b/sette/sette.sh @@ -40,10 +40,10 @@ export USER_INPUT='yes' # Default: yes => request user input on decisions # # Check that git branch is usable export DETACHED_HEAD="no" -git branch --show-current >& /dev/null +git -C ${MAIN_DIR} branch --show-current >& /dev/null if [[ $? == 0 ]] ; then # subdirectory below NEMO_VALIDATION_DIR defaults to branchname - export SETTE_SUB_VAL="$(git branch --show-current)" + export SETTE_SUB_VAL="$(git -C ${MAIN_DIR} branch --show-current)" if [ -z $SETTE_SUB_VAL ] ; then # Probabably on a detached HEAD (possibly testing an old commit). # Verify this and try to recover original commit @@ -76,6 +76,7 @@ if [ $# -gt 0 ]; then y) dry_run=1 echo "";; b) NEMO_DEBUG="-b" + echo "-b: Nemo will be compiled with DEBUG options" echo "";; r) NO_REPORT=1 echo "";; diff --git a/src/ABL/abl.F90 b/src/ABL/abl.F90 index ccffa088d..c89f0c1b6 100644 --- a/src/ABL/abl.F90 +++ b/src/ABL/abl.F90 @@ -41,6 +41,8 @@ MODULE abl ! INTEGER , PUBLIC :: nt_n, nt_a !: now / after indices (equal 1 or 2) ! + !! * Substitutions +# include "do_loop_substitute.h90" !!---------------------------------------------------------------------- !! NEMO/OPA 4.0 , NEMO Consortium (2011) !! $Id: abl.F90 4990 2014-12-15 16:42:49Z timgraham $ @@ -55,18 +57,18 @@ CONTAINS INTEGER :: ierr !!---------------------------------------------------------------------- ! - ALLOCATE( u_abl (1:jpi,1:jpj,1:jpka,jptime ), & - & v_abl (1:jpi,1:jpj,1:jpka,jptime ), & - & tq_abl (1:jpi,1:jpj,1:jpka,jptime,jptq), & - & tke_abl (1:jpi,1:jpj,1:jpka,jptime ), & - & avm_abl (1:jpi,1:jpj,1:jpka ), & - & avt_abl (1:jpi,1:jpj,1:jpka ), & - & mxld_abl(1:jpi,1:jpj,1:jpka ), & - & mxlm_abl(1:jpi,1:jpj,1:jpka ), & - & fft_abl (1:jpi,1:jpj ), & - & pblh (1:jpi,1:jpj ), & - & msk_abl (1:jpi,1:jpj ), & - & rest_eq (1:jpi,1:jpj ), & + ALLOCATE( u_abl (A2D(0),1:jpka,jptime ), & + & v_abl (A2D(0),1:jpka,jptime ), & + & tq_abl (A2D(0),1:jpka,jptime,jptq), & + & tke_abl (A2D(0),1:jpka,jptime ), & + & avm_abl (A2D(0),1:jpka ), & + & avt_abl (A2D(0),1:jpka ), & + & mxld_abl(A2D(0),1:jpka ), & + & mxlm_abl(A2D(0),1:jpka ), & + & fft_abl (A2D(0) ), & + & pblh (A2D(1) ), & ! needed for optional smoothing + & msk_abl (A2D(0) ), & + & rest_eq (A2D(0) ), & & e3t_abl (1:jpka), e3w_abl(1:jpka) , & & ght_abl (1:jpka), ghw_abl(1:jpka) , STAT=ierr ) ! diff --git a/src/ABL/ablmod.F90 b/src/ABL/ablmod.F90 index ffaddbbcd..84a2a99d3 100644 --- a/src/ABL/ablmod.F90 +++ b/src/ABL/ablmod.F90 @@ -67,54 +67,50 @@ CONTAINS !! - Finalize flux computation in psen, pevp, pwndm, ptaui, ptauj, ptaum !! !!---------------------------------------------------------------------- - INTEGER , INTENT(in ) :: kt ! time step index - REAL(wp) , INTENT(in ), DIMENSION(:,: ) :: psst ! sea-surface temperature [Celsius] - REAL(wp) , INTENT(in ), DIMENSION(:,: ) :: pssu ! sea-surface u (U-point) - REAL(wp) , INTENT(in ), DIMENSION(:,: ) :: pssv ! sea-surface v (V-point) - REAL(wp) , INTENT(in ), DIMENSION(:,: ) :: pssq ! sea-surface humidity - REAL(wp) , INTENT(in ), DIMENSION(:,:,:) :: pu_dta ! large-scale windi - REAL(wp) , INTENT(in ), DIMENSION(:,:,:) :: pv_dta ! large-scale windj - REAL(wp) , INTENT(in ), DIMENSION(:,:,:) :: pgu_dta ! large-scale hpgi - REAL(wp) , INTENT(in ), DIMENSION(:,:,:) :: pgv_dta ! large-scale hpgj - REAL(wp) , INTENT(in ), DIMENSION(:,:,:) :: pt_dta ! large-scale pot. temp. - REAL(wp) , INTENT(in ), DIMENSION(:,:,:) :: pq_dta ! large-scale humidity - REAL(wp) , INTENT(in ), DIMENSION(:,: ) :: pslp_dta ! sea-level pressure - REAL(wp) , INTENT(in ), DIMENSION(:,: ) :: pcd_du ! Cd x Du (T-point) - REAL(wp) , INTENT(inout), DIMENSION(:,: ) :: psen ! Ch x Du - REAL(wp) , INTENT(inout), DIMENSION(:,: ) :: pevp ! Ce x Du - REAL(wp) , INTENT(inout), DIMENSION(:,: ) :: pwndm ! ||uwnd - uoce|| - REAL(wp) , INTENT( out), DIMENSION(:,: ) :: plat ! latent heat flux - REAL(wp) , INTENT( out), DIMENSION(:,: ) :: ptaui ! taux - REAL(wp) , INTENT( out), DIMENSION(:,: ) :: ptauj ! tauy - REAL(wp) , INTENT( out), DIMENSION(:,: ) :: ptaum ! ||tau|| - ! + + INTEGER , INTENT(in ) :: kt ! time step index + REAL(wp) , INTENT(in ), DIMENSION(A2D(0)) :: psst ! sea-surface temperature [Celsius] + REAL(wp) , INTENT(in ), DIMENSION(A2D(nn_hls)) :: pssu ! sea-surface u (U-point) + REAL(wp) , INTENT(in ), DIMENSION(A2D(nn_hls)) :: pssv ! sea-surface v (V-point) + REAL(wp) , INTENT(in ), DIMENSION(A2D(0)) :: pssq ! sea-surface humidity + REAL(wp) , INTENT(in ), DIMENSION(A2D(0),1:jpka) :: pu_dta ! large-scale windi + REAL(wp) , INTENT(in ), DIMENSION(A2D(0),1:jpka) :: pv_dta ! large-scale windj + REAL(wp) , INTENT(in ), DIMENSION(A2D(0),1:jpka) :: pgu_dta ! large-scale hpgi + REAL(wp) , INTENT(in ), DIMENSION(A2D(0),1:jpka) :: pgv_dta ! large-scale hpgj + REAL(wp) , INTENT(in ), DIMENSION(A2D(0),1:jpka) :: pt_dta ! large-scale pot. temp. + REAL(wp) , INTENT(in ), DIMENSION(A2D(0),1:jpka) :: pq_dta ! large-scale humidity + REAL(wp) , INTENT(in ), DIMENSION(A2D(0)) :: pslp_dta ! sea-level pressure + REAL(wp) , INTENT(in ), DIMENSION(A2D(0)) :: pcd_du ! Cd x Du (T-point) + REAL(wp) , INTENT(inout), DIMENSION(A2D(0)) :: psen ! Ch x Du + REAL(wp) , INTENT(inout), DIMENSION(A2D(0)) :: pevp ! Ce x Du + REAL(wp) , INTENT(inout), DIMENSION(A2D(0)) :: pwndm ! ||uwnd - uoce|| + REAL(wp) , INTENT( out), DIMENSION(A2D(0)) :: plat ! latent heat flux + REAL(wp) , INTENT( out), DIMENSION(A2D(nn_hls)) :: ptaui ! taux + REAL(wp) , INTENT( out), DIMENSION(A2D(nn_hls)) :: ptauj ! tauy + REAL(wp) , INTENT( out), DIMENSION(A2D(0)) :: ptaum ! ||tau|| #if defined key_si3 - REAL(wp) , INTENT(in ), DIMENSION(:,: ) :: ptm_su ! ice-surface temperature [K] - REAL(wp) , INTENT(in ), DIMENSION(:,: ) :: pssu_ice ! ice-surface u (U-point) - REAL(wp) , INTENT(in ), DIMENSION(:,: ) :: pssv_ice ! ice-surface v (V-point) - REAL(wp) , INTENT(in ), DIMENSION(:,: ) :: pssq_ice ! ice-surface humidity - REAL(wp) , INTENT(in ), DIMENSION(:,: ) :: pcd_du_ice ! Cd x Du over ice (T-point) - REAL(wp) , INTENT(in ), DIMENSION(:,: ) :: psen_ice ! Ch x Du over ice (T-point) - 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 (T-point) - REAL(wp) , INTENT( out), DIMENSION(:,: ) :: ptauj_ice ! ice-surface tauy stress (T-point) + REAL(wp) , INTENT(in ), DIMENSION(A2D(0)) :: ptm_su ! ice-surface temperature [K] + REAL(wp) , INTENT(in ), DIMENSION(A2D(nn_hls)) :: pssu_ice ! ice-surface u (U-point) + REAL(wp) , INTENT(in ), DIMENSION(A2D(nn_hls)) :: pssv_ice ! ice-surface v (V-point) + REAL(wp) , INTENT(in ), DIMENSION(A2D(0)) :: pssq_ice ! ice-surface humidity + REAL(wp) , INTENT(in ), DIMENSION(A2D(0)) :: pcd_du_ice ! Cd x Du over ice (T-point) + REAL(wp) , INTENT(in ), DIMENSION(A2D(0)) :: psen_ice ! Ch x Du over ice (T-point) + REAL(wp) , INTENT(in ), DIMENSION(A2D(0)) :: pevp_ice ! Ce x Du over ice (T-point) + REAL(wp) , INTENT(in ), DIMENSION(A2D(0)) :: pwndm_ice ! ||uwnd - uice|| + REAL(wp) , INTENT(in ), DIMENSION(A2D(0)) :: pfrac_oce ! ocean fraction + REAL(wp) , INTENT( out), DIMENSION(A2D(nn_hls)) :: ptaui_ice ! ice-surface taux stress (T-point) + REAL(wp) , INTENT( out), DIMENSION(A2D(nn_hls)) :: ptauj_ice ! ice-surface tauy stress (T-point) #endif ! - REAL(wp), DIMENSION(1:jpi,1:jpj ) :: zwnd_i, zwnd_j - REAL(wp), DIMENSION(1:jpi,1:jpj ) :: zsspt - REAL(wp), DIMENSION(1:jpi,1:jpj ) :: ztabs, zpre - REAL(wp), DIMENSION(1:jpi ,2:jpka) :: zCF - ! - REAL(wp), DIMENSION(1:jpi ,1:jpka) :: z_elem_a - REAL(wp), DIMENSION(1:jpi ,1:jpka) :: z_elem_b - REAL(wp), DIMENSION(1:jpi ,1:jpka) :: z_elem_c + REAL(wp), DIMENSION(A2D(0)) :: zwnd_i, zwnd_j + REAL(wp), DIMENSION(A2D(0)) :: zsspt, ztabs, zpre + REAL(wp), DIMENSION(A1Di(0),2:jpka) :: zCF + REAL(wp), DIMENSION(A1Di(0),1:jpka) :: z_elem_a, z_elem_b, z_elem_c ! - INTEGER :: ji, jj, jk, jtra, jbak ! dummy loop indices - REAL(wp) :: zztmp, zcff, ztemp, zhumi, zcff1, zztmp1, zztmp2 - REAL(wp) :: zcff2, zfcor, zmsk, zsig, zcffu, zcffv, zzice,zzoce - LOGICAL :: SemiImp_Cor = .TRUE. + INTEGER :: ji, jj, jk, jtra, jbak ! dummy loop indices + REAL(wp) :: zztmp, zcff, ztemp, zhumi, zcff1, zztmp1, zztmp2 + REAL(wp) :: zcff2, zfcor, zmsk, zsig, zcffu, zcffv, zzice,zzoce + LOGICAL :: SemiImp_Cor = .TRUE. ! !!--------------------------------------------------------------------- ! @@ -124,12 +120,12 @@ CONTAINS WRITE(numout,*) '~~~~~~' ENDIF ! - IF( kt == nit000 ) ALLOCATE ( ustar2( 1:jpi, 1:jpj ) ) - IF( kt == nit000 ) ALLOCATE ( zrough( 1:jpi, 1:jpj ) ) + IF( kt == nit000 ) ALLOCATE ( ustar2( A2D(0) ) ) + IF( kt == nit000 ) ALLOCATE ( zrough( A2D(0) ) ) !! Compute ustar squared as Cd || Uatm-Uoce ||^2 !! needed for surface boundary condition of TKE !! pwndm contains | U10m - U_oce | (see blk_oce_1 in sbcblk) - DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) + DO_2D( 0, 0, 0, 0 ) zzoce = pCd_du (ji,jj) * pwndm (ji,jj) #if defined key_si3 zzice = pCd_du_ice(ji,jj) * pwndm_ice(ji,jj) @@ -137,39 +133,40 @@ CONTAINS #else ustar2(ji,jj) = zzoce #endif - !#LB: sorry Cdn_oce is gone: - !zrough(ji,jj) = ght_abl(2) * EXP( - vkarmn / SQRT( MAX( Cdn_oce(ji,jj), 1.e-4 ) ) ) !<-- recover the value of z0 from Cdn_oce + zsspt(ji,jj) = theta_exner( psst(ji,jj)+rt0, pslp_dta(ji,jj) ) ! potential SST [K] END_2D + !#LB: sorry Cdn_oce is gone: + !zrough(:,:) = ght_abl(2) * EXP( - vkarmn / SQRT( MAX( Cdn_oce(:,:), 1.e-4 ) ) ) !<-- recover the value of z0 from Cdn_oce zrough(:,:) = z0_from_Cd( ght_abl(2), pCd_du(:,:) / MAX( pwndm(:,:), 0.5_wp ) ) ! #LB: z0_from_Cd is define in sbc_phy.F90... ! sea surface potential temperature [K] - zsspt(:,:) = theta_exner( psst(:,:)+rt0, pslp_dta(:,:) ) + !zsspt(:,:) = theta_exner( psst(A2D(0))+rt0, pslp_dta(:,:) ) + - ! ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< ! ! 1 *** Advance TKE to time n+1 and compute Avm_abl, Avt_abl, PBLh ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< - CALL abl_zdf_tke( ) !--> Avm_abl, Avt_abl, pblh defined on (1,jpi) x (1,jpj) + CALL abl_zdf_tke( ) !--> Avm_abl, Avt_abl, pblh defined on (A1Di(0)) x (A1Dj(0)) ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< ! ! 2 *** Advance tracers to time n+1 ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< !------------- - DO jj = 1, jpj ! outer loop !--> tq_abl computed on (1:jpi) x (1:jpj) + DO_1Dj(0,0) ! outer loop !--> tq_abl computed on (A1Di(0)) x (A1Dj(0)) !------------- ! Compute matrix elements for interior points DO jk = 3, jpkam1 - DO ji = 1, jpi ! vector opt. + DO_1Di(0,0) z_elem_a( ji, jk ) = - rDt_abl * Avt_abl( ji, jj, jk-1 ) / e3w_abl( jk-1 ) ! lower-diagonal z_elem_c( ji, jk ) = - rDt_abl * Avt_abl( ji, jj, jk ) / e3w_abl( jk ) ! upper-diagonal z_elem_b( ji, jk ) = e3t_abl(jk) - z_elem_a( ji, jk ) - z_elem_c( ji, jk ) ! diagonal - END DO + END_1D END DO ! Boundary conditions - DO ji = 1, jpi ! vector opt. + DO_1Di(0,0) ! Neumann at the bottom z_elem_a( ji, 2 ) = 0._wp z_elem_c( ji, 2 ) = - rDt_abl * Avt_abl( ji, jj, 2 ) / e3w_abl( 2 ) @@ -177,19 +174,18 @@ CONTAINS z_elem_a( ji, jpka ) = - rDt_abl * Avt_abl( ji, jj, jpka ) / e3w_abl( jpka ) z_elem_c( ji, jpka ) = 0._wp z_elem_b( ji, jpka ) = e3t_abl( jpka ) - z_elem_a( ji, jpka ) - END DO + END_1D DO jtra = 1,jptq ! loop on active tracers DO jk = 3, jpkam1 - !DO ji = 2, jpim1 - DO ji = 1,jpi !!GS: to be checked if needed + DO_1Di(0,0) tq_abl( ji, jj, jk, nt_a, jtra ) = e3t_abl(jk) * tq_abl( ji, jj, jk, nt_n, jtra ) ! initialize right-hand-side - END DO + END_1D END DO IF(jtra == jp_ta) THEN - DO ji = 1,jpi ! surface boundary condition for temperature + DO_1Di(0,0) ! surface boundary condition for temperature zztmp1 = psen(ji, jj) zztmp2 = psen(ji, jj) * zsspt(ji, jj) #if defined key_si3 @@ -199,9 +195,9 @@ CONTAINS z_elem_b( ji, 2 ) = e3t_abl( 2 ) - z_elem_c( ji, 2 ) + rDt_abl * zztmp1 tq_abl ( ji, jj, 2, nt_a, jtra ) = e3t_abl( 2 ) * tq_abl ( ji, jj, 2, nt_n, jtra ) + rDt_abl * zztmp2 tq_abl ( ji, jj, jpka, nt_a, jtra ) = e3t_abl( jpka ) * tq_abl ( ji, jj, jpka, nt_n, jtra ) - END DO + END_1D ELSE ! jp_qa - DO ji = 1,jpi ! surface boundary condition for humidity + DO_1Di(0,0) ! surface boundary condition for humidity zztmp1 = pevp(ji, jj) zztmp2 = pevp(ji, jj) * pssq(ji, jj) #if defined key_si3 @@ -211,31 +207,31 @@ CONTAINS z_elem_b( ji, 2 ) = e3t_abl( 2 ) - z_elem_c( ji, 2 ) + rDt_abl * zztmp1 tq_abl ( ji, jj, 2 , nt_a, jtra ) = e3t_abl( 2 ) * tq_abl ( ji, jj, 2, nt_n, jtra ) + rDt_abl * zztmp2 tq_abl ( ji, jj, jpka, nt_a, jtra ) = e3t_abl( jpka ) * tq_abl ( ji, jj, jpka, nt_n, jtra ) - END DO + END_1D END IF !! !! Matrix inversion !! ---------------------------------------------------------- - DO ji = 1,jpi + DO_1Di(0,0) zcff = 1._wp / z_elem_b( ji, 2 ) zCF ( ji, 2 ) = - zcff * z_elem_c( ji, 2 ) tq_abl( ji, jj, 2, nt_a, jtra ) = zcff * tq_abl( ji, jj, 2, nt_a, jtra ) - END DO + END_1D DO jk = 3, jpka - DO ji = 1,jpi + DO_1Di(0,0) zcff = 1._wp / ( z_elem_b( ji, jk ) + z_elem_a( ji, jk ) * zCF( ji, jk-1 ) ) zCF(ji,jk) = - zcff * z_elem_c( ji, jk ) tq_abl(ji,jj,jk,nt_a,jtra) = zcff * ( tq_abl(ji,jj,jk ,nt_a,jtra) & & - z_elem_a(ji, jk) * tq_abl(ji,jj,jk-1,nt_a,jtra) ) - END DO + END_1D END DO !!FL at this point we could check positivity of tq_abl(:,:,:,nt_a,jp_qa) ... test to do ... DO jk = jpkam1,2,-1 - DO ji = 1,jpi + DO_1Di(0,0) tq_abl(ji,jj,jk,nt_a,jtra) = tq_abl(ji,jj,jk,nt_a,jtra) + & & zCF(ji,jk) * tq_abl(ji,jj,jk+1,nt_a,jtra) - END DO + END_1D END DO END DO !<-- loop on tracers @@ -254,7 +250,7 @@ CONTAINS !------------- ! ! Advance u_abl & v_abl to time n+1 - DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) + DO_2D( 0, 0, 0, 0 ) zcff = ( fft_abl(ji,jj) * rDt_abl )*( fft_abl(ji,jj) * rDt_abl ) ! (f dt)**2 u_abl( ji, jj, jk, nt_a ) = e3t_abl(jk) *( & @@ -269,13 +265,13 @@ CONTAINS END_2D ! !------------- - END DO ! end outer loop !<-- u_abl and v_abl are properly updated on (1:jpi) x (1:jpj) + END DO ! end outer loop !<-- u_abl and v_abl are properly updated on (A1Di(0)) x (A1Dj(0)) !------------- ! IF( ln_geos_winds ) THEN - DO jj = 1, jpj ! outer loop + DO_1Dj( 0, 0 ) ! outer loop DO jk = 1, jpka - DO ji = 1, jpi + DO_1Di( 0, 0 ) u_abl( ji, jj, jk, nt_a ) = u_abl( ji, jj, jk, nt_a ) & & - rDt_abl * e3t_abl(jk) * fft_abl(ji , jj) * pgv_dta(ji ,jj ,jk) v_abl( ji, jj, jk, nt_a ) = v_abl( ji, jj, jk, nt_a ) & @@ -286,14 +282,14 @@ CONTAINS END IF ! IF( ln_hpgls_frc ) THEN - DO jj = 1, jpj ! outer loop + DO_1Dj( 0, 0 ) ! outer loop DO jk = 1, jpka - DO ji = 1, jpi + DO_1Di( 0, 0 ) u_abl( ji, jj, jk, nt_a ) = u_abl( ji, jj, jk, nt_a ) - rDt_abl * e3t_abl(jk) * pgu_dta(ji,jj,jk) v_abl( ji, jj, jk, nt_a ) = v_abl( ji, jj, jk, nt_a ) - rDt_abl * e3t_abl(jk) * pgv_dta(ji,jj,jk) - ENDDO + END_1D ENDDO - ENDDO + END_1D END IF ELSE ! SemiImp_Cor = .FALSE. @@ -306,29 +302,25 @@ CONTAINS ! IF( MOD( kt, 2 ) == 0 ) then ! Advance u_abl & v_abl to time n+1 - DO jj = 1, jpj - DO ji = 1, jpi - zcff = fft_abl(ji,jj) * ( v_abl ( ji , jj , jk, nt_n ) - pgv_dta(ji ,jj ,jk) ) - u_abl( ji, jj, jk, nt_a ) = u_abl( ji, jj, jk, nt_n ) + rDt_abl * zcff - zcff = fft_abl(ji,jj) * ( u_abl ( ji , jj , jk, nt_a ) - pgu_dta(ji ,jj ,jk) ) - v_abl( ji, jj, jk, nt_a ) = e3t_abl(jk) *( v_abl( ji, jj, jk, nt_n ) - rDt_abl * zcff ) - u_abl( ji, jj, jk, nt_a ) = e3t_abl(jk) * u_abl( ji, jj, jk, nt_a ) - END DO - END DO + DO_2D( 0, 0, 0, 0 ) + zcff = fft_abl(ji,jj) * ( v_abl ( ji , jj , jk, nt_n ) - pgv_dta(ji ,jj ,jk) ) + u_abl( ji, jj, jk, nt_a ) = u_abl( ji, jj, jk, nt_n ) + rDt_abl * zcff + zcff = fft_abl(ji,jj) * ( u_abl ( ji , jj , jk, nt_a ) - pgu_dta(ji ,jj ,jk) ) + v_abl( ji, jj, jk, nt_a ) = e3t_abl(jk) *( v_abl( ji, jj, jk, nt_n ) - rDt_abl * zcff ) + u_abl( ji, jj, jk, nt_a ) = e3t_abl(jk) * u_abl( ji, jj, jk, nt_a ) + END_2D ELSE - DO jj = 1, jpj - DO ji = 1, jpi - zcff = fft_abl(ji,jj) * ( u_abl ( ji , jj , jk, nt_n ) - pgu_dta(ji ,jj ,jk) ) - v_abl( ji, jj, jk, nt_a ) = v_abl( ji, jj, jk, nt_n ) - rDt_abl * zcff - zcff = fft_abl(ji,jj) * ( v_abl ( ji , jj , jk, nt_a ) - pgv_dta(ji ,jj ,jk) ) - u_abl( ji, jj, jk, nt_a ) = e3t_abl(jk) *( u_abl( ji, jj, jk, nt_n ) + rDt_abl * zcff ) - v_abl( ji, jj, jk, nt_a ) = e3t_abl(jk) * v_abl( ji, jj, jk, nt_a ) - END DO - END DO + DO_2D( 0, 0, 0, 0 ) + zcff = fft_abl(ji,jj) * ( u_abl ( ji , jj , jk, nt_n ) - pgu_dta(ji ,jj ,jk) ) + v_abl( ji, jj, jk, nt_a ) = v_abl( ji, jj, jk, nt_n ) - rDt_abl * zcff + zcff = fft_abl(ji,jj) * ( v_abl ( ji , jj , jk, nt_a ) - pgv_dta(ji ,jj ,jk) ) + u_abl( ji, jj, jk, nt_a ) = e3t_abl(jk) *( u_abl( ji, jj, jk, nt_n ) + rDt_abl * zcff ) + v_abl( ji, jj, jk, nt_a ) = e3t_abl(jk) * v_abl( ji, jj, jk, nt_a ) + END_2D END IF ! !------------- - END DO ! end outer loop !<-- u_abl and v_abl are properly updated on (1:jpi) x (1:jpj) + END DO ! end outer loop !<-- u_abl and v_abl are properly updated on (A1Di(0)) x (A1Dj(0)) !------------- ENDIF ! ln_geos_winds @@ -340,18 +332,19 @@ CONTAINS ! ! Vertical diffusion for u_abl !------------- - DO jj = 1, jpj ! outer loop + DO_1Dj(0,0) ! outer loop !------------- DO jk = 3, jpkam1 - DO ji = 1, jpi + DO_1Di(0,0) z_elem_a( ji, jk ) = - rDt_abl * Avm_abl( ji, jj, jk-1 ) / e3w_abl( jk-1 ) ! lower-diagonal z_elem_c( ji, jk ) = - rDt_abl * Avm_abl( ji, jj, jk ) / e3w_abl( jk ) ! upper-diagonal z_elem_b( ji, jk ) = e3t_abl(jk) - z_elem_a( ji, jk ) - z_elem_c( ji, jk ) ! diagonal - END DO + END_1D END DO - DO ji = 2, jpi ! boundary conditions (Avm_abl and pcd_du must be available at ji=jpi) + DO_1Di(0,0) ! boundary conditions + !++ Surface boundary condition z_elem_a( ji, 2 ) = 0._wp z_elem_c( ji, 2 ) = - rDt_abl * Avm_abl( ji, jj, 2 ) / e3w_abl( 2 ) @@ -381,29 +374,29 @@ CONTAINS u_abl ( ji, jj, jpka, nt_a ) = e3t_abl( jpka ) * pu_dta(ji,jj,jk) !ENDIF - END DO + END_1D !! !! Matrix inversion !! ---------------------------------------------------------- - DO ji = 1, jpi !DO ji = 2, jpi !!GS: TBI + DO_1Di(0,0) zcff = 1._wp / z_elem_b( ji, 2 ) zCF (ji, 2 ) = - zcff * z_elem_c( ji, 2 ) u_abl (ji,jj,2,nt_a) = zcff * u_abl ( ji, jj, 2, nt_a ) - END DO + END_1D DO jk = 3, jpka - DO ji = 2, jpi + DO_1Di(0,0) zcff = 1._wp / ( z_elem_b( ji, jk ) + z_elem_a( ji, jk ) * zCF (ji, jk-1 ) ) zCF(ji,jk) = - zcff * z_elem_c( ji, jk ) u_abl(ji,jj,jk,nt_a) = zcff * ( u_abl(ji,jj,jk ,nt_a) & & - z_elem_a(ji, jk) * u_abl(ji,jj,jk-1,nt_a) ) - END DO + END_1D END DO DO jk = jpkam1,2,-1 - DO ji = 2, jpi + DO_1Di(0,0) u_abl(ji,jj,jk,nt_a) = u_abl(ji,jj,jk,nt_a) + zCF(ji,jk) * u_abl(ji,jj,jk+1,nt_a) - END DO + END_1D END DO !------------- @@ -413,18 +406,19 @@ CONTAINS ! ! Vertical diffusion for v_abl !------------- - DO jj = 2, jpj ! outer loop + DO_1Dj(0,0) ! outer loop !------------- ! DO jk = 3, jpkam1 - DO ji = 1, jpi + DO_1Di(0,0) z_elem_a( ji, jk ) = - rDt_abl * Avm_abl( ji, jj, jk-1 ) / e3w_abl( jk-1 ) ! lower-diagonal z_elem_c( ji, jk ) = - rDt_abl * Avm_abl( ji, jj, jk ) / e3w_abl( jk ) ! upper-diagonal z_elem_b( ji, jk ) = e3t_abl(jk) - z_elem_a( ji, jk ) - z_elem_c( ji, jk ) ! diagonal - END DO + END_1D END DO - DO ji = 1, jpi ! boundary conditions (Avm_abl and pcd_du must be available at jj=jpj) + DO_1Di(0,0) ! boundary conditions + !++ Surface boundary condition z_elem_a( ji, 2 ) = 0._wp z_elem_c( ji, 2 ) = - rDt_abl * Avm_abl( ji, jj, 2 ) / e3w_abl( 2 ) @@ -454,29 +448,29 @@ CONTAINS v_abl ( ji, jj, jpka, nt_a ) = e3t_abl( jpka ) * pv_dta(ji,jj,jk) !ENDIF - END DO + END_1D !! !! Matrix inversion !! ---------------------------------------------------------- - DO ji = 1, jpi + DO_1Di(0,0) zcff = 1._wp / z_elem_b( ji, 2 ) zCF (ji, 2 ) = - zcff * z_elem_c( ji, 2 ) v_abl (ji,jj,2,nt_a) = zcff * v_abl ( ji, jj, 2, nt_a ) - END DO + END_1D DO jk = 3, jpka - DO ji = 1, jpi + DO_1Di(0,0) zcff = 1._wp / ( z_elem_b( ji, jk ) + z_elem_a( ji, jk ) * zCF (ji, jk-1 ) ) zCF(ji,jk) = - zcff * z_elem_c( ji, jk ) v_abl(ji,jj,jk,nt_a) = zcff * ( v_abl(ji,jj,jk ,nt_a) & & - z_elem_a(ji, jk) * v_abl(ji,jj,jk-1,nt_a) ) - END DO + END_1D END DO DO jk = jpkam1,2,-1 - DO ji = 1, jpi + DO_1Di(0,0) v_abl(ji,jj,jk,nt_a) = v_abl(ji,jj,jk,nt_a) + zCF(ji,jk) * v_abl(ji,jj,jk+1,nt_a) - END DO + END_1D END DO ! !------------- @@ -491,7 +485,7 @@ CONTAINS !------------- DO jk = 2, jpka ! outer loop !------------- - DO_2D( nn_hls-1, nn_hls, nn_hls-1, nn_hls ) + DO_2D( 0, 0, 0, 0) zcff1 = pblh( ji, jj ) zsig = ght_abl(jk) / MAX( jp_pblh_min, MIN( jp_pblh_max, zcff1 ) ) zsig = MIN( jp_bmax , MAX( zsig, jp_bmin) ) @@ -514,7 +508,7 @@ CONTAINS !------------- DO jk = 2, jpka ! outer loop !------------- - DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) + DO_2D( 0, 0, 0, 0) zcff1 = pblh( ji, jj ) zsig = ght_abl(jk) / MAX( jp_pblh_min, MIN( jp_pblh_max, zcff1 ) ) zsig = MIN( jp_bmax , MAX( zsig, jp_bmin) ) @@ -528,21 +522,17 @@ CONTAINS tq_abl( ji, jj, jk, nt_a, jp_qa ) = (1._wp - zcff ) * tq_abl( ji, jj, jk, nt_a, jp_qa ) & & + zcff * pq_dta( ji, jj, jk ) - END_2D !------------- END DO ! end outer loop !------------- ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< - ! ! 6 *** MPI exchanges & IOM outputs + ! ! 6 *** IOM outputs ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< ! - CALL lbc_lnk( 'ablmod', u_abl(:,:,:,nt_a ), 'T', -1._wp, v_abl(:,:,:,nt_a) , 'T', -1._wp ) - !CALL lbc_lnk( 'ablmod', tq_abl(:,:,:,nt_a,jp_ta), 'T', 1._wp, tq_abl(:,:,:,nt_a,jp_qa), 'T', 1._wp, kfillmode = jpfillnothing ) ! ++++ this should not be needed... - ! #if defined key_xios ! 2D & first ABL level - IF ( iom_use("pblh" ) ) CALL iom_put ( "pblh", pblh(:,: ) ) + IF ( iom_use("pblh" ) ) CALL iom_put ( "pblh", pblh(A2D(0) ) ) IF ( iom_use("uz1_abl") ) CALL iom_put ( "uz1_abl", u_abl(:,:,2,nt_a ) ) IF ( iom_use("vz1_abl") ) CALL iom_put ( "vz1_abl", v_abl(:,:,2,nt_a ) ) IF ( iom_use("tz1_abl") ) CALL iom_put ( "tz1_abl", tq_abl(:,:,2,nt_a,jp_ta) ) @@ -588,7 +578,7 @@ CONTAINS ! ! 7 *** Finalize flux computation ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< ! - DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) + DO_2D( 0, 0, 0, 0 ) ztemp = tq_abl( ji, jj, 2, nt_a, jp_ta ) zhumi = tq_abl( ji, jj, 2, nt_a, jp_qa ) zpre( ji, jj ) = pres_temp( zhumi, pslp_dta(ji,jj), ght_abl(2), ptpot=ztemp, pta=ztabs( ji, jj ) ) @@ -599,15 +589,13 @@ CONTAINS rhoa( ji, jj ) = zcff END_2D - DO_2D( nn_hls-1, nn_hls, nn_hls-1, nn_hls ) + DO_2D( 0, 0, 0, 0 ) zwnd_i(ji,jj) = u_abl(ji ,jj,2,nt_a) - 0.5_wp * ( pssu(ji ,jj) + pssu(ji-1,jj) ) * rn_vfac zwnd_j(ji,jj) = v_abl(ji,jj ,2,nt_a) - 0.5_wp * ( pssv(ji,jj ) + pssv(ji,jj-1) ) * rn_vfac END_2D ! - CALL lbc_lnk( 'ablmod', zwnd_i(:,:) , 'T', -1.0_wp, zwnd_j(:,:) , 'T', -1.0_wp ) - ! ! ... scalar wind ( = | U10m - U_oce | ) at T-point (masked) - DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) + DO_2D( 0, 0, 0, 0 ) zcff = SQRT( zwnd_i(ji,jj) * zwnd_i(ji,jj) & & + zwnd_j(ji,jj) * zwnd_j(ji,jj) ) ! * msk_abl(ji,jj) zztmp = rhoa(ji,jj) * pcd_du(ji,jj) @@ -617,13 +605,14 @@ CONTAINS zwnd_i(ji,jj) = zztmp * zwnd_i(ji,jj) zwnd_j(ji,jj) = zztmp * zwnd_j(ji,jj) END_2D + CALL iom_put( "taum_oce", ptaum ) + ! ... 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 + DO_2D( 0, 0, 0, 0 ) + ptaui(ji,jj) = zwnd_i(ji,jj) * msk_abl(ji,jj) + ptauj(ji,jj) = zwnd_j(ji,jj) * msk_abl(ji,jj) END_2D - - CALL iom_put( "taum_oce", ptaum ) + CALL lbc_lnk( 'ablmod', ptaui(:,:) , 'T', -1.0_wp, ptauj(:,:) , 'T', -1.0_wp ) IF(sn_cfctl%l_prtctl) THEN CALL prt_ctl( tab2d_1=ptaui , clinfo1=' abl_stp: utau : ', mask1=tmask, & @@ -635,24 +624,12 @@ CONTAINS ! ------------------------------------------------------------ ! ! Wind stress relative to the moving ice ( U10m - U_ice ) ! ! ------------------------------------------------------------ ! - !DO_2D( 0, 0, 0, 0 ) - ! 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) ) & - ! & * ( 0.5_wp * ( u_abl(ji+1,jj,2,nt_a) + u_abl(ji,jj,2,nt_a) ) - pssu_ice(ji,jj) ) - ! 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) ) & - ! & * ( 0.5_wp * ( v_abl(ji,jj+1,2,nt_a) + v_abl(ji,jj,2,nt_a) ) - pssv_ice(ji,jj) ) - !END_2D - !CALL lbc_lnk( 'ablmod', ptaui_ice, 'U', -1.0_wp, ptauj_ice, 'V', -1.0_wp ) - !! - !IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab2d_1=ptaui_ice , clinfo1=' abl_stp: putaui : ' & - ! & , tab2d_2=ptauj_ice , clinfo2=' pvtaui : ' ) - - ! ------------------------------------------------------------ ! - ! Wind stress relative to the moving ice ( U10m - U_ice ) ! - ! ------------------------------------------------------------ ! - 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 ) + DO_2D(0,0,0,0) + zztmp = rhoa(ji,jj) * pCd_du_ice(ji,jj) + ptaui_ice(ji,jj) = zztmp * ( u_abl(ji,jj,2,nt_a) - 0.5_wp * ( pssu_ice(ji,jj) + pssu_ice(ji-1,jj ) ) * rn_vfac ) + ptauj_ice(ji,jj) = zztmp * ( v_abl(ji,jj,2,nt_a) - 0.5_wp * ( pssv_ice(ji,jj) + pssv_ice(ji ,jj-1) ) * rn_vfac ) END_2D + CALL lbc_lnk( 'ablmod', ptaui_ice(:,:) , 'T', -1.0_wp, ptauj_ice(:,:) , 'T', -1.0_wp ) ! IF(sn_cfctl%l_prtctl) THEN CALL prt_ctl( tab2d_1=ptaui_ice , clinfo1=' abl_stp: utau_ice : ', mask1=tmask, & @@ -690,21 +667,22 @@ CONTAINS !! - avmu, avmv : production of TKE by shear at u and v-points !! (= Kz dz[Ub] * dz[Un] ) !! --------------------------------------------------------------------- - INTEGER :: ji, jj, jk, tind, jbak, jkup, jkdwn - INTEGER, DIMENSION(1:jpi ) :: ikbl - REAL(wp) :: zcff, zcff2, ztken, zesrf, zetop, ziRic, ztv - REAL(wp) :: zdU , zdV , zcff1, zshear, zbuoy, zsig, zustar2 - REAL(wp) :: zdU2, zdV2, zbuoy1, zbuoy2 ! zbuoy for BL89 - REAL(wp) :: zwndi, zwndj - REAL(wp), DIMENSION(1:jpi, 1:jpka) :: zsh2 - REAL(wp), DIMENSION(1:jpi,1:jpj,1:jpka) :: zbn2 - REAL(wp), DIMENSION(1:jpi,1:jpka ) :: zFC, zRH, zCF - REAL(wp), DIMENSION(1:jpi,1:jpka ) :: z_elem_a - REAL(wp), DIMENSION(1:jpi,1:jpka ) :: z_elem_b - REAL(wp), DIMENSION(1:jpi,1:jpka ) :: z_elem_c - LOGICAL :: ln_Patankar = .FALSE. - LOGICAL :: ln_dumpvar = .FALSE. - LOGICAL , DIMENSION(1:jpi ) :: ln_foundl + + INTEGER :: ji, jj, jk, tind, jbak, jkup, jkdwn + INTEGER, DIMENSION(A1Di(0)) :: ikbl + REAL(wp) :: zcff, zcff2, ztken, zesrf, zetop, ziRic, ztv + REAL(wp) :: zdU , zdV , zcff1, zshear, zbuoy, zsig, zustar2 + REAL(wp) :: zdU2, zdV2, zbuoy1, zbuoy2 ! zbuoy for BL89 + REAL(wp) :: zwndi, zwndj + REAL(wp), DIMENSION(A1Di(0),1:jpka) :: zsh2 + REAL(wp), DIMENSION(A2D(0) ,1:jpka) :: zbn2 + REAL(wp), DIMENSION(A1Di(0),1:jpka) :: zFC, zRH, zCF + REAL(wp), DIMENSION(A1Di(0),1:jpka) :: z_elem_a + REAL(wp), DIMENSION(A1Di(0),1:jpka) :: z_elem_b + REAL(wp), DIMENSION(A1Di(0),1:jpka) :: z_elem_c + LOGICAL :: ln_Patankar = .FALSE. + LOGICAL :: ln_dumpvar = .FALSE. + LOGICAL , DIMENSION(A1Di(0)) :: ln_foundl ! tind = nt_n ziRic = 1._wp / rn_Ric @@ -713,33 +691,33 @@ CONTAINS ! ! Advance TKE equation to time n+1 ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< !------------- - DO jj = 1, jpj ! outer loop + DO_1Dj( 0, 0 ) ! outer loop !------------- ! ! Compute vertical shear DO jk = 2, jpkam1 - DO ji = 1, jpi + DO_1Di( 0, 0 ) zcff = 1.0_wp / e3w_abl( jk )**2 zdU = zcff* Avm_abl(ji,jj,jk) * (u_abl( ji, jj, jk+1, tind)-u_abl( ji, jj, jk , tind) )**2 zdV = zcff* Avm_abl(ji,jj,jk) * (v_abl( ji, jj, jk+1, tind)-v_abl( ji, jj, jk , tind) )**2 zsh2(ji,jk) = zdU+zdV !<-- zsh2 = Km ( ( du/dz )^2 + ( dv/dz )^2 ) - END DO + END_1D END DO ! ! Compute brunt-vaisala frequency DO jk = 2, jpkam1 - DO ji = 1,jpi + DO_1Di( 0, 0 ) zcff = grav * itvref / e3w_abl( jk ) zcff1 = tq_abl( ji, jj, jk+1, tind, jp_ta) - tq_abl( ji, jj, jk , tind, jp_ta) zcff2 = tq_abl( ji, jj, jk+1, tind, jp_ta) * tq_abl( ji, jj, jk+1, tind, jp_qa) & & - tq_abl( ji, jj, jk , tind, jp_ta) * tq_abl( ji, jj, jk , tind, jp_qa) zbn2(ji,jj,jk) = zcff * ( zcff1 + rctv0 * zcff2 ) !<-- zbn2 defined on (2,jpi) - END DO + END_1D END DO ! ! Terms for the tridiagonal problem DO jk = 2, jpkam1 - DO ji = 1, jpi + DO_1Di( 0, 0 ) zshear = zsh2( ji, jk ) ! zsh2 is already multiplied by Avm_abl at this point zsh2(ji,jk) = zsh2( ji, jk ) / Avm_abl( ji, jj, jk ) ! reformulate zsh2 as a 'true' vertical shear for PBLH computation zbuoy = - Avt_abl( ji, jj, jk ) * zbn2( ji, jj, jk ) @@ -756,10 +734,10 @@ CONTAINS & - e3w_abl(jk) * rDt_abl * zbuoy tke_abl( ji, jj, jk, nt_a ) = e3w_abl(jk) * ( tke_abl( ji, jj, jk, nt_n ) + rDt_abl * zshear ) ! right-hand-side END IF - END DO + END_1D END DO - DO ji = 1,jpi ! vector opt. + DO_1Di( 0, 0 ) zesrf = MAX( rn_Esfc * ustar2(ji,jj), tke_min ) zetop = tke_min @@ -788,44 +766,41 @@ CONTAINS !! !! Matrix inversion !! ---------------------------------------------------------- - DO ji = 1,jpi + DO_1Di( 0, 0 ) zcff = 1._wp / z_elem_b( ji, 1 ) zCF (ji, 1 ) = - zcff * z_elem_c( ji, 1 ) tke_abl(ji,jj,1,nt_a) = zcff * tke_abl ( ji, jj, 1, nt_a ) - END DO + END_1D DO jk = 2, jpka - DO ji = 1,jpi + DO_1Di( 0, 0 ) zcff = 1._wp / ( z_elem_b( ji, jk ) + z_elem_a( ji, jk ) * zCF(ji, jk-1 ) ) zCF(ji,jk) = - zcff * z_elem_c( ji, jk ) tke_abl(ji,jj,jk,nt_a) = zcff * ( tke_abl(ji,jj,jk ,nt_a) & & - z_elem_a(ji, jk) * tke_abl(ji,jj,jk-1,nt_a) ) - END DO + END_1D END DO DO jk = jpkam1,1,-1 - DO ji = 1,jpi + DO_1Di( 0, 0 ) tke_abl(ji,jj,jk,nt_a) = tke_abl(ji,jj,jk,nt_a) + zCF(ji,jk) * tke_abl(ji,jj,jk+1,nt_a) - END DO + END_1D END DO -!!FL should not be needed because of Patankar procedure - tke_abl(2:jpi,jj,1:jpka,nt_a) = MAX( tke_abl(2:jpi,jj,1:jpka,nt_a), tke_min ) + !!FL should not be needed because of Patankar procedure + tke_abl(A2D(0),1:jpka,nt_a) = MAX( tke_abl(A2D(0),1:jpka,nt_a), tke_min ) !! !! Diagnose PBL height !! ---------------------------------------------------------- - - - ! ! arrays zRH, zFC and zCF are available at this point ! and zFC(:, 1 ) = 0. ! diagnose PBL height based on zsh2 and zbn2 zFC ( : ,1) = 0._wp - ikbl( 1:jpi ) = 0 + ikbl( A1Di(0) ) = 0 DO jk = 2,jpka - DO ji = 1, jpi + DO_1Di( 0, 0 ) zcff = ghw_abl( jk-1 ) zcff1 = zcff / ( zcff + rn_epssfc * pblh ( ji, jj ) ) zcff = ghw_abl( jk ) @@ -837,11 +812,11 @@ CONTAINS - rn_Cek * ( fft_abl( ji, jj ) * fft_abl( ji, jj ) ) ) & & ) IF( ikbl(ji) == 0 .and. zFC( ji, jk ).lt.0._wp ) ikbl(ji)=jk - END DO + END_1D END DO ! ! finalize the computation of the PBL height - DO ji = 1, jpi + DO_1Di( 0, 0 ) jk = ikbl(ji) IF( jk > 2 ) THEN ! linear interpolation to get subgrid value of pblh pblh( ji, jj ) = ( ghw_abl( jk-1 ) * zFC( ji, jk ) & @@ -852,16 +827,16 @@ CONTAINS ELSE pblh( ji, jj ) = ghw_abl(jpka) END IF - END DO + END_1D !------------- END DO !------------- ! ! Optional : could add pblh smoothing if pblh is noisy horizontally ... IF(ln_smth_pblh) THEN - CALL lbc_lnk( 'ablmod', pblh, 'T', 1.0_wp) !, kfillmode = jpfillnothing) - CALL smooth_pblh( pblh, msk_abl ) - CALL lbc_lnk( 'ablmod', pblh, 'T', 1.0_wp) !, kfillmode = jpfillnothing) + CALL lbc_lnk( 'ablmod', pblh, 'T', 1.0_wp ) + CALL smooth_pblh( pblh, tmask(A2D(0),1) ) + CALL lbc_lnk( 'ablmod', pblh, 'T', 1.0_wp ) ENDIF ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< ! ! Diagnostic mixing length computation @@ -872,55 +847,55 @@ CONTAINS CASE ( 0 ) ! Deardroff 80 length-scale bounded by the distance to surface and bottom # define zlup zRH # define zldw zFC - DO jj = 1, jpj ! outer loop + DO_1Dj(0,0) ! outer loop ! - DO ji = 1, jpi + DO_1Di(0,0) mxld_abl( ji, jj, 1 ) = mxl_min mxld_abl( ji, jj, jpka ) = mxl_min mxlm_abl( ji, jj, 1 ) = mxl_min mxlm_abl( ji, jj, jpka ) = mxl_min zldw ( ji, 1 ) = zrough(ji,jj) * rn_Lsfc zlup ( ji, jpka ) = mxl_min - END DO + END_1D ! DO jk = 2, jpkam1 - DO ji = 1, jpi + DO_1Di(0,0) zbuoy = MAX( zbn2(ji, jj, jk), rsmall ) mxlm_abl( ji, jj, jk ) = MAX( mxl_min, & & SQRT( 2._wp * tke_abl( ji, jj, jk, nt_a ) / zbuoy ) ) - END DO + END_1D END DO ! ! Limit mxl DO jk = jpkam1,1,-1 - DO ji = 1, jpi + DO_1Di(0,0) zlup(ji,jk) = MIN( zlup(ji,jk+1) + (ghw_abl(jk+1)-ghw_abl(jk)) , mxlm_abl(ji, jj, jk) ) - END DO + END_1D END DO ! DO jk = 2, jpka - DO ji = 1, jpi + DO_1Di(0,0) zldw(ji,jk) = MIN( zldw(ji,jk-1) + (ghw_abl(jk)-ghw_abl(jk-1)) , mxlm_abl(ji, jj, jk) ) - END DO + END_1D END DO ! ! DO jk = 1, jpka -! DO ji = 1, jpi +! DO_1Di(0,0) ! mxlm_abl( ji, jj, jk ) = SQRT( zldw( ji, jk ) * zlup( ji, jk ) ) ! mxld_abl( ji, jj, jk ) = MIN ( zldw( ji, jk ), zlup( ji, jk ) ) -! END DO +! END_1D ! END DO ! DO jk = 1, jpka - DO ji = 1, jpi + DO_1Di(0,0) ! zcff = 2.*SQRT(2.)*( zldw( ji, jk )**(-2._wp/3._wp) + zlup( ji, jk )**(-2._wp/3._wp) )**(-3._wp/2._wp) zcff = SQRT( zldw( ji, jk ) * zlup( ji, jk ) ) mxlm_abl( ji, jj, jk ) = MAX( zcff, mxl_min ) mxld_abl( ji, jj, jk ) = MAX( MIN( zldw( ji, jk ), zlup( ji, jk ) ), mxl_min ) - END DO + END_1D END DO ! - END DO + END_1D ! outer loop # undef zlup # undef zldw ! @@ -928,18 +903,18 @@ CONTAINS CASE ( 1 ) ! Modified Deardroff 80 length-scale bounded by the distance to surface and bottom # define zlup zRH # define zldw zFC - DO jj = 1, jpj ! outer loop + DO_1Dj( 0, 0 ) ! outer loop ! DO jk = 2, jpkam1 - DO ji = 1,jpi - zcff = 1.0_wp / e3w_abl( jk )**2 + DO_1Di( 0, 0 ) + zcff = 1.0_wp / e3w_abl( jk )**2 zdU = zcff* (u_abl( ji, jj, jk+1, tind)-u_abl( ji, jj, jk , tind) )**2 zdV = zcff* (v_abl( ji, jj, jk+1, tind)-v_abl( ji, jj, jk , tind) )**2 zsh2(ji,jk) = SQRT(zdU+zdV) !<-- zsh2 = SQRT ( ( du/dz )^2 + ( dv/dz )^2 ) - ENDDO - ENDDO - ! - DO ji = 1, jpi + END_1D + END DO + ! + DO_1Di( 0, 0 ) zcff = zrough(ji,jj) * rn_Lsfc mxld_abl ( ji, jj, 1 ) = zcff mxld_abl ( ji, jj, jpka ) = mxl_min @@ -947,42 +922,42 @@ CONTAINS mxlm_abl ( ji, jj, jpka ) = mxl_min zldw ( ji, 1 ) = zcff zlup ( ji, jpka ) = mxl_min - END DO + END_1D ! DO jk = 2, jpkam1 - DO ji = 1, jpi + DO_1Di( 0, 0 ) zbuoy = MAX( zbn2(ji, jj, jk), rsmall ) zcff = 2.0_wp*SQRT(tke_abl( ji, jj, jk, nt_a )) / ( rn_Rod*zsh2(ji,jk) & & + SQRT(rn_Rod*rn_Rod*zsh2(ji,jk)*zsh2(ji,jk)+2.0_wp*zbuoy ) ) mxlm_abl( ji, jj, jk ) = MAX( mxl_min, zcff ) - END DO + END_1D END DO ! ! Limit mxl DO jk = jpkam1,1,-1 - DO ji = 1, jpi + DO_1Di( 0, 0 ) zlup(ji,jk) = MIN( zlup(ji,jk+1) + (ghw_abl(jk+1)-ghw_abl(jk)) , mxlm_abl(ji, jj, jk) ) - END DO + END_1D END DO ! DO jk = 2, jpka - DO ji = 1, jpi + DO_1Di( 0, 0 ) zldw(ji,jk) = MIN( zldw(ji,jk-1) + (ghw_abl(jk)-ghw_abl(jk-1)) , mxlm_abl(ji, jj, jk) ) - END DO + END_1D END DO ! DO jk = 1, jpka - DO ji = 1, jpi + DO_1Di( 0, 0 ) !mxlm_abl( ji, jj, jk ) = SQRT( zldw( ji, jk ) * zlup( ji, jk ) ) !zcff = 2.*SQRT(2.)*( zldw( ji, jk )**(-2._wp/3._wp) + zlup( ji, jk )**(-2._wp/3._wp) )**(-3._wp/2._wp) zcff = SQRT( zldw( ji, jk ) * zlup( ji, jk ) ) mxlm_abl( ji, jj, jk ) = MAX( zcff, mxl_min ) !mxld_abl( ji, jj, jk ) = MIN( zldw( ji, jk ), zlup( ji, jk ) ) mxld_abl( ji, jj, jk ) = MAX( MIN( zldw( ji, jk ), zlup( ji, jk ) ), mxl_min ) - END DO + END_1D END DO - ! - END DO + ! + END_1D ! outer loop # undef zlup # undef zldw ! @@ -992,35 +967,35 @@ CONTAINS # define zldw zFC ! zCF is used for matrix inversion ! - DO jj = 1, jpj ! outer loop + DO_1Dj( 0, 0 ) ! outer loop - DO ji = 1, jpi + DO_1Di( 0, 0 ) zcff = zrough(ji,jj) * rn_Lsfc zlup( ji, 1 ) = zcff zldw( ji, 1 ) = zcff zlup( ji, jpka ) = mxl_min zldw( ji, jpka ) = mxl_min - END DO + END_1D DO jk = 2,jpka-1 - DO ji = 1, jpi + DO_1Di( 0, 0 ) zlup(ji,jk) = ghw_abl(jpka) - ghw_abl(jk) zldw(ji,jk) = ghw_abl(jk ) - ghw_abl( 1) - END DO + END_1D END DO !! !! BL89 search for lup !! ---------------------------------------------------------- DO jk=2,jpka-1 ! - DO ji = 1, jpi + DO_1Di( 0, 0 ) zCF(ji,1:jpka) = 0._wp zCF(ji, jk ) = - tke_abl( ji, jj, jk, nt_a ) ln_foundl(ji ) = .false. - END DO + END_1D ! DO jkup=jk+1,jpka-1 - DO ji = 1, jpi + DO_1Di( 0, 0 ) zbuoy1 = MAX( zbn2(ji,jj,jkup ), rsmall ) zbuoy2 = MAX( zbn2(ji,jj,jkup-1), rsmall ) zCF (ji,jkup) = zCF (ji,jkup-1) + 0.5_wp * e3t_abl(jkup) * & @@ -1035,7 +1010,7 @@ CONTAINS zlup(ji,jk) = ghw_abl(jkup ) - ghw_abl(jk) ln_foundl(ji) = .true. END IF - END DO + END_1D END DO ! END DO @@ -1044,14 +1019,14 @@ CONTAINS !! ---------------------------------------------------------- DO jk=2,jpka-1 ! - DO ji = 1, jpi + DO_1Di( 0, 0 ) zCF(ji,1:jpka) = 0._wp zCF(ji, jk ) = - tke_abl( ji, jj, jk, nt_a ) ln_foundl(ji ) = .false. - END DO + END_1D ! DO jkdwn=jk-1,1,-1 - DO ji = 1, jpi + DO_1Di( 0, 0 ) zbuoy1 = MAX( zbn2(ji,jj,jkdwn+1), rsmall ) zbuoy2 = MAX( zbn2(ji,jj,jkdwn ), rsmall ) zCF (ji,jkdwn) = zCF (ji,jkdwn+1) + 0.5_wp * e3t_abl(jkdwn+1) & @@ -1066,21 +1041,21 @@ CONTAINS zldw(ji,jk) = ghw_abl(jk) - ghw_abl(jkdwn ) ln_foundl(ji) = .true. END IF - END DO + END_1D END DO ! END DO DO jk = 1, jpka - DO ji = 1, jpi + DO_1Di( 0, 0 ) !zcff = 2.*SQRT(2.)*( zldw( ji, jk )**(-2._wp/3._wp) + zlup( ji, jk )**(-2._wp/3._wp) )**(-3._wp/2._wp) zcff = SQRT( zldw( ji, jk ) * zlup( ji, jk ) ) mxlm_abl( ji, jj, jk ) = MAX( zcff, mxl_min ) mxld_abl( ji, jj, jk ) = MAX( MIN( zldw( ji, jk ), zlup( ji, jk ) ), mxl_min ) - END DO + END_1D END DO - END DO + END_1D ! outer loop # undef zlup # undef zldw ! @@ -1090,46 +1065,46 @@ CONTAINS # define zldw zFC ! zCF is used for matrix inversion ! - DO jj = 1, jpj ! outer loop + DO_1Dj( 0, 0 ) ! outer loop ! DO jk = 2, jpkam1 - DO ji = 1,jpi + DO_1Di( 0, 0 ) zcff = 1.0_wp / e3w_abl( jk )**2 zdU = zcff* (u_abl( ji, jj, jk+1, tind)-u_abl( ji, jj, jk , tind) )**2 zdV = zcff* (v_abl( ji, jj, jk+1, tind)-v_abl( ji, jj, jk , tind) )**2 zsh2(ji,jk) = SQRT(zdU+zdV) !<-- zsh2 = SQRT ( ( du/dz )^2 + ( dv/dz )^2 ) - ENDDO - ENDDO + END_1D + ENDDO zsh2(:, 1) = zsh2( :, 2) zsh2(:, jpka) = zsh2( :, jpkam1) - DO ji = 1, jpi - zcff = zrough(ji,jj) * rn_Lsfc - zlup( ji, 1 ) = zcff - zldw( ji, 1 ) = zcff - zlup( ji, jpka ) = mxl_min - zldw( ji, jpka ) = mxl_min - END DO + DO_1Di( 0, 0 ) + zcff = zrough(ji,jj) * rn_Lsfc + zlup( ji, 1 ) = zcff + zldw( ji, 1 ) = zcff + zlup( ji, jpka ) = mxl_min + zldw( ji, jpka ) = mxl_min + END_1D DO jk = 2,jpka-1 - DO ji = 1, jpi + DO_1Di( 0, 0 ) zlup(ji,jk) = ghw_abl(jpka) - ghw_abl(jk) zldw(ji,jk) = ghw_abl(jk ) - ghw_abl( 1) - END DO + END_1D END DO !! !! BL89 search for lup !! ---------------------------------------------------------- DO jk=2,jpka-1 ! - DO ji = 1, jpi + DO_1Di( 0, 0 ) zCF(ji,1:jpka) = 0._wp zCF(ji, jk ) = - tke_abl( ji, jj, jk, nt_a ) ln_foundl(ji ) = .false. - END DO + END_1D ! DO jkup=jk+1,jpka-1 - DO ji = 1, jpi + DO_1Di( 0, 0 ) zbuoy1 = MAX( zbn2(ji,jj,jkup ), rsmall ) zbuoy2 = MAX( zbn2(ji,jj,jkup-1), rsmall ) zCF (ji,jkup) = zCF (ji,jkup-1) + 0.5_wp * e3t_abl(jkup) * & @@ -1148,7 +1123,7 @@ CONTAINS zlup(ji,jk) = ghw_abl(jkup ) - ghw_abl(jk) ln_foundl(ji) = .true. END IF - END DO + END_1D END DO ! END DO @@ -1157,14 +1132,14 @@ CONTAINS !! ---------------------------------------------------------- DO jk=2,jpka-1 ! - DO ji = 1, jpi + DO_1Di( 0, 0 ) zCF(ji,1:jpka) = 0._wp zCF(ji, jk ) = - tke_abl( ji, jj, jk, nt_a ) ln_foundl(ji ) = .false. - END DO + END_1D ! DO jkdwn=jk-1,1,-1 - DO ji = 1, jpi + DO_1Di( 0, 0 ) zbuoy1 = MAX( zbn2(ji,jj,jkdwn+1), rsmall ) zbuoy2 = MAX( zbn2(ji,jj,jkdwn ), rsmall ) zCF (ji,jkdwn) = zCF (ji,jkdwn+1) + 0.5_wp * e3t_abl(jkdwn+1) & @@ -1183,18 +1158,18 @@ CONTAINS zldw(ji,jk) = ghw_abl(jk) - ghw_abl(jkdwn ) ln_foundl(ji) = .true. END IF - END DO + END_1D END DO ! END DO DO jk = 1, jpka - DO ji = 1, jpi + DO_1Di( 0, 0 ) !zcff = 2.*SQRT(2.)*( zldw( ji, jk )**(-2._wp/3._wp) + zlup( ji, jk )**(-2._wp/3._wp) )**(-3._wp/2._wp) zcff = SQRT( zldw( ji, jk ) * zlup( ji, jk ) ) mxlm_abl( ji, jj, jk ) = MAX( zcff, mxl_min ) mxld_abl( ji, jj, jk ) = MAX( MIN( zldw( ji, jk ), zlup( ji, jk ) ), mxl_min ) - END DO + END_1D END DO END DO @@ -1208,10 +1183,10 @@ CONTAINS ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< !------------- - DO jj = 1, jpj ! outer loop + DO_1Dj( 0, 0 ) ! outer loop !------------- DO jk = 1, jpka - DO ji = 1, jpi ! vector opt. + DO_1Di( 0, 0 ) zcff = MAX( rn_phimax, rn_Ric * mxlm_abl( ji, jj, jk ) * mxld_abl( ji, jj, jk ) & & * MAX( zbn2(ji, jj, jk), rsmall ) / tke_abl( ji, jj, jk, nt_a ) ) zcff2 = 1. / ( 1. + zcff ) !<-- phi_z(z) @@ -1219,10 +1194,10 @@ CONTAINS !!FL: MAX function probably useless because of the definition of mxl_min Avm_abl( ji, jj, jk ) = MAX( rn_Cm * zcff , avm_bak ) Avt_abl( ji, jj, jk ) = MAX( rn_Ct * zcff * zcff2 , avt_bak ) - END DO + END_1D END DO !------------- - END DO + END_1D ! outer loop !------------- !--------------------------------------------------------------------------------------------------- @@ -1240,42 +1215,43 @@ CONTAINS !! ** Purpose : 2D Hanning filter on atmospheric PBL height !! !! --------------------------------------------------------------------- - REAL(wp), DIMENSION(jpi,jpj), INTENT(in ) :: msk - REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: pvar2d - INTEGER :: ji,jj - REAL(wp) :: smth_a, smth_b - REAL(wp), DIMENSION(jpi,jpj) :: zdX,zdY,zFX,zFY - REAL(wp) :: zumsk,zvmsk - !! + + REAL(wp), DIMENSION(A2D(1)), INTENT(in ) :: msk + REAL(wp), DIMENSION(A2D(1)), INTENT(inout) :: pvar2d + INTEGER :: ji,jj + REAL(wp) :: smth_a, smth_b + REAL(wp), DIMENSION(A2D(1)) :: zdX,zdY,zFX,zFY + REAL(wp) :: zumsk,zvmsk + !!========================================================= !! !! Hanning filter smth_a = 1._wp / 8._wp smth_b = 1._wp / 4._wp ! - DO_2D( nn_hls, nn_hls-1, nn_hls, nn_hls ) + DO_2D( 1, 0, 1, 1 ) zumsk = msk(ji,jj) * msk(ji+1,jj) zdX ( ji, jj ) = ( pvar2d( ji+1,jj ) - pvar2d( ji ,jj ) ) * zumsk END_2D - DO_2D( nn_hls, nn_hls, nn_hls, nn_hls-1 ) + DO_2D( 1, 1, 1, 0 ) zvmsk = msk(ji,jj) * msk(ji,jj+1) zdY ( ji, jj ) = ( pvar2d( ji, jj+1 ) - pvar2d( ji ,jj ) ) * zvmsk END_2D - DO_2D( nn_hls-1, nn_hls-1, nn_hls, nn_hls-1 ) + DO_2D( 0, 0, 1, 0 ) zFY ( ji, jj ) = zdY ( ji, jj ) & & + smth_a* ( (zdX ( ji, jj+1 ) - zdX( ji-1, jj+1 )) & & - (zdX ( ji, jj ) - zdX( ji-1, jj )) ) END_2D - DO_2D( nn_hls, nn_hls-1, nn_hls-1, nn_hls-1 ) + DO_2D( 1, 0, 0, 0 ) zFX( ji, jj ) = zdX( ji, jj ) & & + smth_a*( (zdY( ji+1, jj ) - zdY( ji+1, jj-1)) & & - (zdY( ji , jj ) - zdY( ji , jj-1)) ) END_2D - DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) + DO_2D( 0, 0, 0, 0 ) pvar2d( ji ,jj ) = pvar2d( ji ,jj ) & & + msk(ji,jj) * smth_b * ( & & zFX( ji, jj ) - zFX( ji-1, jj ) & diff --git a/src/ABL/ablrst.F90 b/src/ABL/ablrst.F90 index 925234b3e..de4172974 100644 --- a/src/ABL/ablrst.F90 +++ b/src/ABL/ablrst.F90 @@ -202,12 +202,12 @@ CONTAINS ! --- mandatory fields --- ! CALL iom_get( numrar, jpdom_auto, 'u_abl', u_abl(:,:,:,nt_n ), cd_type = 'T', psgn = -1._wp ) CALL iom_get( numrar, jpdom_auto, 'v_abl', v_abl(:,:,:,nt_n ), cd_type = 'T', psgn = -1._wp ) - CALL iom_get( numrar, jpdom_auto, 't_abl', tq_abl(:,:,:,nt_n,jp_ta), kfill = jpfillcopy ) + CALL iom_get( numrar, jpdom_auto, 't_abl', tq_abl(:,:,:,nt_n,jp_ta) ) !, kfill = jpfillcopy ) CALL iom_get( numrar, jpdom_auto, 'q_abl', tq_abl(:,:,:,nt_n,jp_qa) ) CALL iom_get( numrar, jpdom_auto, 'tke_abl', tke_abl(:,:,:,nt_n ) ) - CALL iom_get( numrar, jpdom_auto, 'avm_abl', avm_abl(:,:,: ), kfill = jpfillcopy ) + CALL iom_get( numrar, jpdom_auto, 'avm_abl', avm_abl(:,:,: ) ) !, kfill = jpfillcopy ) CALL iom_get( numrar, jpdom_auto, 'avt_abl', avt_abl(:,:,: ) ) - CALL iom_get( numrar, jpdom_auto,'mxld_abl',mxld_abl(:,:,: ), kfill = jpfillcopy ) + CALL iom_get( numrar, jpdom_auto,'mxld_abl',mxld_abl(:,:,: ) ) !, kfill = jpfillcopy ) CALL iom_get( numrar, jpdom_auto, 'pblh', pblh(:,: ), kfill = jpfillcopy ) IF(.NOT.lrxios) CALL iom_delay_rst( 'READ', 'ABL', numrar ) ! read only abl delayed global communication variables diff --git a/src/ABL/sbcabl.F90 b/src/ABL/sbcabl.F90 index d0eac5483..cc34c4a19 100644 --- a/src/ABL/sbcabl.F90 +++ b/src/ABL/sbcabl.F90 @@ -261,7 +261,7 @@ CONTAINS rest_eq(:,:) = 1._wp END IF ! T-mask - msk_abl(:,:) = tmask(:,:,1) + msk_abl(:,:) = tmask(A2D(0),1) !!------------------------------------------------------------------------------------------- @@ -320,9 +320,11 @@ CONTAINS !!--------------------------------------------------------------------- INTEGER , INTENT(in) :: kt ! ocean time step !! - REAL(wp), DIMENSION(jpi,jpj) :: zssq, zcd_du, zsen, zlat, zevp + !REAL(wp), DIMENSION(jpi,jpj) :: zssq, zcd_du, zsen, zlat, zevp + REAL(wp), DIMENSION(A2D(0)) :: zssq, zcd_du, zsen, zlat, zevp #if defined key_si3 - REAL(wp), DIMENSION(jpi,jpj) :: zssqi, zcd_dui, zseni, zevpi + !REAL(wp), DIMENSION(jpi,jpj) :: zssqi, zcd_dui, zseni, zevpi + REAL(wp), DIMENSION(A2D(0)) :: zssqi, zcd_dui, zseni, zevpi #endif INTEGER :: jbak, jbak_dta, ji, jj !!--------------------------------------------------------------------- @@ -357,7 +359,7 @@ CONTAINS !! 3 - Advance ABL variables from now (n) to after (n+1) !!------------------------------------------------------------------------------------------- - CALL abl_stp( kt, tsk_m, ssu_m, ssv_m, zssq, & ! <<= in + CALL abl_stp( kt, tsk_m(A2D(0)), ssu_m, ssv_m, zssq, & ! <<= in & sf(jp_wndi)%fnow(:,:,:), sf(jp_wndj)%fnow(:,:,:), & ! <<= in & sf(jp_tair)%fnow(:,:,:), sf(jp_humi)%fnow(:,:,:), & ! <<= in & sf(jp_slp )%fnow(:,:,1), & ! <<= in @@ -366,7 +368,7 @@ CONTAINS & zlat, wndm, utau, vtau, taum & ! =>> out #if defined key_si3 & , tm_su, u_ice, v_ice, zssqi, zcd_dui & ! <<= in - & , zseni, zevpi, wndm_ice, ato_i & ! <<= in + & , zseni, zevpi, wndm_ice, ato_i(A2D(0)) & ! <<= in & , utau_ice, vtau_ice & ! =>> out #endif & ) diff --git a/src/OCE/SBC/sbcblk.F90 b/src/OCE/SBC/sbcblk.F90 index ba5255c90..118951c9a 100644 --- a/src/OCE/SBC/sbcblk.F90 +++ b/src/OCE/SBC/sbcblk.F90 @@ -559,6 +559,7 @@ CONTAINS ! Potential temperature of air at z=rn_zqt (most reanalysis products provide absolute temp., not potential temp.) IF( ln_tair_pot ) THEN ! temperature read into file is already potential temperature, do nothing... + IF((kt==nit000).AND.lwp) WRITE(numout,*) ' *** sbc_blk() => air temperature already converted to POTENTIAL!' theta_air_zt(:,:) = sf(jp_tair )%fnow(:,:,1) ELSE ! temperature read into file is ABSOLUTE temperature (that's the case for ECMWF products for example...) diff --git a/src/OCE/ZDF/zdfgls.F90 b/src/OCE/ZDF/zdfgls.F90 index f44ccc814..798fdaf49 100644 --- a/src/OCE/ZDF/zdfgls.F90 +++ b/src/OCE/ZDF/zdfgls.F90 @@ -973,7 +973,8 @@ CONTAINS CALL ctl_stop( 'zdf_gls_init: wrong value for nn_mxlice, should be 0,1,2,3 ') END SELECT IF ( (nn_mxlice>0).AND.(nn_ice==0) ) THEN - CALL ctl_stop( 'zdf_gls_init: with no ice at all, nn_mxlice must be 0 ') + CALL ctl_warn( 'zdf_gls_init: with no ice at all, nn_mxlice is set to 0 ') + nn_mxlice = 0 ELSEIF ( (nn_mxlice>1).AND.(nn_ice==1) ) THEN CALL ctl_stop( 'zdf_gls_init: with no ice model, nn_mxlice must be 0 or 1') ENDIF diff --git a/src/OCE/ZDF/zdftke.F90 b/src/OCE/ZDF/zdftke.F90 index 60371fbcd..0661ea22d 100644 --- a/src/OCE/ZDF/zdftke.F90 +++ b/src/OCE/ZDF/zdftke.F90 @@ -761,7 +761,8 @@ CONTAINS CALL ctl_stop( 'zdf_tke_init: wrong value for nn_mxlice, should be 0,1,2,3 or 4') END SELECT IF ( (nn_mxlice>0).AND.(nn_ice==0) ) THEN - CALL ctl_stop( 'zdf_tke_init: with no ice at all, nn_mxlice must be 0 ') + CALL ctl_warn( 'zdf_tke_init: with no ice at all, nn_mxlice is set to 0 ') + nn_mxlice = 0 ELSEIF ( (nn_mxlice>1).AND.(nn_ice==1) ) THEN CALL ctl_stop( 'zdf_tke_init: with no ice model, nn_mxlice must be 0 or 1') ENDIF -- GitLab