New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 12928 for NEMO/branches/2019/dev_r11078_OSMOSIS_IMMERSE_Nurser/src/OCE/ZDF/zdfiwm.F90 – NEMO

Ignore:
Timestamp:
2020-05-14T21:46:00+02:00 (4 years ago)
Author:
smueller
Message:

Synchronizing with /NEMO/trunk@12925 (ticket #2170)

Location:
NEMO/branches/2019/dev_r11078_OSMOSIS_IMMERSE_Nurser
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2019/dev_r11078_OSMOSIS_IMMERSE_Nurser

    • Property svn:externals
      •  

        old new  
        66^/vendors/FCM@HEAD            ext/FCM 
        77^/vendors/IOIPSL@HEAD         ext/IOIPSL 
         8 
         9# SETTE 
         10^/utils/CI/sette@HEAD         sette 
  • NEMO/branches/2019/dev_r11078_OSMOSIS_IMMERSE_Nurser/src/OCE/ZDF/zdfiwm.F90

    r12178 r12928  
    2323   USE phycst         ! physical constants 
    2424   ! 
     25   USE fldread        ! field read 
    2526   USE prtctl         ! Print control 
    2627   USE in_out_manager ! I/O manager 
     
    4950 
    5051   !! * Substitutions 
    51 #  include "vectopt_loop_substitute.h90" 
     52#  include "do_loop_substitute.h90" 
    5253   !!---------------------------------------------------------------------- 
    5354   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    6970 
    7071 
    71    SUBROUTINE zdf_iwm( kt, p_avm, p_avt, p_avs ) 
     72   SUBROUTINE zdf_iwm( kt, Kmm, p_avm, p_avt, p_avs ) 
    7273      !!---------------------------------------------------------------------- 
    7374      !!                  ***  ROUTINE zdf_iwm  *** 
     
    8788      !!              This is divided into three components: 
    8889      !!                 1. Bottom-intensified low-mode dissipation at critical slopes 
    89       !!                     zemx_iwm(z) = ( ecri_iwm / rau0 ) * EXP( -(H-z)/hcri_iwm ) 
     90      !!                     zemx_iwm(z) = ( ecri_iwm / rho0 ) * EXP( -(H-z)/hcri_iwm ) 
    9091      !!                                   / ( 1. - EXP( - H/hcri_iwm ) ) * hcri_iwm 
    9192      !!              where hcri_iwm is the characteristic length scale of the bottom  
    9293      !!              intensification, ecri_iwm a map of available power, and H the ocean depth. 
    9394      !!                 2. Pycnocline-intensified low-mode dissipation 
    94       !!                     zemx_iwm(z) = ( epyc_iwm / rau0 ) * ( sqrt(rn2(z))^nn_zpyc ) 
     95      !!                     zemx_iwm(z) = ( epyc_iwm / rho0 ) * ( sqrt(rn2(z))^nn_zpyc ) 
    9596      !!                                   / SUM( sqrt(rn2(z))^nn_zpyc * e3w(z) ) 
    9697      !!              where epyc_iwm is a map of available power, and nn_zpyc 
     
    9899      !!              energy dissipation. 
    99100      !!                 3. WKB-height dependent high mode dissipation 
    100       !!                     zemx_iwm(z) = ( ebot_iwm / rau0 ) * rn2(z) * EXP(-z_wkb(z)/hbot_iwm) 
     101      !!                     zemx_iwm(z) = ( ebot_iwm / rho0 ) * rn2(z) * EXP(-z_wkb(z)/hbot_iwm) 
    101102      !!                                   / SUM( rn2(z) * EXP(-z_wkb(z)/hbot_iwm) * e3w(z) ) 
    102103      !!              where hbot_iwm is the characteristic length scale of the WKB bottom  
     
    118119      !!---------------------------------------------------------------------- 
    119120      INTEGER                    , INTENT(in   ) ::   kt             ! ocean time step 
     121      INTEGER                    , INTENT(in   ) ::   Kmm            ! time level index 
    120122      REAL(wp), DIMENSION(:,:,:) , INTENT(inout) ::   p_avm          ! momentum Kz (w-points) 
    121123      REAL(wp), DIMENSION(:,:,:) , INTENT(inout) ::   p_avt, p_avs   ! tracer   Kz (w-points) 
     
    148150      !                       !* Critical slope mixing: distribute energy over the time-varying ocean depth, 
    149151      !                                                 using an exponential decay from the seafloor. 
    150       DO jj = 1, jpj                ! part independent of the level 
    151          DO ji = 1, jpi 
    152             zhdep(ji,jj) = gdepw_0(ji,jj,mbkt(ji,jj)+1)       ! depth of the ocean 
    153             zfact(ji,jj) = rau0 * (  1._wp - EXP( -zhdep(ji,jj) / hcri_iwm(ji,jj) )  ) 
    154             IF( zfact(ji,jj) /= 0._wp )   zfact(ji,jj) = ecri_iwm(ji,jj) / zfact(ji,jj) 
    155          END DO 
    156       END DO 
    157 !!gm gde3w ==>>>  check for ssh taken into account.... seem OK gde3w_n=gdept_n - sshn 
    158       DO jk = 2, jpkm1              ! complete with the level-dependent part 
    159          DO jj = 1, jpj              
    160             DO ji = 1, jpi 
    161                IF ( zfact(ji,jj) == 0._wp .OR. wmask(ji,jj,jk) == 0._wp ) THEN   ! optimization 
    162                   zemx_iwm(ji,jj,jk) = 0._wp 
    163                ELSE 
    164                   zemx_iwm(ji,jj,jk) = zfact(ji,jj) * (  EXP( ( gde3w_n(ji,jj,jk  ) - zhdep(ji,jj) ) / hcri_iwm(ji,jj) )     & 
    165                        &                               - EXP( ( gde3w_n(ji,jj,jk-1) - zhdep(ji,jj) ) / hcri_iwm(ji,jj) ) )   & 
    166                        &                            / ( gde3w_n(ji,jj,jk) - gde3w_n(ji,jj,jk-1) ) 
    167                ENDIF 
    168             END DO 
    169          END DO 
    170 !!gm delta(gde3w_n) = e3t_n  !!  Please verify the grid-point position w versus t-point 
     152      DO_2D_11_11 
     153         zhdep(ji,jj) = gdepw_0(ji,jj,mbkt(ji,jj)+1)       ! depth of the ocean 
     154         zfact(ji,jj) = rho0 * (  1._wp - EXP( -zhdep(ji,jj) / hcri_iwm(ji,jj) )  ) 
     155         IF( zfact(ji,jj) /= 0._wp )   zfact(ji,jj) = ecri_iwm(ji,jj) / zfact(ji,jj) 
     156      END_2D 
     157!!gm gde3w ==>>>  check for ssh taken into account.... seem OK gde3w_n=gdept(:,:,:,Kmm) - ssh(:,:,Kmm) 
     158      DO_3D_11_11( 2, jpkm1 ) 
     159         IF ( zfact(ji,jj) == 0._wp .OR. wmask(ji,jj,jk) == 0._wp ) THEN   ! optimization 
     160            zemx_iwm(ji,jj,jk) = 0._wp 
     161         ELSE 
     162            zemx_iwm(ji,jj,jk) = zfact(ji,jj) * (  EXP( ( gde3w(ji,jj,jk  ) - zhdep(ji,jj) ) / hcri_iwm(ji,jj) )     & 
     163                 &                               - EXP( ( gde3w(ji,jj,jk-1) - zhdep(ji,jj) ) / hcri_iwm(ji,jj) ) )   & 
     164                 &                            / ( gde3w(ji,jj,jk) - gde3w(ji,jj,jk-1) ) 
     165         ENDIF 
     166      END_3D 
     167!!gm delta(gde3w) = e3t(:,:,:,Kmm)  !!  Please verify the grid-point position w versus t-point 
    171168!!gm it seems to me that only 1/hcri_iwm  is used ==>  compute it one for all 
    172169 
    173       END DO 
    174170 
    175171      !                        !* Pycnocline-intensified mixing: distribute energy over the time-varying  
     
    182178         zfact(:,:) = 0._wp 
    183179         DO jk = 2, jpkm1              ! part independent of the level 
    184             zfact(:,:) = zfact(:,:) + e3w_n(:,:,jk) * SQRT(  MAX( 0._wp, rn2(:,:,jk) )  ) * wmask(:,:,jk) 
    185          END DO 
    186          ! 
    187          DO jj = 1, jpj 
    188             DO ji = 1, jpi 
    189                IF( zfact(ji,jj) /= 0 )   zfact(ji,jj) = epyc_iwm(ji,jj) / ( rau0 * zfact(ji,jj) ) 
    190             END DO 
    191          END DO 
     180            zfact(:,:) = zfact(:,:) + e3w(:,:,jk,Kmm) * SQRT(  MAX( 0._wp, rn2(:,:,jk) )  ) * wmask(:,:,jk) 
     181         END DO 
     182         ! 
     183         DO_2D_11_11 
     184            IF( zfact(ji,jj) /= 0 )   zfact(ji,jj) = epyc_iwm(ji,jj) / ( rho0 * zfact(ji,jj) ) 
     185         END_2D 
    192186         ! 
    193187         DO jk = 2, jpkm1              ! complete with the level-dependent part 
     
    199193         zfact(:,:) = 0._wp 
    200194         DO jk = 2, jpkm1              ! part independent of the level 
    201             zfact(:,:) = zfact(:,:) + e3w_n(:,:,jk) * MAX( 0._wp, rn2(:,:,jk) ) * wmask(:,:,jk) 
    202          END DO 
    203          ! 
    204          DO jj= 1, jpj 
    205             DO ji = 1, jpi 
    206                IF( zfact(ji,jj) /= 0 )   zfact(ji,jj) = epyc_iwm(ji,jj) / ( rau0 * zfact(ji,jj) ) 
    207             END DO 
    208          END DO 
     195            zfact(:,:) = zfact(:,:) + e3w(:,:,jk,Kmm) * MAX( 0._wp, rn2(:,:,jk) ) * wmask(:,:,jk) 
     196         END DO 
     197         ! 
     198         DO_2D_11_11 
     199            IF( zfact(ji,jj) /= 0 )   zfact(ji,jj) = epyc_iwm(ji,jj) / ( rho0 * zfact(ji,jj) ) 
     200         END_2D 
    209201         ! 
    210202         DO jk = 2, jpkm1              ! complete with the level-dependent part 
     
    220212      zfact(:,:)   = 0._wp 
    221213      DO jk = 2, jpkm1 
    222          zfact(:,:) = zfact(:,:) + e3w_n(:,:,jk) * SQRT(  MAX( 0._wp, rn2(:,:,jk) )  ) * wmask(:,:,jk) 
     214         zfact(:,:) = zfact(:,:) + e3w(:,:,jk,Kmm) * SQRT(  MAX( 0._wp, rn2(:,:,jk) )  ) * wmask(:,:,jk) 
    223215         zwkb(:,:,jk) = zfact(:,:) 
    224216      END DO 
    225217!!gm even better: 
    226218!      DO jk = 2, jpkm1 
    227 !         zwkb(:,:) = zwkb(:,:) + e3w_n(:,:,jk) * SQRT(  MAX( 0._wp, rn2(:,:,jk) )  ) 
     219!         zwkb(:,:) = zwkb(:,:) + e3w(:,:,jk,Kmm) * SQRT(  MAX( 0._wp, rn2(:,:,jk) )  ) 
    228220!      END DO 
    229221!      zfact(:,:) = zwkb(:,:,jpkm1) 
     
    231223!!gm 
    232224      ! 
    233       DO jk = 2, jpkm1 
    234          DO jj = 1, jpj 
    235             DO ji = 1, jpi 
    236                IF( zfact(ji,jj) /= 0 )   zwkb(ji,jj,jk) = zhdep(ji,jj) * ( zfact(ji,jj) - zwkb(ji,jj,jk) )   & 
    237                   &                                     * wmask(ji,jj,jk) / zfact(ji,jj) 
    238             END DO 
    239          END DO 
    240       END DO 
     225      DO_3D_11_11( 2, jpkm1 ) 
     226         IF( zfact(ji,jj) /= 0 )   zwkb(ji,jj,jk) = zhdep(ji,jj) * ( zfact(ji,jj) - zwkb(ji,jj,jk) )   & 
     227            &                                     * wmask(ji,jj,jk) / zfact(ji,jj) 
     228      END_3D 
    241229      zwkb(:,:,1) = zhdep(:,:) * wmask(:,:,1) 
    242230      ! 
    243       DO jk = 2, jpkm1 
    244          DO jj = 1, jpj 
    245             DO ji = 1, jpi 
    246                IF ( rn2(ji,jj,jk) <= 0._wp .OR. wmask(ji,jj,jk) == 0._wp ) THEN   ! optimization 
    247                   zweight(ji,jj,jk) = 0._wp 
    248                ELSE 
    249                   zweight(ji,jj,jk) = rn2(ji,jj,jk) * hbot_iwm(ji,jj)    & 
    250                      &   * (  EXP( -zwkb(ji,jj,jk) / hbot_iwm(ji,jj) ) - EXP( -zwkb(ji,jj,jk-1) / hbot_iwm(ji,jj) )  ) 
    251                ENDIF 
    252             END DO 
    253          END DO 
    254       END DO 
     231      DO_3D_11_11( 2, jpkm1 ) 
     232         IF ( rn2(ji,jj,jk) <= 0._wp .OR. wmask(ji,jj,jk) == 0._wp ) THEN   ! optimization 
     233            zweight(ji,jj,jk) = 0._wp 
     234         ELSE 
     235            zweight(ji,jj,jk) = rn2(ji,jj,jk) * hbot_iwm(ji,jj)    & 
     236               &   * (  EXP( -zwkb(ji,jj,jk) / hbot_iwm(ji,jj) ) - EXP( -zwkb(ji,jj,jk-1) / hbot_iwm(ji,jj) )  ) 
     237         ENDIF 
     238      END_3D 
    255239      ! 
    256240      zfact(:,:) = 0._wp 
     
    259243      END DO 
    260244      ! 
    261       DO jj = 1, jpj 
    262          DO ji = 1, jpi 
    263             IF( zfact(ji,jj) /= 0 )   zfact(ji,jj) = ebot_iwm(ji,jj) / ( rau0 * zfact(ji,jj) ) 
    264          END DO 
    265       END DO 
     245      DO_2D_11_11 
     246         IF( zfact(ji,jj) /= 0 )   zfact(ji,jj) = ebot_iwm(ji,jj) / ( rho0 * zfact(ji,jj) ) 
     247      END_2D 
    266248      ! 
    267249      DO jk = 2, jpkm1              ! complete with the level-dependent part 
    268250         zemx_iwm(:,:,jk) = zemx_iwm(:,:,jk) + zweight(:,:,jk) * zfact(:,:) * wmask(:,:,jk)   & 
    269             &                                / ( gde3w_n(:,:,jk) - gde3w_n(:,:,jk-1) ) 
    270 !!gm  use of e3t_n just above? 
     251            &                                / ( gde3w(:,:,jk) - gde3w(:,:,jk-1) ) 
     252!!gm  use of e3t(:,:,:,Kmm) just above? 
    271253      END DO 
    272254      ! 
    273255!!gm  this is to be replaced by just a constant value znu=1.e-6 m2/s 
    274256      ! Calculate molecular kinematic viscosity 
    275       znu_t(:,:,:) = 1.e-4_wp * (  17.91_wp - 0.53810_wp * tsn(:,:,:,jp_tem) + 0.00694_wp * tsn(:,:,:,jp_tem) * tsn(:,:,:,jp_tem)  & 
    276          &                                  + 0.02305_wp * tsn(:,:,:,jp_sal)  ) * tmask(:,:,:) * r1_rau0 
     257      znu_t(:,:,:) = 1.e-4_wp * (  17.91_wp - 0.53810_wp * ts(:,:,:,jp_tem,Kmm) + 0.00694_wp * ts(:,:,:,jp_tem,Kmm) * ts(:,:,:,jp_tem,Kmm)  & 
     258         &                                  + 0.02305_wp * ts(:,:,:,jp_sal,Kmm)  ) * tmask(:,:,:) * r1_rho0 
    277259      DO jk = 2, jpkm1 
    278260         znu_w(:,:,jk) = 0.5_wp * ( znu_t(:,:,jk-1) + znu_t(:,:,jk) ) * wmask(:,:,jk) 
     
    291273      ! 
    292274      IF( ln_mevar ) THEN              ! Variable mixing efficiency case : modify zav_wave in the 
    293          DO jk = 2, jpkm1              ! energetic (Reb > 480) and buoyancy-controlled (Reb <10.224 ) regimes 
    294             DO jj = 1, jpj 
    295                DO ji = 1, jpi 
    296                   IF( zReb(ji,jj,jk) > 480.00_wp ) THEN 
    297                      zav_wave(ji,jj,jk) = 3.6515_wp * znu_w(ji,jj,jk) * SQRT( zReb(ji,jj,jk) ) 
    298                   ELSEIF( zReb(ji,jj,jk) < 10.224_wp ) THEN 
    299                      zav_wave(ji,jj,jk) = 0.052125_wp * znu_w(ji,jj,jk) * zReb(ji,jj,jk) * SQRT( zReb(ji,jj,jk) ) 
    300                   ENDIF 
    301                END DO 
    302             END DO 
    303          END DO 
     275         DO_3D_11_11( 2, jpkm1 ) 
     276            IF( zReb(ji,jj,jk) > 480.00_wp ) THEN 
     277               zav_wave(ji,jj,jk) = 3.6515_wp * znu_w(ji,jj,jk) * SQRT( zReb(ji,jj,jk) ) 
     278            ELSEIF( zReb(ji,jj,jk) < 10.224_wp ) THEN 
     279               zav_wave(ji,jj,jk) = 0.052125_wp * znu_w(ji,jj,jk) * zReb(ji,jj,jk) * SQRT( zReb(ji,jj,jk) ) 
     280            ENDIF 
     281         END_3D 
    304282      ENDIF 
    305283      ! 
     
    311289         zztmp = 0._wp 
    312290!!gm used of glosum 3D.... 
    313          DO jk = 2, jpkm1 
    314             DO jj = 1, jpj 
    315                DO ji = 1, jpi 
    316                   zztmp = zztmp + e3w_n(ji,jj,jk) * e1e2t(ji,jj)   & 
    317                      &          * MAX( 0._wp, rn2(ji,jj,jk) ) * zav_wave(ji,jj,jk) * wmask(ji,jj,jk) * tmask_i(ji,jj) 
    318                END DO 
    319             END DO 
    320          END DO 
     291         DO_3D_11_11( 2, jpkm1 ) 
     292            zztmp = zztmp + e3w(ji,jj,jk,Kmm) * e1e2t(ji,jj)   & 
     293               &          * MAX( 0._wp, rn2(ji,jj,jk) ) * zav_wave(ji,jj,jk) * wmask(ji,jj,jk) * tmask_i(ji,jj) 
     294         END_3D 
    321295         CALL mpp_sum( 'zdfiwm', zztmp ) 
    322          zztmp = rau0 * zztmp ! Global integral of rauo * Kz * N^2 = power contributing to mixing  
     296         zztmp = rho0 * zztmp ! Global integral of rauo * Kz * N^2 = power contributing to mixing  
    323297         ! 
    324298         IF(lwp) THEN 
     
    337311      IF( ln_tsdiff ) THEN          !* Option for differential mixing of salinity and temperature 
    338312         ztmp1 = 0.505_wp + 0.495_wp * TANH( 0.92_wp * ( LOG10( 1.e-20_wp ) - 0.60_wp ) ) 
    339          DO jk = 2, jpkm1              ! Calculate S/T diffusivity ratio as a function of Reb 
    340             DO jj = 1, jpj 
    341                DO ji = 1, jpi 
    342                   ztmp2 = zReb(ji,jj,jk) * 5._wp * r1_6 
    343                   IF ( ztmp2 > 1.e-20_wp .AND. wmask(ji,jj,jk) == 1._wp ) THEN 
    344                      zav_ratio(ji,jj,jk) = 0.505_wp + 0.495_wp * TANH( 0.92_wp * ( LOG10(ztmp2) - 0.60_wp ) ) 
    345                   ELSE 
    346                      zav_ratio(ji,jj,jk) = ztmp1 * wmask(ji,jj,jk) 
    347                   ENDIF 
    348                END DO 
    349             END DO 
    350          END DO 
     313         DO_3D_11_11( 2, jpkm1 ) 
     314            ztmp2 = zReb(ji,jj,jk) * 5._wp * r1_6 
     315            IF ( ztmp2 > 1.e-20_wp .AND. wmask(ji,jj,jk) == 1._wp ) THEN 
     316               zav_ratio(ji,jj,jk) = 0.505_wp + 0.495_wp * TANH( 0.92_wp * ( LOG10(ztmp2) - 0.60_wp ) ) 
     317            ELSE 
     318               zav_ratio(ji,jj,jk) = ztmp1 * wmask(ji,jj,jk) 
     319            ENDIF 
     320         END_3D 
    351321         CALL iom_put( "av_ratio", zav_ratio ) 
    352322         DO jk = 2, jpkm1           !* update momentum & tracer diffusivity with wave-driven mixing 
     
    368338                                    !* output useful diagnostics: Kz*N^2 ,  
    369339!!gm Kz*N2 should take into account the ratio avs/avt if it is used.... (see diaar5) 
    370                                     !  vertical integral of rau0 * Kz * N^2 , energy density (zemx_iwm) 
     340                                    !  vertical integral of rho0 * Kz * N^2 , energy density (zemx_iwm) 
    371341      IF( iom_use("bflx_iwm") .OR. iom_use("pcmap_iwm") ) THEN 
    372342         ALLOCATE( z2d(jpi,jpj) , z3d(jpi,jpj,jpk) ) 
     
    374344         z2d(:,:) = 0._wp 
    375345         DO jk = 2, jpkm1 
    376             z2d(:,:) = z2d(:,:) + e3w_n(:,:,jk) * z3d(:,:,jk) * wmask(:,:,jk) 
    377          END DO 
    378          z2d(:,:) = rau0 * z2d(:,:) 
     346            z2d(:,:) = z2d(:,:) + e3w(:,:,jk,Kmm) * z3d(:,:,jk) * wmask(:,:,jk) 
     347         END DO 
     348         z2d(:,:) = rho0 * z2d(:,:) 
    379349         CALL iom_put( "bflx_iwm", z3d ) 
    380350         CALL iom_put( "pcmap_iwm", z2d ) 
     
    383353      CALL iom_put( "emix_iwm", zemx_iwm ) 
    384354       
    385       IF(ln_ctl)   CALL prt_ctl(tab3d_1=zav_wave , clinfo1=' iwm - av_wave: ', tab3d_2=avt, clinfo2=' avt: ', kdim=jpk) 
     355      IF(sn_cfctl%l_prtctl)   CALL prt_ctl(tab3d_1=zav_wave , clinfo1=' iwm - av_wave: ', tab3d_2=avt, clinfo2=' avt: ', kdim=jpk) 
    386356      ! 
    387357   END SUBROUTINE zdf_iwm 
     
    414384      !!              de Lavergne et al. in prep., 2017 
    415385      !!---------------------------------------------------------------------- 
    416       INTEGER  ::   ji, jj, jk   ! dummy loop indices 
    417       INTEGER  ::   inum         ! local integer 
     386      INTEGER  ::   ifpr               ! dummy loop indices 
     387      INTEGER  ::   inum               ! local integer 
    418388      INTEGER  ::   ios 
    419389      REAL(wp) ::   zbot, zpyc, zcri   ! local scalars 
    420       !! 
    421       NAMELIST/namzdf_iwm/ nn_zpyc, ln_mevar, ln_tsdiff 
    422       !!---------------------------------------------------------------------- 
    423       ! 
    424       REWIND( numnam_ref )              ! Namelist namzdf_iwm in reference namelist : Wave-driven mixing 
     390      ! 
     391      CHARACTER(len=256)            ::   cn_dir                 ! Root directory for location of ssr files 
     392      INTEGER, PARAMETER            ::   jpiwm  = 5             ! maximum number of files to read 
     393      INTEGER, PARAMETER            ::   jp_mpb = 1 
     394      INTEGER, PARAMETER            ::   jp_mpp = 2 
     395      INTEGER, PARAMETER            ::   jp_mpc = 3 
     396      INTEGER, PARAMETER            ::   jp_dsb = 4 
     397      INTEGER, PARAMETER            ::   jp_dsc = 5 
     398      ! 
     399      TYPE(FLD_N), DIMENSION(jpiwm) ::   slf_iwm                ! array of namelist informations 
     400      TYPE(FLD_N)                   ::   sn_mpb, sn_mpp, sn_mpc ! informations about Mixing Power field to be read 
     401      TYPE(FLD_N)                   ::   sn_dsb, sn_dsc         ! informations about Decay Scale field to be read 
     402      TYPE(FLD  ), DIMENSION(jpiwm) ::   sf_iwm                 ! structure of input fields (file informations, fields read) 
     403      ! 
     404      NAMELIST/namzdf_iwm/ nn_zpyc, ln_mevar, ln_tsdiff, & 
     405         &                 cn_dir, sn_mpb, sn_mpp, sn_mpc, sn_dsb, sn_dsc 
     406      !!---------------------------------------------------------------------- 
     407      ! 
    425408      READ  ( numnam_ref, namzdf_iwm, IOSTAT = ios, ERR = 901) 
    426409901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namzdf_iwm in reference namelist' ) 
    427410      ! 
    428       REWIND( numnam_cfg )              ! Namelist namzdf_iwm in configuration namelist : Wave-driven mixing 
    429411      READ  ( numnam_cfg, namzdf_iwm, IOSTAT = ios, ERR = 902 ) 
    430412902   IF( ios >  0 )   CALL ctl_nam ( ios , 'namzdf_iwm in configuration namelist' ) 
     
    456438      IF( zdf_iwm_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'zdf_iwm_init : unable to allocate iwm arrays' ) 
    457439      ! 
     440      ! store namelist information in an array 
     441      slf_iwm(jp_mpb) = sn_mpb ; slf_iwm(jp_mpp) = sn_mpp ; slf_iwm(jp_mpc) = sn_mpc 
     442      slf_iwm(jp_dsb) = sn_dsb ; slf_iwm(jp_dsc) = sn_dsc 
     443      ! 
     444      DO ifpr= 1, jpiwm 
     445         ALLOCATE( sf_iwm(ifpr)%fnow(jpi,jpj,1)   ) 
     446         IF( slf_iwm(ifpr)%ln_tint )ALLOCATE( sf_iwm(ifpr)%fdta(jpi,jpj,1,2) ) 
     447      END DO 
     448 
     449      ! fill sf_iwm with sf_iwm and control print 
     450      CALL fld_fill( sf_iwm, slf_iwm , cn_dir, 'zdfiwm_init', 'iwm input file', 'namiwm' ) 
     451 
     452      !                             ! hard-coded default definition (to be defined in namelist ?) 
     453      sf_iwm(jp_mpb)%fnow(:,:,1) = 1.e-6 
     454      sf_iwm(jp_mpp)%fnow(:,:,1) = 1.e-6 
     455      sf_iwm(jp_mpc)%fnow(:,:,1) = 1.e-10 
     456      sf_iwm(jp_dsb)%fnow(:,:,1) = 100. 
     457      sf_iwm(jp_dsc)%fnow(:,:,1) = 100. 
     458 
    458459      !                             ! read necessary fields 
    459       CALL iom_open('mixing_power_bot',inum)       ! energy flux for high-mode wave breaking [W/m2] 
    460       CALL iom_get  (inum, jpdom_data, 'field', ebot_iwm, 1 )  
    461       CALL iom_close(inum) 
    462       ! 
    463       CALL iom_open('mixing_power_pyc',inum)       ! energy flux for pynocline-intensified wave breaking [W/m2] 
    464       CALL iom_get  (inum, jpdom_data, 'field', epyc_iwm, 1 ) 
    465       CALL iom_close(inum) 
    466       ! 
    467       CALL iom_open('mixing_power_cri',inum)       ! energy flux for critical slope wave breaking [W/m2] 
    468       CALL iom_get  (inum, jpdom_data, 'field', ecri_iwm, 1 ) 
    469       CALL iom_close(inum) 
    470       ! 
    471       CALL iom_open('decay_scale_bot',inum)        ! spatially variable decay scale for high-mode wave breaking [m] 
    472       CALL iom_get  (inum, jpdom_data, 'field', hbot_iwm, 1 ) 
    473       CALL iom_close(inum) 
    474       ! 
    475       CALL iom_open('decay_scale_cri',inum)        ! spatially variable decay scale for critical slope wave breaking [m] 
    476       CALL iom_get  (inum, jpdom_data, 'field', hcri_iwm, 1 ) 
    477       CALL iom_close(inum) 
    478  
    479       ebot_iwm(:,:) = ebot_iwm(:,:) * ssmask(:,:) 
    480       epyc_iwm(:,:) = epyc_iwm(:,:) * ssmask(:,:) 
    481       ecri_iwm(:,:) = ecri_iwm(:,:) * ssmask(:,:) 
     460      CALL fld_read( nit000, 1, sf_iwm ) 
     461 
     462      ebot_iwm(:,:) = sf_iwm(1)%fnow(:,:,1) * ssmask(:,:) ! energy flux for high-mode wave breaking [W/m2] 
     463      epyc_iwm(:,:) = sf_iwm(2)%fnow(:,:,1) * ssmask(:,:) ! energy flux for pynocline-intensified wave breaking [W/m2] 
     464      ecri_iwm(:,:) = sf_iwm(3)%fnow(:,:,1) * ssmask(:,:) ! energy flux for critical slope wave breaking [W/m2] 
     465      hbot_iwm(:,:) = sf_iwm(4)%fnow(:,:,1)               ! spatially variable decay scale for high-mode wave breaking [m] 
     466      hcri_iwm(:,:) = sf_iwm(5)%fnow(:,:,1)               ! spatially variable decay scale for critical slope wave breaking [m] 
    482467 
    483468      zbot = glob_sum( 'zdfiwm', e1e2t(:,:) * ebot_iwm(:,:) ) 
    484469      zpyc = glob_sum( 'zdfiwm', e1e2t(:,:) * epyc_iwm(:,:) ) 
    485470      zcri = glob_sum( 'zdfiwm', e1e2t(:,:) * ecri_iwm(:,:) ) 
     471 
    486472      IF(lwp) THEN 
    487473         WRITE(numout,*) '      High-mode wave-breaking energy:             ', zbot * 1.e-12_wp, 'TW' 
Note: See TracChangeset for help on using the changeset viewer.