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 3995 for branches/2013/dev_r3948_NOC_FK/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_mle.F90 – NEMO

Ignore:
Timestamp:
2013-07-31T14:27:58+02:00 (11 years ago)
Author:
gm
Message:

dev_r3406_CNRS_LIM3: fix a bug on zpsim_v and remove unnecessary print

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2013/dev_r3948_NOC_FK/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_mle.F90

    r3959 r3995  
    106106      CALL wrk_alloc( jpi, jpj, jpk, zpsi_uw, zpsi_vw) 
    107107      CALL wrk_alloc( jpi, jpj, inml_mle) 
    108 !!gm   Caution HERE : * compute the nmln from the 10m density to deal with the diurnal cycle 
    109 !!gm                  * check that nmld is the number of T-level not w-levels... 
    110  
    111       DO jk = jpkm1, nlb10, -1         ! from the bottom to nlb10 
     108      ! 
     109      !                                      !==  MLD used for MLE  ==! 
     110      !                                                ! compute from the 10m density to deal with the diurnal cycle 
     111      inml_mle(:,:) = mbkt(:,:) + 1                    ! init. to number of ocean w-level (T-level + 1) 
     112      DO jk = jpkm1, nlb10, -1                         ! from the bottom to nlb10 (10m) 
     113         DO jj = 1, jpj 
     114            DO ji = 1, jpi                             ! index of the w-level at the ML based 
     115               IF( rhop(ji,jj,jk) > rhop(ji,jj,nla10) + rn_rho_c_mle )   inml_mle(ji,jj) = jk      ! Mixed layer 
     116            END DO 
     117         END DO 
     118      END DO 
     119      ikmax = MAXVAL( inml_mle(:,:) )                  ! max level of the computation 
     120      ! 
     121      ! 
     122      zmld(:,:) = 0._wp                      !==   Horizontal shape of the MLE  ==! 
     123      zbm (:,:) = 0._wp 
     124      zn2 (:,:) = 0._wp 
     125      DO jk = 1, ikmax                                 ! MLD and mean buoyancy and N2 over the mixed layer 
    112126         DO jj = 1, jpj 
    113127            DO ji = 1, jpi 
    114                IF( rhop(ji,jj,jk) > rhop(ji,jj,nla10) + rn_rho_c_mle )   inml_mle(ji,jj) = jk      ! Mixed layer 
     128               zc = fse3t(ji,jj,jk) * REAL( MIN( MAX( 0, inml_mle(ji,jj)-jk ) , 1  )  )    ! zc being 0 outside the ML t-points 
     129               zmld(ji,jj) = zmld(ji,jj) + zc 
     130               zbm (ji,jj) = zbm (ji,jj) + zc * (rau0 - rhop(ji,jj,jk) ) * r1_rau0 
     131               zn2 (ji,jj) = zn2 (ji,jj) + zc * (rn2(ji,jj,jk)+rn2(ji,jj,jk+1))*0.5_wp 
    115132            END DO 
    116133         END DO 
    117134      END DO 
    118       ikmax = MAXVAL( inml_mle(:,:) )                  ! max level of the computation 
    119       ! 
    120       zbm(:,:) = 0._wp   ;   zmld(:,:) = 0._wp; zn2(:,:) = 0._wp   ! Horizontal shape of the MLE 
    121       ! 
    122  
    123  
    124    IF(lwp) WRITE(numout,*) 'mle 0' 
    125  
    126            ilocs = MAXLOC( inml_mle(:,:) - 1 ) 
    127            ii = ilocs(1) 
    128            ij = ilocs(2) 
    129             ! ii=120 
    130             ! ij=110 
    131              !  inml_mle(ii,ij) = Max( 4, inml_mle(ii,ij) ) 
    132             WRITE(numout,*) 'mle zmld max', ii,ij, inml_mle(ii,ij), 'ikmax', ikmax 
    133  
    134  
    135       DO jk = 1, ikmax                             ! MLD and mean buoyancy and N2 over the mixed layer 
    136    IF(lwp) WRITE(numout,*) 'mle zc',jk,' c ', REAL( MIN( MAX( 0, inml_mle(ii,ij)-jk ) , 1  )  ) 
    137          DO jj = 1, jpj 
    138             DO ji = 1, jpi 
    139                zc = fse3t(ji,jj,jk) * REAL( MIN( MAX( 0, inml_mle(ji,jj)-jk ) , 1  )  )    ! e3t being 0 outside the ML 
    140                zmld(ji,jj) = zmld(ji,jj) + zc 
    141                zbm (ji,jj) = zbm (ji,jj) + zc * (rau0 - rhop(ji,jj,jk) ) / rau0 
    142                zn2 (ji,jj) = zn2 (ji,jj) + zc * (rn2(ji,jj,jk)+rn2(ji,jj,jk+1))*0.5_wp 
    143             END DO 
    144          END DO 
    145       END DO 
    146  
    147       SELECT CASE( nn_mld_uv )                           ! MLD at u- & v-pts 
     135 
     136      SELECT CASE( nn_mld_uv )                         ! MLD at u- & v-pts 
    148137      CASE ( 0 )                                               != min of the 2 neighbour MLDs 
    149138         DO jj = 1, jpjm1 
     
    168157         END DO 
    169158      END SELECT 
    170  
     159      !                                                ! convert density into buoyancy 
    171160      zbm(:,:) = + grav * zbm(:,:) / MAX( fse3t(:,:,1), zmld(:,:) ) 
    172  
    173  
    174    write(numout,*) 'n2: max ', MAXVAL(zn2(:,:) , MASK = tmask(:,:,1) ==1. ),    & 
    175          &            ' min ', MINVAL(zn2(:,:) , MASK = tmask(:,:,1) ==1. ) 
    176    write(numout,*) 'h: max ', MAXVAL(zhu(1:jpim1,1:jpjm1) ), ' min ', MINVAL(zhu(1:jpim1,1:jpjm1) ) 
    177    write(numout,*) 'h: max ', MAXVAL(zhv(1:jpim1,1:jpjm1) ), ' min ', MINVAL(zhv(1:jpim1,1:jpjm1) ) 
    178  
    179       IF(lwp) write(numout,*) 'mle: max zbm', MAXVAL(-zbm(:,:) , MASK = tmask(:,:,1) ==1. ) / grav,    & 
    180          &                    ' min '       , MINVAL(-zbm(:,:) , MASK = tmask(:,:,1) ==1. ) / grav 
    181  
    182  
    183       ! 
    184       ! Magnitude of the MLE stream function: 
     161      ! 
     162      ! 
     163      !                                      !==  Magnitude of the MLE stream function  ==! 
    185164      ! 
    186165      !                 di[bm]  Ds 
     
    188167      !                  e1u   Lf fu                      and the e2u for the "transport" 
    189168      !                                                      (not *e3u as divided by e3u at the end) 
    190 ! MLD criteria set in zdfmxl: rho_c = 0.01 kg/m3 
    191       zpsim_u(:,:) = 0._wp 
    192       zpsim_v(:,:) = 0._wp 
    193       ! 
    194  
     169      ! 
    195170      IF( nn_mle == 0 ) THEN           ! Fox-Kemper et al. 2010 formulation 
    196171         DO jj = 1, jpjm1 
     
    202177               zpsim_v(ji,jj) = rn_ce * zhv(ji,jj) * zhv(ji,jj)  * e1v(ji,jj)                                            & 
    203178                  &           * ( zbm(ji,jj+1) - zbm(ji,jj) ) * MIN( 111.e3_wp          , e2v(ji,jj)                 )   & 
    204                   &           / (         e1u(ji,jj)          * MAX( rn_lf * rfv(ji,jj) , SQRT( rb_c * zhv(ji,jj) ) )   ) 
     179                  &           / (         e2v(ji,jj)          * MAX( rn_lf * rfv(ji,jj) , SQRT( rb_c * zhv(ji,jj) ) )   ) 
    205180            END DO 
    206181         END DO 
     
    212187                  &                  * ( zbm(ji+1,jj) - zbm(ji,jj) ) * MIN( 111.e3_wp , e1u(ji,jj) ) 
    213188                  ! 
    214                zpsim_v(ji,jj) = rc_f *   zhv(ji,jj)   * zhv(ji,jj)   * e1v(ji,jj) / e1u(ji,jj)          & 
     189               zpsim_v(ji,jj) = rc_f *   zhv(ji,jj)   * zhv(ji,jj)   * e1v(ji,jj) / e2v(ji,jj)          & 
    215190                  &                  * ( zbm(ji,jj+1) - zbm(ji,jj) ) * MIN( 111.e3_wp , e2v(ji,jj) ) 
    216191            END DO 
     
    227202      ENDIF 
    228203      ! 
    229       ! 
    230       !                                        ! structure function value at uw- and vw-points 
     204      !                                      !==  structure function value at uw- and vw-points  ==! 
    231205      zhu(:,:) = 1._wp / zhu(:,:)                   ! hu --> 1/hu 
    232206      zhv(:,:) = 1._wp / zhv(:,:) 
     
    249223         END DO 
    250224      END DO 
    251  
    252 !!!gm start 
    253    write(numout,*) 'psi max u ', MAXVAL(ABS(zpsim_u(1:jpim1,1:jpjm1)*umask(1:jpim1,1:jpjm1,1) )),    & 
    254          &                ' v ', MAXVAL(ABS(zpsim_v(1:jpim1,1:jpjm1)*vmask(1:jpim1,1:jpjm1,1) )) 
    255  
    256             ilocs = MAXLOC( ABS(zpsim_v(1:jpim1,1:jpjm1)*vmask(1:jpim1,1:jpjm1,1) )) !, inml_mle(1:jpim1,1:jpjm1 ) 
    257             ii = ilocs(1) 
    258             ij = ilocs(2) 
    259                WRITE(numout,*)        'hu, hv'      , 1./zhu(ii,ij), 1./zhv(ii,ij) 
    260                WRITE(numout,*) 
    261                WRITE(numout,*)        'psim at i,j = ', ii, ij, zpsim_u(ii,ij) 
    262 !           ji=ii 
    263 !           jj=ij 
    264  
    265       DO jk = 2, ikmax                                ! start from 2 : surface value = 0 
    266                zcuw = 1.e0 - ( fsdepw(ji+1,jj,jk) + fsdepw(ji,jj,jk) ) * zhu(ji,jj) 
    267                zcvw = 1.e0 - ( fsdepw(ji,jj+1,jk) + fsdepw(ji,jj,jk) ) * zhv(ji,jj) 
    268                zcuw = zcuw * zcuw 
    269                zcvw = zcvw * zcvw 
    270                zmuw = MAX(  0._wp , ( 1._wp - zcuw ) * ( 1._wp + r5_21 * zcuw )  ) 
    271                zmvw = MAX(  0._wp , ( 1._wp - zcvw ) * ( 1._wp + r5_21 * zcvw )  ) 
    272                ! 
    273                WRITE(numout,*)        jk, zmuw, zmvw 
    274       END DO 
    275  
    276    IF(lwp) WRITE(numout,*) 'mle d' 
    277       !!!gm 
    278       IF(lwp) write(numout,*) 'max uw', MAXVAL(zpsi_uw,MASK = umask==1.), 'vw', MAXVAL(zpsi_vw,MASK = vmask==1.) 
    279       IF(lwp) write(numout,*) 'min uw', MinVAL(zpsi_uw,MASK = umask==1.), 'vw', MinVAL(zpsi_vw,MASK = vmask==1.) 
    280       !!!gm 
    281  
    282       IF(lwp) write(numout,*) 'rho i i1', rhop(ii,ij,1), rhop(ii+1,ij,1) 
    283       IF(lwp) write(numout,*) 'zbm i i1', zbm(ii,ij), zbm(ii+1,ij) 
    284       IF(lwp) write(numout,*) 'Dbm i i1', zbm(ii+1,ij)-zbm(ii,ij) 
    285       IF(lwp) write(numout,*) 'zmld i i1', zmld(ii,ij), zmld(ii+1,ij) 
    286       IF(lwp) write(numout,*) 'Dmld i i1', zmld(ii+1,ij)-zmld(ii,ij) 
    287       IF(lwp) write(numout,*) 'hu   i   ', zhu(ii,ij) 
    288       IF(lwp) write(numout,*) 'mld  i i1', inml_mle(ii,ij), inml_mle(ii+1,ij) 
    289  
    290       IF(lwp) write(numout,*) 'rho j j1', rhop(ii,ij,1), rhop(ii,ij+1,1) 
    291       IF(lwp) write(numout,*) 'zbm j j1', zbm(ii,ij), zbm(ii,ij+1) 
    292       IF(lwp) write(numout,*) 'Dbm j j1', zbm(ii,ij+1)-zbm(ii,ij) 
    293       IF(lwp) write(numout,*) 'zmld j j1', zmld(ii,ij), zmld(ii,ij+1) 
    294       IF(lwp) write(numout,*) 'Dmld j j1', zmld(ii,ij+1)-zmld(ii,ij) 
    295       IF(lwp) write(numout,*) 'hv   i   ', zhv(ii,ij) 
    296       IF(lwp) write(numout,*) 'mld  i i1', inml_mle(ii,ij), inml_mle(ii,ij+1) 
    297  
    298       WRITE(numout,*) 'mle zpsi: i,j,k',ii,ij, 'mask', tmask(ii,ij,1), 'mld', zmld(ii,ij) 
    299       DO jk = 1, jpk 
    300          WRITE(numout,*) jk, zpsi_uw(ii,ij,jk), ' u - v', zpsi_vw(ii,ij,jk) 
    301       END DO 
    302 !!!gm end 
    303  
    304  
    305  
    306  
    307  
     225      ! 
     226      !                                      !==  transport increased by the MLE induced transport ==! 
    308227      DO jk = 1, ikmax 
    309228         DO jj = 1, jpjm1                          ! CAUTION pu,pv must be defined at row/column i=1 / j=1 
     
    321240      END DO 
    322241 
    323 !!!gm start 
    324       IF(lwp) THEN 
    325          write(numout,*) jk, 'velocity min/max' 
    326          Do jk=1,jpk 
    327          write(numout,*) jk, 'uw', MAXVAL(ABS((zpsi_uw(:,:,jk)-zpsi_uw(:,:,jk+1))/( e2u(:,:)*fse3u(:,:,jk)))),  & 
    328             &                'vw', MAXVAL(ABS((zpsi_vw(:,:,jk)-zpsi_vw(:,:,jk+1))/( e1v(:,:)*fse3v(:,:,jk)))) 
    329          end do 
    330       endif 
    331 !!!gm end 
    332  
    333       IF( cdtype == 'TRA') THEN 
     242      IF( cdtype == 'TRA') THEN              !==  outputs  ==! 
    334243         ! 
    335244         zLf_NH(:,:) = SQRT( rb_c * zmld(:,:) ) * r1_ft(:,:)      ! Lf = N H / f 
    336          ! 
     245         CALL iom_put( "Lf_NHpf" , zLf_NH  )    ! Lf = N H / f 
     246         ! 
     247         zpsim_u(jpi,:) = 0._wp   ;   zpsim_u(:,jpj) = 0._wp      ! set a value where it is undefined 
     248         zpsim_v(jpi,:) = 0._wp   ;   zpsim_v(:,jpj) = 0._wp      ! before the output 
    337249         CALL iom_put( "psiu_mle", zpsi_uw )    ! i-mle streamfunction 
    338250         CALL iom_put( "psiv_mle", zpsi_vw )    ! j-mle streamfunction 
    339          CALL iom_put( "Lf_NHpf" , zLf_NH  )    ! Lf = N H / f 
    340251      ENDIF 
    341252      CALL wrk_dealloc( jpi, jpj, zpsim_u, zpsim_v, zmld, zbm, zhu, zhv, zn2, zLf_NH, zLf_MH) 
Note: See TracChangeset for help on using the changeset viewer.