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 5111 for branches/2015/dev_r5094_UKMO_ISFCLEAN/NEMOGCM/NEMO/OPA_SRC – NEMO

Ignore:
Timestamp:
2015-03-02T16:30:57+01:00 (9 years ago)
Author:
mathiot
Message:

add some missing if ln_isfcav, test of option compatibility with ln_isfcav + small documentation

Location:
branches/2015/dev_r5094_UKMO_ISFCLEAN/NEMOGCM/NEMO/OPA_SRC
Files:
18 edited

Legend:

Unmodified
Added
Removed
  • branches/2015/dev_r5094_UKMO_ISFCLEAN/NEMOGCM/NEMO/OPA_SRC/DIA/diaar5.F90

    r4990 r5111  
    105105      END DO 
    106106      IF( .NOT.lk_vvl ) THEN 
    107          DO ji=1,jpi 
    108             DO jj=1,jpj 
    109                zbotpres(ji,jj) = zbotpres(ji,jj) + sshn(ji,jj) * zrhd(ji,jj,mikt(ji,jj)) + riceload(ji,jj) 
    110             END DO 
    111          END DO 
     107         IF ( ln_isfcav ) THEN 
     108            DO ji=1,jpi 
     109               DO jj=1,jpj 
     110                  zbotpres(ji,jj) = zbotpres(ji,jj) + sshn(ji,jj) * zrhd(ji,jj,mikt(ji,jj)) + riceload(ji,jj) 
     111               END DO 
     112            END DO 
     113         ELSE 
     114            zbotpres(:,:) = zbotpres(:,:) + sshn(:,:) * zrhd(:,:,1) 
     115         END IF 
    112116      END IF 
    113117      !                                          
     
    127131      END DO 
    128132      IF( .NOT.lk_vvl ) THEN 
    129          DO ji=1,jpi 
    130             DO jj=1,jpj 
    131                zbotpres(ji,jj) = zbotpres(ji,jj) + sshn(ji,jj) * zrhd(ji,jj,mikt(ji,jj)) + riceload(ji,jj) 
    132             END DO 
    133          END DO 
     133         IF ( ln_isfcav ) THEN 
     134            DO ji=1,jpi 
     135               DO jj=1,jpj 
     136                  zbotpres(ji,jj) = zbotpres(ji,jj) + sshn(ji,jj) * zrhd(ji,jj,mikt(ji,jj)) + riceload(ji,jj) 
     137               END DO 
     138            END DO 
     139         ELSE 
     140            zbotpres(:,:) = zbotpres(:,:) + sshn(:,:) * zrhd(:,:,1) 
     141         END IF 
    134142      END IF 
    135143      !     
     
    157165      END DO 
    158166      IF( .NOT.lk_vvl ) THEN 
    159          DO ji=1,jpi 
    160             DO jj=1,jpj 
    161                ztemp = ztemp + zarea_ssh(ji,jj) * tsn(ji,jj,mikt(ji,jj),jp_tem)  
    162                zsal  = zsal  + zarea_ssh(ji,jj) * tsn(ji,jj,mikt(ji,jj),jp_sal)  
    163             END DO 
    164          END DO 
     167         IF ( ln_isfcav ) THEN 
     168            DO ji=1,jpi 
     169               DO jj=1,jpj 
     170                  ztemp = ztemp + zarea_ssh(ji,jj) * tsn(ji,jj,mikt(ji,jj),jp_tem)  
     171                  zsal  = zsal  + zarea_ssh(ji,jj) * tsn(ji,jj,mikt(ji,jj),jp_sal)  
     172               END DO 
     173            END DO 
     174         ELSE 
     175            ztemp = ztemp + zarea_ssh(:,:) * tsn(:,:,1,jp_tem)  
     176            zsal  = zsal  + zarea_ssh(:,:) * tsn(:,:,1,jp_sal)  
     177         END IF 
    165178      ENDIF 
    166179      IF( lk_mpp ) THEN   
  • branches/2015/dev_r5094_UKMO_ISFCLEAN/NEMOGCM/NEMO/OPA_SRC/DIA/diahsb.F90

    r4990 r5111  
    9696      z_frc_trd_t =           glob_sum( sbc_tsc(:,:,jp_tem) * surf(:,:) )                               ! heat fluxes 
    9797      z_frc_trd_s =           glob_sum( sbc_tsc(:,:,jp_sal) * surf(:,:) )                               ! salt fluxes 
    98       ! Add runoff heat & salt input 
     98      ! Add runoff    heat & salt input 
    9999      IF( ln_rnf    )   z_frc_trd_t = z_frc_trd_t + glob_sum( rnf_tsc(:,:,jp_tem) * surf(:,:) ) 
    100100      IF( ln_rnf_sal)   z_frc_trd_s = z_frc_trd_s + glob_sum( rnf_tsc(:,:,jp_sal) * surf(:,:) ) 
    101       ! Add geothermal ice shelf 
     101      ! Add ice shelf heat & salt input 
    102102      IF( nn_isf .GE. 1 )  THEN 
    103103          z_frc_trd_t = z_frc_trd_t & 
     
    112112      ! 
    113113      IF( .NOT. lk_vvl ) THEN 
    114          z2d0=0.0_wp ; z2d1=0.0_wp 
    115          DO ji=1,jpi 
    116             DO jj=1,jpj 
    117               z2d0(ji,jj) = surf(ji,jj) * wn(ji,jj,mikt(ji,jj)) * tsb(ji,jj,mikt(ji,jj),jp_tem) 
    118               z2d1(ji,jj) = surf(ji,jj) * wn(ji,jj,mikt(ji,jj)) * tsb(ji,jj,mikt(ji,jj),jp_sal) 
     114         IF ( ln_isfcav ) THEN 
     115            DO ji=1,jpi 
     116               DO jj=1,jpj 
     117                  z2d0(ji,jj) = surf(ji,jj) * wn(ji,jj,mikt(ji,jj)) * tsb(ji,jj,mikt(ji,jj),jp_tem) 
     118                  z2d1(ji,jj) = surf(ji,jj) * wn(ji,jj,mikt(ji,jj)) * tsb(ji,jj,mikt(ji,jj),jp_sal) 
     119               ENDDO 
    119120            ENDDO 
    120          ENDDO 
     121         ELSE 
     122            z2d0(:,:) = surf(:,:) * wn(:,:,1) * tsb(:,:,1,jp_tem) 
     123            z2d1(:,:) = surf(:,:) * wn(:,:,1) * tsb(:,:,1,jp_sal) 
     124         END IF 
    121125         z_wn_trd_t = - glob_sum( z2d0 )  
    122126         z_wn_trd_s = - glob_sum( z2d1 ) 
     
    144148      ! heat & salt content variation (associated with ssh) 
    145149      IF( .NOT. lk_vvl ) THEN 
    146          z2d0 = 0._wp   ;   z2d1 = 0._wp 
    147          DO ji = 1, jpi 
    148             DO jj = 1, jpj 
    149               z2d0(ji,jj) = surf(ji,jj) * ( tsn(ji,jj,mikt(ji,jj),jp_tem) * sshn(ji,jj) - ssh_hc_loc_ini(ji,jj) )  
    150               z2d1(ji,jj) = surf(ji,jj) * ( tsn(ji,jj,mikt(ji,jj),jp_sal) * sshn(ji,jj) - ssh_sc_loc_ini(ji,jj) )  
     150         IF ( ln_isfcav ) THEN 
     151            DO ji = 1, jpi 
     152               DO jj = 1, jpj 
     153                  z2d0(ji,jj) = surf(ji,jj) * ( tsn(ji,jj,mikt(ji,jj),jp_tem) * sshn(ji,jj) - ssh_hc_loc_ini(ji,jj) )  
     154                  z2d1(ji,jj) = surf(ji,jj) * ( tsn(ji,jj,mikt(ji,jj),jp_sal) * sshn(ji,jj) - ssh_sc_loc_ini(ji,jj) )  
     155               END DO 
    151156            END DO 
    152          END DO 
     157         ELSE 
     158            z2d0(:,:) = surf(:,:) * ( tsn(:,:,1,jp_tem) * sshn(:,:) - ssh_hc_loc_ini(:,:) )  
     159            z2d1(:,:) = surf(:,:) * ( tsn(:,:,1,jp_sal) * sshn(:,:) - ssh_sc_loc_ini(:,:) )  
     160         END IF 
    153161         z_ssh_hc = glob_sum( z2d0 )  
    154162         z_ssh_sc = glob_sum( z2d1 )  
     
    277285          frc_s = 0._wp                                           ! salt content   -    -   -    -         
    278286          IF( .NOT. lk_vvl ) THEN 
    279              DO ji=1,jpi 
    280                 DO jj=1,jpj 
    281                    ssh_hc_loc_ini(ji,jj) = tsn(ji,jj,mikt(ji,jj),jp_tem) * sshn(ji,jj)   ! initial heat content in ssh 
    282                    ssh_sc_loc_ini(ji,jj) = tsn(ji,jj,mikt(ji,jj),jp_sal) * sshn(ji,jj)   ! initial salt content in ssh 
     287             IF ( ln_isfcav ) THEN 
     288                DO ji=1,jpi 
     289                   DO jj=1,jpj 
     290                      ssh_hc_loc_ini(ji,jj) = tsn(ji,jj,mikt(ji,jj),jp_tem) * sshn(ji,jj)   ! initial heat content in ssh 
     291                      ssh_sc_loc_ini(ji,jj) = tsn(ji,jj,mikt(ji,jj),jp_sal) * sshn(ji,jj)   ! initial salt content in ssh 
     292                   ENDDO 
    283293                ENDDO 
    284              ENDDO 
     294             ELSE 
     295                ssh_hc_loc_ini(:,:) = tsn(:,:,1,jp_tem) * sshn(:,:)   ! initial heat content in ssh 
     296                ssh_sc_loc_ini(:,:) = tsn(:,:,1,jp_sal) * sshn(:,:)   ! initial salt content in ssh 
     297             END IF 
    285298             frc_wn_t = 0._wp                                       ! initial heat content misfit due to free surface 
    286299             frc_wn_s = 0._wp                                       ! initial salt content misfit due to free surface 
  • branches/2015/dev_r5094_UKMO_ISFCLEAN/NEMOGCM/NEMO/OPA_SRC/DOM/domzgr.F90

    r5098 r5111  
    298298      ENDIF 
    299299 
     300!      IF ( ln_isfcav ) THEN 
    300301! need to be like this to compute the pressure gradient with ISF. If not, level beneath the ISF are not aligned (sum(e3t) /= depth) 
    301302! define e3t_0 and e3w_0 as the differences between gdept and gdepw respectively 
    302       DO jk = 1, jpkm1 
    303          e3t_1d(jk) = gdepw_1d(jk+1)-gdepw_1d(jk)  
    304       END DO 
    305       e3t_1d(jpk) = e3t_1d(jpk-1)   ! we don't care because this level is masked in NEMO 
    306  
    307       DO jk = 2, jpk 
    308          e3w_1d(jk) = gdept_1d(jk) - gdept_1d(jk-1)  
    309       END DO 
    310       e3w_1d(1  ) = 2._wp * (gdept_1d(1) - gdepw_1d(1))  
     303         DO jk = 1, jpkm1 
     304            e3t_1d(jk) = gdepw_1d(jk+1)-gdepw_1d(jk)  
     305         END DO 
     306         e3t_1d(jpk) = e3t_1d(jpk-1)   ! we don't care because this level is masked in NEMO 
     307 
     308         DO jk = 2, jpk 
     309            e3w_1d(jk) = gdept_1d(jk) - gdept_1d(jk-1)  
     310         END DO 
     311         e3w_1d(1  ) = 2._wp * (gdept_1d(1) - gdepw_1d(1))  
     312!      END IF 
    311313 
    312314!!gm BUG in s-coordinate this does not work! 
     
    473475         ! (ISF) TODO build ice draft netcdf file for isomip and build the corresponding part of code 
    474476         IF( cp_cfg == "isomip" .AND. ln_isfcav ) THEN  
    475            !  
    476            risfdep(:,:)=200.e0  
    477            misfdep(:,:)=1  
    478            ij0 = 1 ; ij1 = 40  
    479            DO jj = mj0(ij0), mj1(ij1)  
    480               risfdep(:,jj)=700.0_wp-(gphit(:,jj)+80.0_wp)*125.0_wp  
    481                 END DO  
     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  
    482483            WHERE( bathy(:,:) <= 0._wp )  risfdep(:,:) = 0._wp  
    483            !  
     484         !  
    484485         ELSEIF ( cp_cfg == "isomip2" .AND. ln_isfcav ) THEN 
    485486         !  
     
    536537            CALL iom_get  ( inum, jpdom_data, 'Bathymetry', bathy ) 
    537538            CALL iom_close( inum ) 
    538             !   
     539            !                                                 
    539540            risfdep(:,:)=0._wp          
    540541            misfdep(:,:)=1              
     
    963964      !!---------------------------------------------------------------------- 
    964965      !! 
    965       INTEGER  ::   ji, jj, jk, jl   ! dummy loop indices 
     966      INTEGER  ::   ji, jj, jk       ! dummy loop indices 
    966967      INTEGER  ::   ik, it           ! temporary integers 
    967       INTEGER  ::   id, jd, nprocd 
    968968      LOGICAL  ::   ll_print         ! Allow  control print for debugging 
    969969      REAL(wp) ::   ze3tp , ze3wp    ! Last ocean level thickness at T- and W-points 
    970970      REAL(wp) ::   zdepwp, zdepth   ! Ajusted ocean depth to avoid too small e3t 
    971       REAL(wp) ::   zmax, zmin       ! Maximum and minimum depth 
     971      REAL(wp) ::   zmax             ! Maximum depth 
    972972      REAL(wp) ::   zdiff            ! temporary scalar 
    973973      REAL(wp) ::   zrefdep          ! temporary scalar 
     
    10461046                     &                             * ((gdept_1d(     ik  ) - gdepw_1d(ik) )   & 
    10471047                     &                             / ( gdepw_1d(     ik+1) - gdepw_1d(ik) )) 
    1048                   e3t_0(ji,jj,ik) = e3t_1d (ik) * ( gdepw_0 (ji,jj,ik+1) - gdepw_1d(ik) )   &  
    1049                      &                          / ( gdepw_1d(      ik+1) - gdepw_1d(ik) )  
     1048                  e3t_0  (ji,jj,ik) = e3t_1d (ik) * ( gdepw_0 (ji,jj,ik+1) - gdepw_1d(ik) )   &  
     1049                     &                             / ( gdepw_1d(      ik+1) - gdepw_1d(ik) )  
    10501050                  e3w_0(ji,jj,ik) = 0.5_wp * ( gdepw_0(ji,jj,ik+1) + gdepw_1d(ik+1) - 2._wp * gdepw_1d(ik) )   & 
    10511051                     &                     * ( e3w_1d(ik) / ( gdepw_1d(ik+1) - gdepw_1d(ik) ) ) 
     
    12021202      
    12031203      ! Compute gdep3w_0 (vertical sum of e3w) 
    1204       IF ( ln_isfcav ) THEN 
     1204      IF ( ln_isfcav ) THEN ! if cavity 
    12051205         WHERE (misfdep == 0) misfdep = 1 
    12061206         DO jj = 1,jpj 
     
    13351335            misfdep(jpi,:) = misfdep(  2  ,:)  
    13361336         ENDIF 
    1337   
     1337 
    13381338         IF( nperio == 1 .OR. nperio  ==  4 .OR. nperio  ==  6 ) THEN 
    13391339            mbathy( 1 ,:) = mbathy(jpim1,:)             ! local domain is cyclic east-west 
    13401340            mbathy(jpi,:) = mbathy(  2  ,:) 
    13411341         ENDIF 
    1342   
     1342 
    13431343         ! split last cell if possible (only where water column is 2 cell or less) 
    13441344         DO jk = jpkm1, 1, -1 
     
    13581358            END WHERE 
    13591359         END DO 
    1360   
     1360 
    13611361  
    13621362 ! Case where bathy and risfdep compatible but not the level variable mbathy/misfdep because of partial cell condition 
     
    16391639               IF( zmbathy(ji,jj) .LT. misfdep(ji  ,jj+1) ) ibtestjp1 = 0 
    16401640               ibtest=MAX(ibtestim1, ibtestip1, ibtestjm1, ibtestjp1) 
    1641                IF( ibtest == 0 ) THEN 
     1641               IF( ibtest == 0 .AND. misfdep(ji,jj) .GE. 2) THEN 
    16421642                  mbathy(ji,jj) = 0 ; bathy(ji,jj) = 0.0_wp ; misfdep(ji,jj) = 0 ; risfdep(ji,jj) = 0.0_wp ; 
    16431643               END IF 
  • branches/2015/dev_r5094_UKMO_ISFCLEAN/NEMOGCM/NEMO/OPA_SRC/DOM/istate.F90

    r5104 r5111  
    140140            &            CALL zps_hde    ( nit000, jpts, tsb, gtsu, gtsv,  &    ! Partial steps: before horizontal gradient 
    141141            &                                            rhd, gru , grv    )  ! of t, s, rd at the last ocean level 
    142          IF( ln_zps .AND. ln_isfcav)                                       & 
     142         IF( ln_zps .AND.       ln_isfcav)                                 & 
    143143            &            CALL zps_hde_isf( nit000, jpts, tsb, gtsu, gtsv,  &    ! Partial steps for top cell (ISF) 
    144144            &                                            rhd, gru , grv , aru , arv , gzu , gzv , ge3ru , ge3rv ,   & 
  • branches/2015/dev_r5094_UKMO_ISFCLEAN/NEMOGCM/NEMO/OPA_SRC/DYN/dynadv.F90

    r4990 r5111  
    127127      IF( ln_dynzad_zts .AND. .NOT. ln_dynadv_vec )   & 
    128128          CALL ctl_stop( 'Sub timestepping of vertical advection requires vector form; set ln_dynadv_vec = .TRUE.' ) 
     129      IF( ln_dynzad_zts .AND. ln_isfcav )   & 
     130          CALL ctl_stop( 'Sub timestepping of vertical advection does not work with ln_isfcav = .TRUE.' ) 
    129131 
    130132      !                               ! Set nadv 
  • branches/2015/dev_r5094_UKMO_ISFCLEAN/NEMOGCM/NEMO/OPA_SRC/DYN/dynhpg.F90

    r5098 r5111  
    167167                           & the standard jacobian formulation hpg_sco or & 
    168168                           & the pressure jacobian formulation hpg_prj') 
     169 
     170      IF( ln_hpg_isf .AND. .NOT. ln_isfcav )   & 
     171           &   CALL ctl_stop( ' hpg_isf not available is ln_isfcav = false ' ) 
    169172      ! 
    170173      !                               ! Set nhpg from ln_hpg_... flags 
     
    187190      IF( ( .NOT. ln_hpg_isf ) .AND. ln_isfcav )   & 
    188191          &  CALL ctl_stop( 'Only hpg_isf has been corrected to work with ice shelf cavity.' ) 
     192      !  
     193      ! initialisation of ice load 
     194      riceload(:,:)=0.0 
    189195      ! 
    190196   END SUBROUTINE dyn_hpg_init 
  • branches/2015/dev_r5094_UKMO_ISFCLEAN/NEMOGCM/NEMO/OPA_SRC/DYN/dynzad.F90

    r5098 r5111  
    9090         DO jj = 2, jpjm1                 ! vertical momentum advection at w-point 
    9191            DO ji = fs_2, fs_jpim1        ! vector opt. 
    92                zwuw(ji,jj,jk) = ( zww(ji+1,jj  ) + zww(ji,jj) ) * ( un(ji,jj,jk-1)-un(ji,jj,jk) ) !* wumask(ji,jj,jk) 
    93                zwvw(ji,jj,jk) = ( zww(ji  ,jj+1) + zww(ji,jj) ) * ( vn(ji,jj,jk-1)-vn(ji,jj,jk) ) !* wvmask(ji,jj,jk) 
     92               zwuw(ji,jj,jk) = ( zww(ji+1,jj  ) + zww(ji,jj) ) * ( un(ji,jj,jk-1)-un(ji,jj,jk) ) 
     93               zwvw(ji,jj,jk) = ( zww(ji  ,jj+1) + zww(ji,jj) ) * ( vn(ji,jj,jk-1)-vn(ji,jj,jk) ) 
    9494            END DO   
    9595         END DO    
    9696      END DO 
    97       DO jj = 2, jpjm1              ! Surface and bottom values set to zero 
    98          DO ji = fs_2, fs_jpim1           ! vector opt. 
    99             zwuw(ji,jj, 1:miku(ji,jj) ) = 0._wp 
    100             zwvw(ji,jj, 1:mikv(ji,jj) ) = 0._wp 
    101             zwuw(ji,jj,jpk) = 0._wp 
    102             zwvw(ji,jj,jpk) = 0._wp 
    103          END DO   
    104       END DO 
     97      ! 
     98      ! Surface and bottom advective fluxes set to zero 
     99      IF ( ln_isfcav ) THEN 
     100         DO jj = 2, jpjm1 
     101            DO ji = fs_2, fs_jpim1           ! vector opt. 
     102               zwuw(ji,jj, 1:miku(ji,jj) ) = 0._wp 
     103               zwvw(ji,jj, 1:mikv(ji,jj) ) = 0._wp 
     104               zwuw(ji,jj,jpk) = 0._wp 
     105               zwvw(ji,jj,jpk) = 0._wp 
     106            END DO 
     107         END DO 
     108      ELSE 
     109         DO jj = 2, jpjm1         
     110            DO ji = fs_2, fs_jpim1           ! vector opt. 
     111               zwuw(ji,jj, 1 ) = 0._wp 
     112               zwvw(ji,jj, 1 ) = 0._wp 
     113               zwuw(ji,jj,jpk) = 0._wp 
     114               zwvw(ji,jj,jpk) = 0._wp 
     115            END DO   
     116         END DO 
     117      END IF 
    105118 
    106119      DO jk = 1, jpkm1              ! Vertical momentum advection at u- and v-points 
     
    196209         END DO 
    197210      END DO 
    198  
    199       DO jj = 2, jpjm1                    ! Surface and bottom advective fluxes set to zero 
     211      ! 
     212      ! Surface and bottom advective fluxes set to zero 
     213      DO jj = 2, jpjm1         
    200214         DO ji = fs_2, fs_jpim1           ! vector opt. 
    201             zwuw(ji,jj, 1:miku(ji,jj) ) = 0._wp 
    202             zwvw(ji,jj, 1:mikv(ji,jj) ) = 0._wp 
     215            zwuw(ji,jj, 1 ) = 0._wp 
     216            zwvw(ji,jj, 1 ) = 0._wp 
    203217            zwuw(ji,jj,jpk) = 0._wp 
    204218            zwvw(ji,jj,jpk) = 0._wp 
  • branches/2015/dev_r5094_UKMO_ISFCLEAN/NEMOGCM/NEMO/OPA_SRC/LDF/ldfslp.F90

    r5104 r5111  
    157157         ! 
    158158         !==   Local vertical density gradient at T-point   == !   (evaluated from N^2) 
    159          ! surface initialisation  
    160          zdzr(:,:,1) = 0._wp  
    161159         ! interior value 
    162160         DO jk = 2, jpkm1 
     
    169167               &                 * ( pn2(:,:,jk) + pn2(:,:,jk+1) ) * ( 1._wp - 0.5_wp * tmask(:,:,jk+1) ) 
    170168         END DO 
     169         ! surface initialisation  
     170         zdzr(:,:,1) = 0._wp  
    171171         IF ( ln_isfcav ) THEN 
    172          ! if isf need to overwrite the interior value at at the first ocean point 
     172            ! if isf need to overwrite the interior value at at the first ocean point 
    173173            DO jj = 1, jpjm1 
    174174               DO ji = 1, jpim1 
     
    185185         ! ===========================      | vslp = d/dj( prd ) / d/dz( prd ) 
    186186         ! 
    187          DO jj = 2, jpjm1 
    188             DO ji = fs_2, fs_jpim1   ! vector opt. 
    189                IF (miku(ji,jj) .GT. miku(ji+1,jj)) zhmlpu(ji,jj) = hmlpt(ji  ,jj) 
    190                IF (miku(ji,jj) .LT. miku(ji+1,jj)) zhmlpu(ji,jj) = hmlpt(ji+1,jj) 
    191                IF (miku(ji,jj) .EQ. miku(ji+1,jj)) zhmlpu(ji,jj) = MAX(hmlpt(ji  ,jj), hmlpt(ji+1,jj)) 
    192                IF (mikv(ji,jj) .GT. miku(ji,jj+1)) zhmlpv(ji,jj) = hmlpt(ji  ,jj) 
    193                IF (mikv(ji,jj) .LT. miku(ji,jj+1)) zhmlpv(ji,jj) = hmlpt(ji,jj+1) 
    194                IF (mikv(ji,jj) .EQ. miku(ji,jj+1)) zhmlpv(ji,jj) = MAX(hmlpt(ji,jj), hmlpt(ji,jj+1)) 
     187         IF ( ln_isfcav ) THEN 
     188            DO jj = 2, jpjm1 
     189               DO ji = fs_2, fs_jpim1   ! vector opt. 
     190                  IF (miku(ji,jj) .GT. miku(ji+1,jj)) zhmlpu(ji,jj) = MAX(hmlpt(ji  ,jj  ),                   5._wp) 
     191                  IF (miku(ji,jj) .LT. miku(ji+1,jj)) zhmlpu(ji,jj) = MAX(hmlpt(ji+1,jj  ),                   5._wp) 
     192                  IF (miku(ji,jj) .EQ. miku(ji+1,jj)) zhmlpu(ji,jj) = MAX(hmlpt(ji  ,jj  ), hmlpt(ji+1,jj  ), 5._wp) 
     193                  IF (mikv(ji,jj) .GT. miku(ji,jj+1)) zhmlpv(ji,jj) = MAX(hmlpt(ji  ,jj  ),                   5._wp) 
     194                  IF (mikv(ji,jj) .LT. miku(ji,jj+1)) zhmlpv(ji,jj) = MAX(hmlpt(ji  ,jj+1),                   5._wp) 
     195                  IF (mikv(ji,jj) .EQ. miku(ji,jj+1)) zhmlpv(ji,jj) = MAX(hmlpt(ji  ,jj  ), hmlpt(ji  ,jj+1), 5._wp) 
     196               ENDDO 
    195197            ENDDO 
    196          ENDDO 
     198         ELSE 
     199            DO jj = 2, jpjm1 
     200               DO ji = fs_2, fs_jpim1   ! vector opt. 
     201                  zhmlpu(ji,jj) = MAX(hmlpt(ji,jj), hmlpt(ji+1,jj  ), 5._wp) 
     202                  zhmlpv(ji,jj) = MAX(hmlpt(ji,jj), hmlpt(ji  ,jj+1), 5._wp) 
     203               ENDDO 
     204            ENDDO 
     205         END IF 
    197206         DO jk = 2, jpkm1                            !* Slopes at u and v points 
    198207            DO jj = 2, jpjm1 
     
    208217                  zbv = MIN(  zbv, -100._wp* ABS( zav ) , -7.e+3_wp/fse3v(ji,jj,jk)* ABS( zav )  ) 
    209218                  !                                      ! uslp and vslp output in zwz and zww, resp. 
    210                   zfi = MAX( omlmask(ji,jj,jk), omlmask(ji+1,jj,jk) )  
    211                   zfj = MAX( omlmask(ji,jj,jk), omlmask(ji,jj+1,jk) )  
     219                  zfi = MAX( omlmask(ji,jj,jk), omlmask(ji+1,jj  ,jk) ) 
     220                  zfj = MAX( omlmask(ji,jj,jk), omlmask(ji  ,jj+1,jk) ) 
    212221                  ! thickness of water column between surface and level k at u/v point 
    213                   zdepu = 0.5_wp * (( fsdept(ji,jj,jk) + fsdept(ji+1,jj  ,jk) )                   & 
    214                              - 2 * MAX( risfdep(ji,jj), risfdep(ji+1,jj  ) )  & 
    215                              - fse3u(ji,jj,miku(ji,jj))                                         ) 
    216                   zdepv = 0.5_wp * (( fsdept(ji,jj,jk) + fsdept(ji  ,jj+1,jk) )                   & 
    217                              - 2 * MAX( risfdep(ji,jj), risfdep(ji,jj+1) ) & 
    218                              - fse3v(ji,jj,mikv(ji,jj))                                         ) 
    219                   zwz(ji,jj,jk) = ( 1. - zfi) * zau / ( zbu - zeps )                                              & 
    220                      &                 + zfi  * uslpml(ji,jj)                                                     & 
    221                      &                        * zdepu / MAX( zhmlpu(ji,jj), 5._wp ) 
     222                  zdepu = 0.5_wp * ( ( fsdept(ji,jj,jk) + fsdept(ji+1,jj  ,jk) )                              & 
     223                                   - 2 * MAX( risfdep(ji,jj), risfdep(ji+1,jj  ) ) - fse3u(ji,jj,miku(ji,jj)) ) 
     224                  zdepv = 0.5_wp * ( ( fsdept(ji,jj,jk) + fsdept(ji  ,jj+1,jk) )                              & 
     225                                   - 2 * MAX( risfdep(ji,jj), risfdep(ji  ,jj+1) ) - fse3v(ji,jj,mikv(ji,jj)) ) 
     226                  ! 
     227                  zwz(ji,jj,jk) = ( 1. - zfi) * zau / ( zbu - zeps )                                          & 
     228                     &                 + zfi  * uslpml(ji,jj) * zdepu / zhmlpu(ji,jj) 
    222229                  zwz(ji,jj,jk) = zwz(ji,jj,jk) * wumask(ji,jj,jk) 
    223                   zww(ji,jj,jk) = ( 1. - zfj) * zav / ( zbv - zeps )                                              & 
    224                      &                 + zfj  * vslpml(ji,jj)                                                     & 
    225                      &                        * zdepv / MAX( zhmlpv(ji,jj), 5._wp )  
     230                  zww(ji,jj,jk) = ( 1. - zfj) * zav / ( zbv - zeps )                                          & 
     231                     &                 + zfj  * vslpml(ji,jj) * zdepv / zhmlpv(ji,jj)  
    226232                  zww(ji,jj,jk) = zww(ji,jj,jk) * wvmask(ji,jj,jk) 
    227233                   
     
    276282                  uslp(ji,jj,jk) = uslp(ji,jj,jk) * ( umask(ji,jj+1,jk) + umask(ji,jj-1,jk  ) ) * 0.5_wp   & 
    277283                     &                            * ( umask(ji,jj  ,jk) + umask(ji,jj  ,jk+1) ) * 0.5_wp   & 
    278                      &                            *   umask(ji,jj,jk-1) !* umask(ji,jj,jk) * umask(ji,jj,jk+1) 
     284                     &                            *   umask(ji,jj,jk-1) 
    279285                  vslp(ji,jj,jk) = vslp(ji,jj,jk) * ( vmask(ji+1,jj,jk) + vmask(ji-1,jj,jk  ) ) * 0.5_wp   & 
    280286                     &                            * ( vmask(ji  ,jj,jk) + vmask(ji  ,jj,jk+1) ) * 0.5_wp   & 
    281                      &                            *   vmask(ji,jj,jk-1) !* vmask(ji,jj,jk) * vmask(ji,jj,jk+1) 
     287                     &                            *   vmask(ji,jj,jk-1) 
    282288               END DO 
    283289            END DO 
     
    292298               DO ji = fs_2, fs_jpim1   ! vector opt. 
    293299                  !                                  !* Local vertical density gradient evaluated from N^2 
    294                   zbw = zm1_2g * pn2 (ji,jj,jk) * ( prd (ji,jj,jk) + prd (ji,jj,jk-1) + 2. ) * tmask(ji,jj,jk) * tmask(ji,jj,jk-1) 
     300                  zbw = zm1_2g * pn2 (ji,jj,jk) * ( prd (ji,jj,jk) + prd (ji,jj,jk-1) + 2. ) * wmask(ji,jj,jk) 
    295301                  !                                  !* Slopes at w point 
    296302                  !                                        ! i- & j-gradient of density at w-points 
     
    308314                  zbj = MIN( zbw , -100._wp* ABS( zaj ) , -7.e+3_wp/fse3w(ji,jj,jk)* ABS( zaj )  ) 
    309315                  !                                        ! wslpi and wslpj with ML flattening (output in zwz and zww, resp.) 
    310                   zfk = MAX( omlmask(ji,jj,jk), omlmask(ji,jj,jk-1) )    ! zfk=1 in the ML otherwise zfk=0 
     316                  zfk = MAX( omlmask(ji,jj,jk), omlmask(ji,jj,jk-1) )   ! zfk=1 in the ML otherwise zfk=0 
    311317                  zck = ( fsdepw(ji,jj,jk) - fsdepw(ji,jj,mikt(ji,jj) ) ) / MAX( hmlp(ji,jj), 10._wp ) 
    312318                  zwz(ji,jj,jk) = (  zai / ( zbi - zeps ) * ( 1._wp - zfk ) & 
     
    746752            DO ji = 1, jpi 
    747753               ik = nmln(ji,jj) - 1 
    748                IF( jk <= ik .AND. jk >= mikt(ji,jj) ) THEN   ;   omlmask(ji,jj,jk) = 1._wp 
    749                ELSE                  ;   omlmask(ji,jj,jk) = 0._wp 
     754               IF( jk <= ik .AND. jk >= mikt(ji,jj) ) THEN 
     755                  omlmask(ji,jj,jk) = 1._wp 
     756               ELSE 
     757                  omlmask(ji,jj,jk) = 0._wp 
    750758               ENDIF 
    751759            END DO 
  • branches/2015/dev_r5094_UKMO_ISFCLEAN/NEMOGCM/NEMO/OPA_SRC/SBC/sbcfwb.F90

    r5098 r5111  
    8989         ! 
    9090         IF( kn_fwb == 3 .AND. nn_sssr /= 2 )   CALL ctl_stop( 'sbc_fwb: nn_fwb = 3 requires nn_sssr = 2, we stop ' ) 
    91          ! 
    92          area = glob_sum( e1e2t(:,:) )           ! interior global domain surface 
     91         IF( kn_fwb == 3 .AND. ln_isfcav    )   CALL ctl_stop( 'sbc_fwb: nn_fwb = 3 with ln_isfcav = .TRUE. not working, we stop ' ) 
     92         ! 
     93         area = glob_sum( e1e2t(:,:) * tmask(:,:,1))           ! interior global domain surface 
     94         ! isf cavities are excluded because it can feedback to the melting with generation of inhibition of plumes 
     95         ! and in case of no melt, it can generate HSSW. 
    9396         ! 
    9497#if ! defined key_lim2 &&  ! defined key_lim3 && ! defined key_cice 
     
    107110            z_fwf = glob_sum( e1e2t(:,:) * ( emp(:,:) - rnf(:,:) + rdivisf * fwfisf(:,:) -  snwice_fmass(:,:) ) ) / area   ! sum over the global domain 
    108111            zcoef = z_fwf * rcp 
    109             emp(:,:) = emp(:,:) - z_fwf  
    110             qns(:,:) = qns(:,:) + zcoef * sst_m(:,:) ! account for change to the heat budget due to fw correction 
     112            emp(:,:) = emp(:,:) - z_fwf              * tmask(:,:,1) 
     113            qns(:,:) = qns(:,:) + zcoef * sst_m(:,:) * tmask(:,:,1) ! account for change to the heat budget due to fw correction 
    111114         ENDIF 
    112115         ! 
     
    139142         IF( MOD( kt-1, kn_fsbc ) == 0 ) THEN         ! correct the freshwater fluxes 
    140143            zcoef = fwfold * rcp 
    141             emp(:,:) = emp(:,:) + fwfold 
    142             qns(:,:) = qns(:,:) - zcoef * sst_m(:,:) ! account for change to the heat budget due to fw correction 
     144            emp(:,:) = emp(:,:) + fwfold             * tmask(:,:,1) 
     145            qns(:,:) = qns(:,:) - zcoef * sst_m(:,:) * tmask(:,:,1) ! account for change to the heat budget due to fw correction 
    143146         ENDIF 
    144147         ! 
  • branches/2015/dev_r5094_UKMO_ISFCLEAN/NEMOGCM/NEMO/OPA_SRC/SBC/sbcssm.F90

    r4990 r5111  
    6161      !!--------------------------------------------------------------------- 
    6262       
    63       !                                        !* first wet T-, U-, V- ocean level (ISF) variables (T, S, depth, velocity) 
     63      !                                        !* surface T-, U-, V- ocean level variables (T, S, depth, velocity) 
    6464      DO jj = 1, jpj 
    6565         DO ji = 1, jpi 
    66             zub(ji,jj)        = ub (ji,jj,miku(ji,jj)) 
    67             zvb(ji,jj)        = vb (ji,jj,mikv(ji,jj)) 
    6866            zts(ji,jj,jp_tem) = tsn(ji,jj,mikt(ji,jj),jp_tem) 
    6967            zts(ji,jj,jp_sal) = tsn(ji,jj,mikt(ji,jj),jp_sal) 
    7068         END DO 
    7169      END DO 
     70      zub(:,:)        = ub (:,:,1       ) 
     71      zvb(:,:)        = vb (:,:,1       ) 
    7272      ! 
    7373      IF( lk_vvl ) THEN 
    74          DO jj = 1, jpj 
    75             DO ji = 1, jpi 
    76                zdep(ji,jj) = fse3t_n(ji,jj,mikt(ji,jj)) 
    77             END DO 
    78          END DO 
     74         zdep(:,:) = fse3t_n(:,:,1) 
    7975      ENDIF 
    8076      !                                                   ! ---------------------------------------- ! 
  • branches/2015/dev_r5094_UKMO_ISFCLEAN/NEMOGCM/NEMO/OPA_SRC/TRA/traadv.F90

    r4990 r5111  
    206206      IF( lk_esopa         )   ioptio =          1 
    207207 
    208       IF( ( ln_traadv_muscl .OR. ln_traadv_muscl2 .OR. ln_traadv_ubs .OR. ln_traadv_qck ) .AND. nn_isf .NE. 0 )  & 
     208      IF( ( ln_traadv_muscl .OR. ln_traadv_muscl2 .OR. ln_traadv_ubs .OR. ln_traadv_qck .OR. ln_traadv_tvd_zts ) .AND. ln_isfcav )  & 
    209209      &   CALL ctl_stop( 'Only traadv_cen2 and traadv_tvd is compatible with ice shelf cavity') 
    210210 
  • branches/2015/dev_r5094_UKMO_ISFCLEAN/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_tvd.F90

    r5098 r5111  
    147147         ! Surface value 
    148148         IF( lk_vvl ) THEN    
    149             DO jj = 1, jpj 
    150                DO ji = 1, jpi 
    151                   zwz(ji,jj, mikt(ji,jj) ) = 0.e0          ! volume variable 
    152                END DO 
    153             END DO 
     149            IF ( ln_isfcav ) THEN 
     150               DO jj = 1, jpj 
     151                  DO ji = 1, jpi 
     152                     zwz(ji,jj, mikt(ji,jj) ) = 0.e0          ! volume variable 
     153                  END DO 
     154               END DO 
     155            ELSE 
     156               zwz(:,:,1) = 0.e0          ! volume variable 
     157            END IF 
    154158         ELSE                 
    155             DO jj = 1, jpj 
    156                DO ji = 1, jpi 
    157                   zwz(ji,jj, mikt(ji,jj) ) = pwn(ji,jj,mikt(ji,jj)) * ptb(ji,jj,mikt(ji,jj),jn)   ! linear free surface  
    158                END DO 
    159             END DO    
     159            IF ( ln_isfcav ) THEN 
     160               DO jj = 1, jpj 
     161                  DO ji = 1, jpi 
     162                     zwz(ji,jj, mikt(ji,jj) ) = pwn(ji,jj,mikt(ji,jj)) * ptb(ji,jj,mikt(ji,jj),jn)   ! linear free surface  
     163                  END DO 
     164               END DO    
     165            ELSE 
     166               zwz(:,:,1) = pwn(:,:,1) * ptb(:,:,1,jn)   ! linear free surface 
     167            END IF 
    160168         ENDIF 
    161169 
     
    212220         END DO 
    213221         ! surface value 
    214          DO jj = 1, jpj 
    215             DO ji = 1, jpi 
    216                zwz(ji,jj,mikt(ji,jj)) = 0.e0 
    217             END DO 
    218          END DO 
     222         IF ( ln_isfcav ) THEN 
     223            DO jj = 1, jpj 
     224               DO ji = 1, jpi 
     225                  zwz(ji,jj,mikt(ji,jj)) = 0.e0 
     226               END DO 
     227            END DO 
     228         ELSE 
     229            zwz(:,:,1) = 0.e0 
     230         END IF 
    219231         CALL lbc_lnk( zwx, 'U', -1. )   ;   CALL lbc_lnk( zwy, 'V', -1. )         ! Lateral bondary conditions 
    220232         CALL lbc_lnk( zwz, 'W',  1. ) 
     
    360372 
    361373         ! upstream tracer flux in the k direction 
    362          ! Surface value 
    363          IF( lk_vvl ) THEN   ;   zwz(:,:, 1 ) = 0._wp                        ! volume variable 
    364          ELSE                ;   zwz(:,:, 1 ) = pwn(:,:,1) * ptb(:,:,1,jn)   ! linear free surface  
    365          ENDIF 
    366374         ! Interior value 
    367375         DO jk = 2, jpkm1 
     
    374382            END DO 
    375383         END DO 
     384         ! Surface value 
     385         IF( lk_vvl ) THEN 
     386            IF ( ln_isfcav ) THEN 
     387               DO jj = 1, jpj 
     388                  DO ji = 1, jpi 
     389                     zwz(ji,jj, mikt(ji,jj) ) = 0.e0          ! volume variable +    isf 
     390                  END DO 
     391               END DO 
     392            ELSE 
     393               zwz(:,:,1) = 0.e0                              ! volume variable + no isf 
     394            END IF 
     395         ELSE 
     396            IF ( ln_isfcav ) THEN 
     397               DO jj = 1, jpj 
     398                  DO ji = 1, jpi 
     399                     zwz(ji,jj, mikt(ji,jj) ) = pwn(ji,jj,mikt(ji,jj)) * ptb(ji,jj,mikt(ji,jj),jn)   ! linear free surface +    isf 
     400                  END DO 
     401               END DO 
     402            ELSE 
     403               zwz(:,:,1) = pwn(:,:,1) * ptb(:,:,1,jn)                                               ! linear free surface + no isf 
     404            END IF 
     405         ENDIF 
    376406 
    377407         ! total advective trend 
     
    583613 
    584614      DO jk = 1, jpkm1 
     615         ikm1 = MAX(jk-1,1) 
     616         z2dtt = p2dt(jk) 
    585617         DO jj = 2, jpjm1 
    586618            DO ji = fs_2, fs_jpim1   ! vector opt. 
    587                ikm1 = MAX(jk-1,1) 
    588                z2dtt = p2dt(jk) 
    589                 
     619 
    590620               ! search maximum in neighbourhood 
    591621               zup = MAX(  zbup(ji  ,jj  ,jk  ),   & 
  • branches/2015/dev_r5094_UKMO_ISFCLEAN/NEMOGCM/NEMO/OPA_SRC/TRA/traldf.F90

    r4990 r5111  
    290290      IF(lwp) WRITE(numout,*) '              homogeneous ocean T = ', zt0, ' S = ',zs0 
    291291 
     292      ! Initialisation of gtui/gtvi in case of no cavity 
     293      IF ( .NOT. ln_isfcav ) THEN 
     294         gtui(:,:,:) = 0.0_wp 
     295         gtvi(:,:,:) = 0.0_wp 
     296      END IF 
    292297      !                                        ! T & S profile (to be coded +namelist parameter 
    293298 
  • branches/2015/dev_r5094_UKMO_ISFCLEAN/NEMOGCM/NEMO/OPA_SRC/TRA/traldf_iso.F90

    r5098 r5111  
    187187            END DO 
    188188         END DO 
    189             ! surface boundary condition: zdkt(jk=1)=zdkt(jk=2) 
     189         ! surface boundary condition: zdkt(jk=1)=zdkt(jk=2) 
    190190         zdk1t(:,:,1) = ( ptb(:,:,1,jn  ) - ptb(:,:,2,jn) ) * wmask(:,:,2) 
    191191         zdkt (:,:,1) = zdk1t(:,:,1) 
     
    200200         END IF 
    201201 
    202             ! 2. Horizontal fluxes 
    203             ! --------------------    
     202         ! 2. Horizontal fluxes 
     203         ! --------------------    
    204204         DO jk = 1, jpkm1 
    205205            DO jj = 1 , jpjm1 
  • branches/2015/dev_r5094_UKMO_ISFCLEAN/NEMOGCM/NEMO/OPA_SRC/TRA/zpshde.F90

    r5098 r5111  
    8282      !!          di(rho) = rd~ - rd(i,j,k)   or   rd(i+1,j,k) - rd~ 
    8383      !! 
    84       !! ** Action  : compute for top and bottom interfaces 
     84      !! ** Action  : compute for top interfaces 
    8585      !!              - pgtu, pgtv: horizontal gradient of tracer at u- & v-points 
    8686      !!              - pgru, pgrv: horizontal gradient of rho (if present) at u- & v-points 
     
    9191      REAL(wp), DIMENSION(jpi,jpj,    kjpt), INTENT(  out)           ::  pgtu, pgtv  ! hor. grad. of ptra at u- & v-pts  
    9292      REAL(wp), DIMENSION(jpi,jpj,jpk     ), INTENT(in   ), OPTIONAL ::  prd         ! 3D density anomaly fields 
    93       REAL(wp), DIMENSION(jpi,jpj         ), INTENT(  out), OPTIONAL ::  pgru, pgrv      ! hor. grad of prd at u- & v-pts (bottom) 
     93      REAL(wp), DIMENSION(jpi,jpj         ), INTENT(  out), OPTIONAL ::  pgru, pgrv  ! hor. grad of prd at u- & v-pts (bottom) 
    9494      ! 
    9595      INTEGER  ::   ji, jj, jn      ! Dummy loop indices 
    96       INTEGER  ::   iku, ikv, ikum1, ikvm1,ikup1, ikvp1   ! partial step level (ocean bottom level) at u- and v-points 
    97       REAL(wp) ::  ze3wu, ze3wv, zmaxu, zmaxv, zdzwu, zdzwv, zdzwuip1, zdzwvjp1  ! temporary scalars 
     96      INTEGER  ::   iku, ikv, ikum1, ikvm1   ! partial step level (ocean bottom level) at u- and v-points 
     97      REAL(wp) ::  ze3wu, ze3wv, zmaxu, zmaxv  ! temporary scalars 
    9898      REAL(wp), DIMENSION(jpi,jpj)      ::  zri, zrj, zhi, zhj   ! NB: 3rd dim=1 to use eos 
    9999      REAL(wp), DIMENSION(jpi,jpj,kjpt) ::  zti, ztj             !  
     
    172172            END DO 
    173173         END DO 
    174           
     174 
    175175         ! Compute interpolated rd from zti, ztj for the 2 cases at the depth of the partial 
    176176         ! step and store it in  zri, zrj for each  case 
     
    193193            END DO 
    194194         END DO 
    195          CALL lbc_lnk( pgru   , 'U', -1. )   ;   CALL lbc_lnk( pgrv  , 'V', -1. )   ! Lateral boundary conditions 
     195         CALL lbc_lnk( pgru , 'U', -1. )   ;   CALL lbc_lnk( pgrv , 'V', -1. )   ! Lateral boundary conditions 
    196196         ! 
    197197      END IF 
  • branches/2015/dev_r5094_UKMO_ISFCLEAN/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfbfr.F90

    r4990 r5111  
    120120                  zbfrt(ji,jj) = MAX(bfrcoef2d(ji,jj), ztmp) 
    121121                  zbfrt(ji,jj) = MIN(zbfrt(ji,jj), rn_bfri2_max) 
    122 ! (ISF) 
    123                   ikbt = mikt(ji,jj) 
    124 ! JC: possible WAD implementation should modify line below if layers vanish 
    125                   ztmp = tmask(ji,jj,ikbt) * ( vkarmn / LOG( 0.5_wp * fse3t_n(ji,jj,ikbt) / rn_bfrz0 ))**2._wp 
    126                   ztfrt(ji,jj) = MAX(tfrcoef2d(ji,jj), ztmp) 
    127                   ztfrt(ji,jj) = MIN(ztfrt(ji,jj), rn_tfri2_max) 
    128  
    129122               END DO 
    130123            END DO 
     124! (ISF) 
     125            IF ( ln_isfcav ) THEN 
     126               DO jj = 1, jpj 
     127                  DO ji = 1, jpi 
     128                     ikbt = mikt(ji,jj) 
     129! JC: possible WAD implementation should modify line below if layers vanish 
     130                     ztmp = (1-tmask(ji,jj,1)) * ( vkarmn / LOG( 0.5_wp * fse3t_n(ji,jj,ikbt) / rn_bfrz0 ))**2._wp 
     131                     ztfrt(ji,jj) = MAX(tfrcoef2d(ji,jj), ztmp) 
     132                     ztfrt(ji,jj) = MIN(ztfrt(ji,jj), rn_tfri2_max) 
     133                  END DO 
     134               END DO 
     135            END IF 
    131136         !    
    132137         ELSE 
     
    152157               ! 
    153158               ! in case of 2 cell water column, we assume each cell feels the top and bottom friction 
    154                IF ( miku(ji,jj) + 2 .GE. mbku(ji,jj) ) THEN 
    155                   bfrua(ji,jj) = - 0.5_wp * ( ( zbfrt(ji,jj) + zbfrt(ji+1,jj  ) )   & 
    156                                &            + ( ztfrt(ji,jj) + ztfrt(ji+1,jj  ) ) ) & 
    157                                &          * zecu * (1._wp - umask(ji,jj,1)) 
    158                END IF 
    159                IF ( mikv(ji,jj) + 2 .GE. mbkv(ji,jj) ) THEN 
    160                   bfrva(ji,jj) = - 0.5_wp * ( ( zbfrt(ji,jj) + zbfrt(ji  ,jj+1) )   & 
    161                                &            + ( ztfrt(ji,jj) + ztfrt(ji  ,jj+1) ) ) & 
    162                                &          * zecv * (1._wp - vmask(ji,jj,1)) 
    163                END IF 
    164                ! (ISF) ======================================================================== 
    165                ikbu = miku(ji,jj)         ! ocean bottom level at u- and v-points  
    166                ikbv = mikv(ji,jj)         ! (deepest ocean u- and v-points) 
    167                ! 
    168                zvu  = 0.25 * (  vn(ji,jj  ,ikbu) + vn(ji+1,jj  ,ikbu)     & 
    169                   &           + vn(ji,jj-1,ikbu) + vn(ji+1,jj-1,ikbu)  ) 
    170                zuv  = 0.25 * (  un(ji,jj  ,ikbv) + un(ji-1,jj  ,ikbv)     & 
    171                   &           + un(ji,jj+1,ikbv) + un(ji-1,jj+1,ikbv)  ) 
    172                ! 
    173                zecu = SQRT(  un(ji,jj,ikbu) * un(ji,jj,ikbu) + zvu*zvu + rn_bfeb2 ) 
    174                zecv = SQRT(  vn(ji,jj,ikbv) * vn(ji,jj,ikbv) + zuv*zuv + rn_bfeb2 ) 
    175                ! 
    176                tfrua(ji,jj) = - 0.5_wp * ( ztfrt(ji,jj) + ztfrt(ji+1,jj  ) ) * zecu * (1._wp - umask(ji,jj,1)) 
    177                tfrva(ji,jj) = - 0.5_wp * ( ztfrt(ji,jj) + ztfrt(ji  ,jj+1) ) * zecv * (1._wp - vmask(ji,jj,1)) 
    178                ! (ISF) END ==================================================================== 
    179                ! in case of 2 cell water column, we assume each cell feels the top and bottom friction 
    180                IF ( miku(ji,jj) + 2 .GE. mbku(ji,jj) ) THEN 
    181                   tfrua(ji,jj) = - 0.5_wp * ( ( ztfrt(ji,jj) + ztfrt(ji+1,jj  ) )   & 
    182                                &            + ( zbfrt(ji,jj) + zbfrt(ji+1,jj  ) ) ) & 
    183                                &          * zecu * (1._wp - umask(ji,jj,1)) 
    184                END IF 
    185                IF ( mikv(ji,jj) + 2 .GE. mbkv(ji,jj) ) THEN 
    186                   tfrva(ji,jj) = - 0.5_wp * ( ( ztfrt(ji,jj) + ztfrt(ji  ,jj+1) )   & 
    187                                &            + ( zbfrt(ji,jj) + zbfrt(ji  ,jj+1) ) ) & 
    188                                &          * zecv * (1._wp - vmask(ji,jj,1)) 
     159               IF ( ln_isfcav ) THEN 
     160                  IF ( miku(ji,jj) + 1 .GE. mbku(ji,jj) ) THEN 
     161                     bfrua(ji,jj) = - 0.5_wp * ( ( zbfrt(ji,jj) + zbfrt(ji+1,jj  ) )   & 
     162                                  &            + ( ztfrt(ji,jj) + ztfrt(ji+1,jj  ) ) ) & 
     163                                  &          * zecu * (1._wp - umask(ji,jj,1)) 
     164                  END IF 
     165                  IF ( mikv(ji,jj) + 1 .GE. mbkv(ji,jj) ) THEN 
     166                     bfrva(ji,jj) = - 0.5_wp * ( ( zbfrt(ji,jj) + zbfrt(ji  ,jj+1) )   & 
     167                                  &            + ( ztfrt(ji,jj) + ztfrt(ji  ,jj+1) ) ) & 
     168                                  &          * zecv * (1._wp - vmask(ji,jj,1)) 
     169                  END IF 
    189170               END IF 
    190171            END DO 
    191172         END DO 
     173         IF ( ln_isfcav ) THEN 
     174            DO jj = 2, jpjm1 
     175               DO ji = 2, jpim1 
     176                  ! (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                  ! 
     180                  zvu  = 0.25 * (  vn(ji,jj  ,ikbu) + vn(ji+1,jj  ,ikbu)     & 
     181                     &           + vn(ji,jj-1,ikbu) + vn(ji+1,jj-1,ikbu)  ) 
     182                  zuv  = 0.25 * (  un(ji,jj  ,ikbv) + un(ji-1,jj  ,ikbv)     & 
     183                     &           + un(ji,jj+1,ikbv) + un(ji-1,jj+1,ikbv)  ) 
     184              ! 
     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              ! 
     188                  tfrua(ji,jj) = - 0.5_wp * ( ztfrt(ji,jj) + ztfrt(ji+1,jj  ) ) * zecu * (1._wp - umask(ji,jj,1)) 
     189                  tfrva(ji,jj) = - 0.5_wp * ( ztfrt(ji,jj) + ztfrt(ji  ,jj+1) ) * zecv * (1._wp - vmask(ji,jj,1)) 
     190              ! (ISF) END ==================================================================== 
     191              ! in case of 2 cell water column, we assume each cell feels the top and bottom friction 
     192                  IF ( miku(ji,jj) + 1 .GE. mbku(ji,jj) ) THEN 
     193                     tfrua(ji,jj) = - 0.5_wp * ( ( ztfrt(ji,jj) + ztfrt(ji+1,jj  ) )   & 
     194                                  &            + ( zbfrt(ji,jj) + zbfrt(ji+1,jj  ) ) ) & 
     195                                  &          * zecu * (1._wp - umask(ji,jj,1)) 
     196                  END IF 
     197                  IF ( mikv(ji,jj) + 1 .GE. mbkv(ji,jj) ) THEN 
     198                     tfrva(ji,jj) = - 0.5_wp * ( ( ztfrt(ji,jj) + ztfrt(ji  ,jj+1) )   & 
     199                                  &            + ( zbfrt(ji,jj) + zbfrt(ji  ,jj+1) ) ) & 
     200                                  &          * zecv * (1._wp - vmask(ji,jj,1)) 
     201                  END IF 
     202               END DO 
     203            END DO 
     204         END IF 
    192205         ! 
    193206         CALL lbc_lnk( bfrua, 'U', 1. )   ;   CALL lbc_lnk( bfrva, 'V', 1. )      ! Lateral boundary condition 
  • branches/2015/dev_r5094_UKMO_ISFCLEAN/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfini.F90

    r4990 r5111  
    1414   !!---------------------------------------------------------------------- 
    1515   USE par_oce         ! mesh and scale factors 
    16    USE sbc_oce         ! surface module (only for nn_isf in the option compatibility test) 
    1716   USE ldftra_oce      ! ocean active tracers: lateral physics 
    1817   USE ldfdyn_oce      ! ocean dynamics lateral physics 
     
    118117      IF( ioptio == 0 .OR. ioptio > 1 .AND. .NOT. lk_esopa )   & 
    119118         &   CALL ctl_stop( ' one and only one vertical diffusion option has to be defined ' ) 
    120       IF( ( lk_zdfric .OR. lk_zdfgls .OR. lk_zdfkpp ) .AND. nn_isf .NE. 0 )   & 
     119      IF( ( lk_zdfric .OR. lk_zdfgls .OR. lk_zdfkpp ) .AND. ln_isfcav )   & 
    121120         &   CALL ctl_stop( ' only zdfcst and zdftke were tested with ice shelves cavities ' ) 
    122121      ! 
  • branches/2015/dev_r5094_UKMO_ISFCLEAN/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftke.F90

    r5104 r5111  
    336336               avmu(ji,jj,jk) = avmu(ji,jj,jk) * ( un(ji,jj,jk-1) - un(ji,jj,jk) )   & 
    337337                  &                            * ( ub(ji,jj,jk-1) - ub(ji,jj,jk) )   &  
    338                   &           / (  fse3uw_n(ji,jj,jk)         & 
    339                   &              * fse3uw_b(ji,jj,jk) ) 
     338                  &                            / (  fse3uw_n(ji,jj,jk)               & 
     339                  &                              *  fse3uw_b(ji,jj,jk) ) 
    340340               avmv(ji,jj,jk) = avmv(ji,jj,jk) * ( vn(ji,jj,jk-1) - vn(ji,jj,jk) )   & 
    341341                  &                            * ( vb(ji,jj,jk-1) - vb(ji,jj,jk) )   & 
Note: See TracChangeset for help on using the changeset viewer.