Changeset 5332


Ignore:
Timestamp:
2015-06-01T18:59:09+02:00 (5 years ago)
Author:
mathiot
Message:

ISOMIP config code removed (mv to unsupported configuration); correction of bug about initialisation of top friction (ISF only), about coastline modification if ISF activated and top and bottom e3uw if ISF.

Location:
trunk/NEMOGCM/NEMO/OPA_SRC
Files:
3 edited

Legend:

Unmodified
Added
Removed
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/domzgr.F90

    r5120 r5332  
    472472         risfdep(:,:)=0.e0 
    473473         misfdep(:,:)=1 
    474          ! 
    475          ! (ISF) TODO build ice draft netcdf file for isomip and build the corresponding part of code 
    476          IF( cp_cfg == "isomip" .AND. ln_isfcav ) THEN  
    477             risfdep(:,:)=200.e0  
    478             misfdep(:,:)=1  
    479             ij0 = 1 ; ij1 = 40  
    480             DO jj = mj0(ij0), mj1(ij1)  
    481                risfdep(:,jj)=700.0_wp-(gphit(:,jj)+80.0_wp)*125.0_wp  
    482             END DO  
    483             WHERE( bathy(:,:) <= 0._wp )  risfdep(:,:) = 0._wp  
    484          !  
    485          ELSEIF ( cp_cfg == "isomip2" .AND. ln_isfcav ) THEN 
    486          !  
    487             risfdep(:,:)=0.e0 
    488             misfdep(:,:)=1 
    489             ij0 = 1 ; ij1 = 40 
    490             DO jj = mj0(ij0), mj1(ij1) 
    491                risfdep(:,jj)=700.0_wp-(gphit(:,jj)+80.0_wp)*125.0_wp 
    492             END DO 
    493             WHERE( bathy(:,:) <= 0._wp )  risfdep(:,:) = 0._wp 
    494          END IF 
    495474         ! 
    496475         DEALLOCATE( idta, zdta ) 
     
    969948      !! 
    970949      INTEGER  ::   ji, jj, jk       ! dummy loop indices 
    971       INTEGER  ::   ik, it          ! temporary integers 
     950      INTEGER  ::   ik, it, ikb, ikt ! temporary integers 
    972951      LOGICAL  ::   ll_print         ! Allow  control print for debugging 
    973952      REAL(wp) ::   ze3tp , ze3wp    ! Last ocean level thickness at T- and W-points 
     
    11521131      IF ( ln_isfcav ) THEN 
    11531132      ! (ISF) define e3uw (adapted for 2 cells in the water column) 
    1154       ! Need to test if the modification of only mikt and mbkt levels is enough 
    1155          DO jk = 2,jpk                           
    1156             DO jj = 1, jpjm1  
    1157                DO ji = 1, fs_jpim1   ! vector opt.  
    1158                   e3uw_0(ji,jj,jk) = MIN( gdept_0(ji,jj,jk), gdept_0(ji+1,jj  ,jk) ) & 
    1159                     &   - MAX( gdept_0(ji,jj,jk-1), gdept_0(ji+1,jj  ,jk-1) ) 
    1160                   e3vw_0(ji,jj,jk) = MIN( gdept_0(ji,jj,jk), gdept_0(ji  ,jj+1,jk) ) & 
    1161                     &   - MAX( gdept_0(ji,jj,jk-1), gdept_0(ji  ,jj+1,jk-1) ) 
    1162                END DO  
    1163             END DO  
     1133         DO jj = 2, jpjm1  
     1134            DO ji = 2, fs_jpim1   ! vector opt.  
     1135               ikb = MAX(mbathy (ji,jj),mbathy (ji+1,jj)) 
     1136               ikt = MAX(misfdep(ji,jj),misfdep(ji+1,jj)) 
     1137               IF (ikb == ikt+1) e3uw_0(ji,jj,ikb) =  MIN( gdept_0(ji,jj,ikb  ), gdept_0(ji+1,jj  ,ikb  ) ) & 
     1138                                       &            - MAX( gdept_0(ji,jj,ikb-1), gdept_0(ji+1,jj  ,ikb-1) ) 
     1139               ikb = MAX(mbathy (ji,jj),mbathy (ji,jj+1)) 
     1140               ikt = MAX(misfdep(ji,jj),misfdep(ji,jj+1)) 
     1141               IF (ikb == ikt+1) e3vw_0(ji,jj,ikb) =  MIN( gdept_0(ji,jj,ikb  ), gdept_0(ji  ,jj+1,ikb  ) ) & 
     1142                                       &            - MAX( gdept_0(ji,jj,ikb-1), gdept_0(ji  ,jj+1,ikb-1) ) 
     1143            END DO 
    11641144         END DO 
    11651145      END IF 
    1166        
     1146 
    11671147      CALL lbc_lnk( e3u_0 , 'U', 1._wp )   ;   CALL lbc_lnk( e3uw_0, 'U', 1._wp )   ! lateral boundary conditions 
    11681148      CALL lbc_lnk( e3v_0 , 'V', 1._wp )   ;   CALL lbc_lnk( e3vw_0, 'V', 1._wp ) 
     
    15381518  
    15391519 ! remove single point "bay" on isf coast line in the ice shelf draft' 
    1540          DO jk = 1, jpk 
     1520         DO jk = 2, jpk 
    15411521            WHERE (misfdep==0) misfdep=jpk 
    15421522            zmask=0 
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/istate.F90

    r5163 r5332  
    110110         ELSEIF( cp_cfg == 'gyre' ) THEN          
    111111            CALL istate_gyre                     ! GYRE  configuration : start from pre-defined T-S fields 
    112         ELSEIF( cp_cfg == 'isomip' .OR. cp_cfg == 'isomip2') THEN 
    113             IF(lwp) WRITE(numout,*) 'Initialization of T+S for ISOMIP domain'  
    114             tsn(:,:,:,jp_tem)=-1.9*tmask(:,:,:)          ! ISOMIP configuration : start from constant T+S fields  
    115             tsn(:,:,:,jp_sal)=34.4*tmask(:,:,:) 
    116             tsb(:,:,:,:)=tsn(:,:,:,:)   
    117112         ELSE                                    ! Initial T-S, U-V fields read in files 
    118113            IF ( ln_tsd_init ) THEN              ! read 3D T and S data at nit000 
  • trunk/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfbfr.F90

    r5120 r5332  
    171171            END DO 
    172172         END DO 
     173         CALL lbc_lnk( bfrua, 'U', 1. )   ;   CALL lbc_lnk( bfrva, 'V', 1. )      ! Lateral boundary condition 
     174 
    173175         IF ( ln_isfcav ) THEN 
    174176            DO jj = 2, jpjm1 
    175177               DO ji = 2, jpim1 
    176178                  ! (ISF) ======================================================================== 
    177                   ikbu = miku(ji,jj)         ! ocean bottom level at u- and v-points  
    178                   ikbv = mikv(ji,jj)         ! (deepest ocean u- and v-points) 
     179                  ikbu = miku(ji,jj)         ! ocean top level at u- and v-points  
     180                  ikbv = mikv(ji,jj)         ! (1st wet ocean u- and v-points) 
    179181                  ! 
    180182                  zvu  = 0.25 * (  vn(ji,jj  ,ikbu) + vn(ji+1,jj  ,ikbu)     & 
     
    183185                     &           + un(ji,jj+1,ikbv) + un(ji-1,jj+1,ikbv)  ) 
    184186              ! 
    185                   zecu = SQRT(  un(ji,jj,ikbu) * un(ji,jj,ikbu) + zvu*zvu + rn_bfeb2 ) 
    186                   zecv = SQRT(  vn(ji,jj,ikbv) * vn(ji,jj,ikbv) + zuv*zuv + rn_bfeb2 ) 
     187                  zecu = SQRT(  un(ji,jj,ikbu) * un(ji,jj,ikbu) + zvu*zvu + rn_tfeb2 ) 
     188                  zecv = SQRT(  vn(ji,jj,ikbv) * vn(ji,jj,ikbv) + zuv*zuv + rn_tfeb2 ) 
    187189              ! 
    188190                  tfrua(ji,jj) = - 0.5_wp * ( ztfrt(ji,jj) + ztfrt(ji+1,jj  ) ) * zecu * (1._wp - umask(ji,jj,1)) 
     
    202204               END DO 
    203205            END DO 
    204          END IF 
    205          ! 
    206          CALL lbc_lnk( bfrua, 'U', 1. )   ;   CALL lbc_lnk( bfrva, 'V', 1. )      ! Lateral boundary condition 
     206            CALL lbc_lnk( tfrua, 'U', 1. )   ;   CALL lbc_lnk( tfrva, 'V', 1. )      ! Lateral boundary condition 
     207         END IF 
     208         ! 
    207209         ! 
    208210         IF(ln_ctl)   CALL prt_ctl( tab2d_1=bfrua, clinfo1=' bfr  - u: ', mask1=umask,        & 
     
    277279            IF(lwp) WRITE(numout,*) '      coef rn_bfri2 enhancement factor                rn_bfrien  = ',rn_bfrien 
    278280         ENDIF 
    279          IF(lwp) WRITE(numout,*) '      top    friction coef.   rn_bfri1  = ', rn_bfri1 
    280          IF( ln_tfr2d ) THEN 
    281             IF(lwp) WRITE(numout,*) '      read coef. enhancement distribution from file   ln_tfr2d  = ', ln_tfr2d 
    282             IF(lwp) WRITE(numout,*) '      coef rn_tfri2 enhancement factor                rn_tfrien  = ',rn_tfrien 
    283          ENDIF 
     281         IF ( ln_isfcav ) THEN 
     282            IF(lwp) WRITE(numout,*) '      top    friction coef.   rn_bfri1  = ', rn_tfri1 
     283            IF( ln_tfr2d ) THEN 
     284               IF(lwp) WRITE(numout,*) '      read coef. enhancement distribution from file   ln_tfr2d  = ', ln_tfr2d 
     285               IF(lwp) WRITE(numout,*) '      coef rn_tfri2 enhancement factor                rn_tfrien  = ',rn_tfrien 
     286            ENDIF 
     287         END IF 
    284288         ! 
    285289         IF(ln_bfr2d) THEN 
     
    295299         bfrua(:,:) = - bfrcoef2d(:,:) 
    296300         bfrva(:,:) = - bfrcoef2d(:,:) 
     301         ! 
     302         IF ( ln_isfcav ) THEN 
     303            IF(ln_tfr2d) THEN 
     304               ! tfr_coef is a coefficient in [0,1] giving the mask where to apply the bfr enhancement 
     305               CALL iom_open('tfr_coef.nc',inum) 
     306               CALL iom_get (inum, jpdom_data, 'tfr_coef',tfrcoef2d,1) ! tfrcoef2d is used as tmp array 
     307               CALL iom_close(inum) 
     308               tfrcoef2d(:,:) = rn_tfri1 * ( 1 + rn_tfrien * tfrcoef2d(:,:) ) 
     309            ELSE 
     310               tfrcoef2d(:,:) = rn_tfri1  ! initialize tfrcoef2d to the namelist variable 
     311            ENDIF 
     312            ! 
     313            tfrua(:,:) = - tfrcoef2d(:,:) 
     314            tfrva(:,:) = - tfrcoef2d(:,:) 
     315         END IF 
    297316         ! 
    298317      CASE( 2 ) 
     
    311330            IF(lwp) WRITE(numout,*) '      coef rn_bfri2 enhancement factor                rn_bfrien  = ',rn_bfrien 
    312331         ENDIF 
    313          IF(lwp) WRITE(numout,*) '      quadratic top    friction' 
    314          IF(lwp) WRITE(numout,*) '      friction coef.   rn_bfri2  = ', rn_tfri2 
    315          IF(lwp) WRITE(numout,*) '      Max. coef. (log case)   rn_tfri2_max  = ', rn_tfri2_max 
    316          IF(lwp) WRITE(numout,*) '      background tke   rn_tfeb2  = ', rn_tfeb2 
    317          IF(lwp) WRITE(numout,*) '      log formulation   ln_tfr2d = ', ln_loglayer 
    318          IF(lwp) WRITE(numout,*) '      bottom roughness  rn_tfrz0 [m] = ', rn_tfrz0 
    319          IF( rn_tfrz0<=0.e0 ) THEN 
    320             WRITE(ctmp1,*) '      bottom roughness must be strictly positive' 
    321             CALL ctl_stop( ctmp1 ) 
    322          ENDIF 
    323          IF( ln_tfr2d ) THEN 
    324             IF(lwp) WRITE(numout,*) '      read coef. enhancement distribution from file   ln_tfr2d  = ', ln_tfr2d 
    325             IF(lwp) WRITE(numout,*) '      coef rn_tfri2 enhancement factor                rn_tfrien  = ',rn_tfrien 
    326          ENDIF 
     332         IF ( ln_isfcav ) THEN 
     333            IF(lwp) WRITE(numout,*) '      quadratic top    friction' 
     334            IF(lwp) WRITE(numout,*) '      friction coef.    rn_tfri2     = ', rn_tfri2 
     335            IF(lwp) WRITE(numout,*) '      Max. coef. (log case)   rn_tfri2_max  = ', rn_tfri2_max 
     336            IF(lwp) WRITE(numout,*) '      background tke    rn_tfeb2     = ', rn_tfeb2 
     337            IF(lwp) WRITE(numout,*) '      log formulation   ln_tfr2d     = ', ln_loglayer 
     338            IF(lwp) WRITE(numout,*) '      top roughness     rn_tfrz0 [m] = ', rn_tfrz0 
     339            IF( rn_tfrz0<=0.e0 ) THEN 
     340               WRITE(ctmp1,*) '      top roughness must be strictly positive' 
     341               CALL ctl_stop( ctmp1 ) 
     342            ENDIF 
     343            IF( ln_tfr2d ) THEN 
     344               IF(lwp) WRITE(numout,*) '      read coef. enhancement distribution from file   ln_tfr2d  = ', ln_tfr2d 
     345               IF(lwp) WRITE(numout,*) '      coef rn_tfri2 enhancement factor                rn_tfrien  = ',rn_tfrien 
     346            ENDIF 
     347         END IF 
    327348         ! 
    328349         IF(ln_bfr2d) THEN 
     
    336357            bfrcoef2d(:,:) = rn_bfri2  ! initialize bfrcoef2d to the namelist variable 
    337358         ENDIF 
     359          
     360         IF ( ln_isfcav ) THEN 
     361            IF(ln_tfr2d) THEN 
     362               ! tfr_coef is a coefficient in [0,1] giving the mask where to apply the bfr enhancement 
     363               CALL iom_open('tfr_coef.nc',inum) 
     364               CALL iom_get (inum, jpdom_data, 'tfr_coef',tfrcoef2d,1) ! tfrcoef2d is used as tmp array 
     365               CALL iom_close(inum) 
     366               ! 
     367               tfrcoef2d(:,:) = rn_tfri2 * ( 1 + rn_tfrien * tfrcoef2d(:,:) ) 
     368            ELSE 
     369               tfrcoef2d(:,:) = rn_tfri2  ! initialize tfrcoef2d to the namelist variable 
     370            ENDIF 
     371         END IF 
    338372         ! 
    339373         IF ( ln_loglayer.AND.(.NOT.lk_vvl) ) THEN ! set "log layer" bottom friction once for all 
     
    346380               END DO 
    347381            END DO 
     382            IF ( ln_isfcav ) THEN 
     383               DO jj = 1, jpj 
     384                  DO ji = 1, jpi 
     385                     ikbt = mikt(ji,jj) 
     386                     ztmp = tmask(ji,jj,ikbt) * ( vkarmn / LOG( 0.5_wp * fse3t_n(ji,jj,ikbt) / rn_tfrz0 ))**2._wp 
     387                     tfrcoef2d(ji,jj) = MAX(tfrcoef2d(ji,jj), ztmp) 
     388                     tfrcoef2d(ji,jj) = MIN(tfrcoef2d(ji,jj), rn_tfri2_max) 
     389                  END DO 
     390               END DO 
     391            END IF 
    348392         ENDIF 
    349393         ! 
     
    398442             zminbfr = MIN(  zminbfr, MIN( zfru, ABS( bfrcoef2d(ji,jj) ) )  ) 
    399443             zmaxbfr = MAX(  zmaxbfr, MIN( zfrv, ABS( bfrcoef2d(ji,jj) ) )  ) 
     444! (ISF) 
     445             IF ( ln_isfcav ) THEN 
     446                ikbu = miku(ji,jj)       ! 1st wet ocean level at u- and v-points 
     447                ikbv = mikv(ji,jj) 
     448                zfru = 0.5 * fse3u(ji,jj,ikbu) / rdt 
     449                zfrv = 0.5 * fse3v(ji,jj,ikbv) / rdt 
     450                IF( ABS( tfrcoef2d(ji,jj) ) > zfru ) THEN 
     451                   IF( ln_ctl ) THEN 
     452                      WRITE(numout,*) 'TFR ', narea, nimpp+ji, njmpp+jj, ikbu 
     453                      WRITE(numout,*) 'TFR ', ABS( tfrcoef2d(ji,jj) ), zfru 
     454                   ENDIF 
     455                   ictu = ictu + 1 
     456                ENDIF 
     457                IF( ABS( tfrcoef2d(ji,jj) ) > zfrv ) THEN 
     458                   IF( ln_ctl ) THEN 
     459                      WRITE(numout,*) 'TFR ', narea, nimpp+ji, njmpp+jj, ikbv 
     460                      WRITE(numout,*) 'TFR ', tfrcoef2d(ji,jj), zfrv 
     461                   ENDIF 
     462                   ictv = ictv + 1 
     463                ENDIF 
     464                zmintfr = MIN(  zmintfr, MIN( zfru, ABS( tfrcoef2d(ji,jj) ) )  ) 
     465                zmaxtfr = MAX(  zmaxtfr, MIN( zfrv, ABS( tfrcoef2d(ji,jj) ) )  ) 
     466             END IF 
     467! END ISF 
    400468         END DO 
    401469      END DO 
     
    405473         CALL mpp_min( zminbfr ) 
    406474         CALL mpp_max( zmaxbfr ) 
     475         IF ( ln_isfcav) CALL mpp_min( zmintfr ) 
     476         IF ( ln_isfcav) CALL mpp_max( zmaxtfr ) 
    407477      ENDIF 
    408478      IF( .NOT.ln_bfrimp) THEN 
    409479      IF( lwp .AND. ictu + ictv > 0 ) THEN 
    410          WRITE(numout,*) ' Bottom friction stability check failed at ', ictu, ' U-points ' 
    411          WRITE(numout,*) ' Bottom friction stability check failed at ', ictv, ' V-points ' 
     480         WRITE(numout,*) ' Bottom/Top friction stability check failed at ', ictu, ' U-points ' 
     481         WRITE(numout,*) ' Bottom/Top friction stability check failed at ', ictv, ' V-points ' 
    412482         WRITE(numout,*) ' Bottom friction coefficient now ranges from: ', zminbfr, ' to ', zmaxbfr 
    413          WRITE(numout,*) ' Bottom friction coefficient now ranges from: ', zmintfr, ' to ', zmaxtfr 
    414          WRITE(numout,*) ' Bottom friction coefficient will be reduced where necessary' 
     483         IF ( ln_isfcav ) WRITE(numout,*) ' Top friction coefficient now ranges from: ', zmintfr, ' to ', zmaxtfr 
     484         WRITE(numout,*) ' Bottom/Top friction coefficient will be reduced where necessary' 
    415485      ENDIF 
    416486      ENDIF 
Note: See TracChangeset for help on using the changeset viewer.