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 4724 for branches/2014 – NEMO

Changeset 4724 for branches/2014


Ignore:
Timestamp:
2014-07-18T17:32:27+02:00 (10 years ago)
Author:
mathiot
Message:

ISF branch: add comments, fix mpp and restar issues, add test to stop if incompatible options and fix mask issue in sbcice and sbcblk.

Location:
branches/2014/dev_r4650_UKMO2_ice_shelves/NEMOGCM/NEMO/OPA_SRC
Files:
21 edited

Legend:

Unmodified
Added
Removed

  • branches/2014/dev_r4650_UKMO2_ice_shelves/NEMOGCM/NEMO/OPA_SRC/DIA/diaar5.F90

    r4313 r4724  
    104104         zbotpres(:,:) = zbotpres(:,:) + fse3t(:,:,jk) * zrhd(:,:,jk) 
    105105      END DO 
    106       IF( .NOT.lk_vvl )   zbotpres(:,:) = zbotpres(:,:) + sshn(:,:) * zrhd(:,:,1) 
     106      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 
     112      END IF 
    107113      !                                          
    108114      zarho = SUM( area(:,:) * zbotpres(:,:) )  
    … …  
    120126         zbotpres(:,:) = zbotpres(:,:) + fse3t(:,:,jk) * zrhd(:,:,jk) 
    121127      END DO 
    122       IF( .NOT.lk_vvl )   zbotpres(:,:) = zbotpres(:,:) + sshn(:,:) * zrhd(:,:,1) 
     128      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 
     134      END IF 
    123135      !     
    124136      zarho = SUM( area(:,:) * zbotpres(:,:) )  
    … …  
    145157      END DO 
    146158      IF( .NOT.lk_vvl ) THEN 
    147          ztemp = ztemp + SUM( zarea_ssh(:,:) * tsn(:,:,1,jp_tem) ) 
    148          zsal  = zsal  + SUM( zarea_ssh(:,:) * tsn(:,:,1,jp_sal) ) 
     159         DO ji=1,jpi 
     160            DO jj=1,jpj 
     161               ztemp = ztemp + SUM( zarea_ssh(ji,jj) * tsn(ji,jj,mikt(ji,jj),jp_tem) ) 
     162               zsal  = zsal  + SUM( zarea_ssh(ji,jj) * tsn(ji,jj,mikt(ji,jj),jp_sal) ) 
     163            END DO 
     164         END DO 
    149165      ENDIF 
    150166      IF( lk_mpp ) THEN   

  • branches/2014/dev_r4650_UKMO2_ice_shelves/NEMOGCM/NEMO/OPA_SRC/DIA/diaharm.F90

    r4624 r4724  
    196196              DO ji = 1,jpi 
    197197                ! Elevation 
    198                 ana_temp(ji,jj,nhc,1) = ana_temp(ji,jj,nhc,1) + ztemp*sshn(ji,jj)           *tmask(ji,jj,1)         
     198                ana_temp(ji,jj,nhc,1) = ana_temp(ji,jj,nhc,1) + ztemp*sshn(ji,jj)           *tmask_i(ji,jj)         
    199199#if defined key_dynspg_ts 
    200                 ana_temp(ji,jj,nhc,2) = ana_temp(ji,jj,nhc,2) + ztemp*un_b(ji,jj)*hur(ji,jj)*umask(ji,jj,1) 
    201                 ana_temp(ji,jj,nhc,3) = ana_temp(ji,jj,nhc,3) + ztemp*vn_b(ji,jj)*hvr(ji,jj)*vmask(ji,jj,1) 
     200                ana_temp(ji,jj,nhc,2) = ana_temp(ji,jj,nhc,2) + ztemp*un_b(ji,jj)*hur(ji,jj)*umask_i(ji,jj) 
     201                ana_temp(ji,jj,nhc,3) = ana_temp(ji,jj,nhc,3) + ztemp*vn_b(ji,jj)*hvr(ji,jj)*vmask_i(ji,jj) 
    202202#endif 
    203203              END DO 
    … …  
    294294               X1 = ana_amp(ji,jj,jh,1) 
    295295               X2 =-ana_amp(ji,jj,jh,2) 
    296                out_eta(ji,jj,jh       ) = X1 * tmask(ji,jj,1) 
    297                out_eta(ji,jj,jh+nb_ana) = X2 * tmask(ji,jj,1) 
     296               out_eta(ji,jj,jh       ) = X1 * tmask_i(ji,jj) 
     297               out_eta(ji,jj,jh+nb_ana) = X2 * tmask_i(ji,jj) 
    298298            ENDDO 
    299299         ENDDO 
    … …  
    328328               X1=ana_amp(ji,jj,jh,1) 
    329329               X2=-ana_amp(ji,jj,jh,2) 
    330                out_u(ji,jj,jh) = X1 * umask(ji,jj,1) 
    331                out_u (ji,jj,nb_ana+jh) = X2 * umask(ji,jj,1) 
     330               out_u(ji,jj,jh) = X1 * umask_i(ji,jj) 
     331               out_u (ji,jj,nb_ana+jh) = X2 * umask_i(ji,jj) 
    332332            ENDDO 
    333333         ENDDO 
    … …  
    362362               X1=ana_amp(ji,jj,jh,1) 
    363363               X2=-ana_amp(ji,jj,jh,2) 
    364                out_v(ji,jj,jh)=X1 * vmask(ji,jj,1) 
    365                out_v(ji,jj,nb_ana+jh)=X2 * vmask(ji,jj,1) 
     364               out_v(ji,jj,jh)=X1 * vmask_i(ji,jj) 
     365               out_v(ji,jj,nb_ana+jh)=X2 * vmask_i(ji,jj) 
    366366            ENDDO 
    367367         ENDDO 

  • branches/2014/dev_r4650_UKMO2_ice_shelves/NEMOGCM/NEMO/OPA_SRC/DIA/diahsb.F90

    r4704 r4724  
    6969      INTEGER, INTENT(in) ::   kt   ! ocean time-step index 
    7070      !! 
    71       INTEGER    ::   jk                          ! dummy loop indice 
     71      INTEGER    ::   jk, ji, jj                          ! dummy loop indice 
    7272      REAL(dp)   ::   zdiff_hc    , zdiff_sc      ! heat and salt content variations 
    7373      REAL(dp)   ::   zdiff_hc1   , zdiff_sc1     ! -   -   -   -   -   -   -   -  
    … …  
    7979      REAL(dp)   ::   z_wn_trd_t , z_wn_trd_s   !    -     - 
    8080      REAL(dp)   ::   z_ssh_hc , z_ssh_sc   !    -     - 
     81      REAL(dp), DIMENSION(:,:), POINTER       ::   z2d0, z2d1 
    8182      !!--------------------------------------------------------------------------- 
    8283      IF( nn_timing == 1 )   CALL timing_start('dia_hsb')       
    83  
     84      CALL wrk_alloc( jpi,jpj, z2d0, z2d1 ) 
     85      tsn(:,:,:,1) = tsn(:,:,:,1) * tmask(:,:,:) ; tsb(:,:,:,1) = tsb(:,:,:,1) * tmask(:,:,:) ; 
     86      tsn(:,:,:,2) = tsn(:,:,:,2) * tmask(:,:,:) ; tsb(:,:,:,2) = tsb(:,:,:,2) * tmask(:,:,:) ; 
    8487      ! ------------------------- ! 
    8588      ! 1 - Trends due to forcing ! 
    8689      ! ------------------------- ! 
    87       z_frc_trd_v = r1_rau0 * glob_sum( - ( emp(:,:) - rnf(:,:) + rdivisf * fwfisf(:,:)) * surf(:,:) ) ! volume fluxes 
    88       z_frc_trd_t =           glob_sum( sbc_tsc(:,:,jp_tem) * surf(:,:) )       ! heat fluxes 
    89       z_frc_trd_s =           glob_sum( sbc_tsc(:,:,jp_sal) * surf(:,:) )       ! salt fluxes 
     90      z_frc_trd_v = r1_rau0 * glob_sum( - ( emp(:,:) - rnf(:,:) + rdivisf * fwfisf(:,:) ) * surf(:,:) ) ! volume fluxes 
     91      z_frc_trd_t =           glob_sum( sbc_tsc(:,:,jp_tem) * surf(:,:) )                               ! heat fluxes 
     92      z_frc_trd_s =           glob_sum( sbc_tsc(:,:,jp_sal) * surf(:,:) )                               ! salt fluxes 
    9093      ! Add runoff heat & salt input 
    9194      IF( ln_rnf    )   z_frc_trd_t = z_frc_trd_t + glob_sum( rnf_tsc(:,:,jp_tem) * surf(:,:) ) 
    9295      IF( ln_rnf_sal)   z_frc_trd_s = z_frc_trd_s + glob_sum( rnf_tsc(:,:,jp_sal) * surf(:,:) ) 
     96      ! Add geothermal ice shelf 
     97      IF( nn_isf .GE. 1 )  THEN 
     98                        z_frc_trd_t = z_frc_trd_t +                      glob_sum( ( risf_tsc(:,:,jp_tem) - rdivisf * fwfisf(:,:) * -1.9 * r1_rau0 ) * surf(:,:) ) 
     99                        z_frc_trd_s = z_frc_trd_s + (1.0_wp - rdivisf) * glob_sum( risf_tsc(:,:,jp_sal) * surf(:,:) ) 
     100      ENDIF 
    93101 
    94102      ! Add penetrative solar radiation 
    … …  
    98106      ! 
    99107      IF( .NOT. lk_vvl ) THEN 
    100          z_wn_trd_t = - glob_sum( surf(:,:) * wn(:,:,1) * tsb(:,:,1,jp_tem) ) 
    101          z_wn_trd_s = - glob_sum( surf(:,:) * wn(:,:,1) * tsb(:,:,1,jp_sal) ) 
     108         z2d0=0.0_wp ; z2d1=0.0_wp 
     109         DO ji=1,jpi 
     110            DO jj=1,jpj 
     111              z2d0(ji,jj) = surf(ji,jj) * wn(ji,jj,mikt(ji,jj)) * tsb(ji,jj,mikt(ji,jj),jp_tem) 
     112              z2d1(ji,jj) = surf(ji,jj) * wn(ji,jj,mikt(ji,jj)) * tsb(ji,jj,mikt(ji,jj),jp_sal) 
     113            ENDDO 
     114         ENDDO 
     115         z_wn_trd_t = - glob_sum( z2d0 )  
     116         z_wn_trd_s = - glob_sum( z2d1 ) 
    102117      ENDIF 
    103118 
    … …  
    123138      ! heat & salt content variation (associated with ssh) 
    124139      IF( .NOT. lk_vvl ) THEN 
    125          z_ssh_hc = glob_sum( surf(:,:) * ( tsn(:,:,1,jp_tem) * sshn(:,:) - ssh_hc_loc_ini(:,:) ) ) 
    126          z_ssh_sc = glob_sum( surf(:,:) * ( tsn(:,:,1,jp_sal) * sshn(:,:) - ssh_sc_loc_ini(:,:) ) ) 
     140         z2d0=0.0_wp ; z2d1=0.0_wp 
     141         DO ji=1,jpi 
     142            DO jj=1,jpj 
     143              z2d0(ji,jj) = surf(ji,jj) * ( tsn(ji,jj,mikt(ji,jj),jp_tem) * sshn(ji,jj) - ssh_hc_loc_ini(ji,jj) )  
     144              z2d1(ji,jj) = surf(ji,jj) * ( tsn(ji,jj,mikt(ji,jj),jp_sal) * sshn(ji,jj) - ssh_sc_loc_ini(ji,jj) )  
     145            ENDDO 
     146         ENDDO 
     147         z_ssh_hc = glob_sum( z2d0 )  
     148         z_ssh_sc = glob_sum( z2d1 )  
    127149      ENDIF 
    128150 
    … …  
    282304     CHARACTER(len=*), INTENT(in) ::   cdrw   ! "READ"/"WRITE" flag 
    283305     ! 
    284      INTEGER ::   jk   !  
     306     INTEGER ::   jk,ji,jj   !  
    285307     INTEGER ::   id1   ! local integers 
    286308     !!---------------------------------------------------------------------- 
    … …  
    322344          frc_s = 0.d0                                           ! salt content   -    -   -    -         
    323345          IF( .NOT. lk_vvl ) THEN 
    324              ssh_hc_loc_ini(:,:) = tsn(:,:,1,jp_tem) * sshn(:,:)   ! initial heat content in ssh 
    325              ssh_sc_loc_ini(:,:) = tsn(:,:,1,jp_sal) * sshn(:,:)   ! initial salt content in ssh 
     346             DO ji=1,jpi 
     347                DO jj=1,jpj 
     348                   ssh_hc_loc_ini(ji,jj) = tsn(ji,jj,mikt(ji,jj),jp_tem) * sshn(ji,jj)   ! initial heat content in ssh 
     349                   ssh_sc_loc_ini(ji,jj) = tsn(ji,jj,mikt(ji,jj),jp_sal) * sshn(ji,jj)   ! initial salt content in ssh 
     350                ENDDO 
     351             ENDDO 
    326352             frc_wn_t = 0.d0                                       ! initial heat content misfit due to free surface 
    327353             frc_wn_s = 0.d0                                       ! initial salt content misfit due to free surface 

  • branches/2014/dev_r4650_UKMO2_ice_shelves/NEMOGCM/NEMO/OPA_SRC/DOM/domvvl.F90

    r4701 r4724  
    312312      !                                                ! --------------------------------------------- ! 
    313313 
    314       z_scale(:,:) = ( ssha(:,:) - sshb(:,:) ) * tmask_i(:,:) / ( ht_0(:,:) + sshn(:,:) + 1. - tmask_i(:,:) ) 
     314      z_scale(:,:) = ( ssha(:,:) - sshb(:,:) ) * lmask(:,:) / ( ht_0(:,:) + sshn(:,:) + 1. - lmask(:,:) ) 
    315315      DO jk = 1, jpkm1 
    316316         ! formally this is the same as fse3t_a = e3t_0*(1+ssha/ht_0) 
    … …  
    472472            zht(:,:)  = zht(:,:) + tilde_e3t_a(:,:,jk) * tmask(:,:,jk) 
    473473         END DO 
    474          z_scale(:,:) =  - zht(:,:) / ( ht_0(:,:) + sshn(:,:) + 1. - tmask(:,:,1) ) 
     474         z_scale(:,:) =  - zht(:,:) / ( ht_0(:,:) + sshn(:,:) + 1. - lmask ) 
    475475         DO jk = 1, jpkm1 
    476476            dtilde_e3t_a(:,:,jk) = dtilde_e3t_a(:,:,jk) + fse3t_n(:,:,jk) * z_scale(:,:) * tmask(:,:,jk) 
    … …  
    998998 
    999999      IF( ioptio /= 1 )   CALL ctl_stop( 'Choose ONE vertical coordinate in namelist nam_vvl' ) 
     1000      IF( .NOT. ln_vvl_zstar .AND. nn_isf .NE. 0)   CALL ctl_stop( 'vvl_ztilde, vvl_layer, vvl_ztilde_as_zstar, vvl_zstar_at_eqtor not tested with ice shelf cavity (only vvl_zstar was tested)' ) 
    10001001 
    10011002      IF(lwp) THEN                   ! Print the choice 

  • branches/2014/dev_r4650_UKMO2_ice_shelves/NEMOGCM/NEMO/OPA_SRC/DOM/domzgr.F90

    r4666 r4724  
    299299      ENDIF 
    300300 
    301 ! need to be like this to compute the pressure gradient with ISF 
     301! need to be like this to compute the pressure gradient with ISF. If not, level beneath the ISF are not aligned (sum(e3t) /= depth) 
    302302! define e3t_0 and e3w_0 as the differences between gdept and gdepw respectively 
    303303      DO jk = 1, jpkm1 
    … …  
    10231023      END WHERE 
    10241024     
    1025  
     1025! basic check for the compatibility of bathy and icedep. I think it should be offline because it is not perfect and cannot solved all the situation 
    10261026      icompt = 0  
     1027! run the bathy check 10 times to be sure all the modif in the bathy or iceshelf draft are compatible together 
    10271028      DO jl = 1, 10  
    10281029        IF( lk_mpp ) THEN 
    … …  
    10581059            bathy(:,:)   = 0._wp 
    10591060         END WHERE 
    1060  
     1061! Case where bathy and icedep compatible but not the level variable mbathy/micedep because of partial cell condition 
    10611062         DO jj = 1, jpj 
    10621063            DO ji = 1, jpi 
    … …  
    10751076         DO jj = 2, jpjm1 
    10761077            DO ji = 2, jpim1 
     1078               ! T point 
    10771079               IF( zmicedep(ji,jj) == zmbathy(ji,jj) .AND. zmbathy(ji,jj) .GT. 1) THEN 
    10781080                  mbathy(ji,jj) = zmbathy(ji,jj) + 1  
    10791081                  bathy(ji,jj)=gdepw_1d(mbathy(ji,jj)) + MIN( e3zps_min, e3t_1d(mbathy(ji,jj))*e3zps_rat ) 
    10801082               ENDIF 
     1083               ! V point 
    10811084               IF( zmicedep(ji,jj+1) == zmbathy(ji,jj) .AND. zmbathy(ji,jj) .GT. 1) THEN 
    10821085                  mbathy(ji,jj) = zmbathy(ji,jj) + 1  
    10831086                  bathy(ji,jj)=gdepw_1d(mbathy(ji,jj)) + MIN( e3zps_min, e3t_1d(mbathy(ji,jj))*e3zps_rat ) 
    10841087               ENDIF 
     1088               ! V point -1 
    10851089               IF( zmicedep(ji,jj-1) == zmbathy(ji,jj) .AND. zmbathy(ji,jj) .GT. 1) THEN 
    10861090                  mbathy(ji,jj) = zmbathy(ji,jj) + 1  
    10871091                  bathy(ji,jj)=gdepw_1d(mbathy(ji,jj)) + MIN( e3zps_min, e3t_1d(mbathy(ji,jj))*e3zps_rat ) 
    10881092               ENDIF 
     1093               ! U point 
    10891094               IF( zmicedep(ji+1,jj) == zmbathy(ji,jj) .AND. zmbathy(ji,jj) .GT. 1) THEN 
    10901095                  mbathy(ji,jj) = zmbathy(ji,jj) + 1  
    10911096                  bathy(ji,jj)=gdepw_1d(mbathy(ji,jj)) + MIN( e3zps_min, e3t_1d(mbathy(ji,jj))*e3zps_rat ) 
    10921097               ENDIF 
     1098               ! U point -1 
    10931099               IF( zmicedep(ji-1,jj) == zmbathy(ji,jj) .AND. zmbathy(ji,jj) .GT. 1) THEN 
    10941100                  mbathy(ji,jj) = zmbathy(ji,jj) + 1  
    … …  
    12171223            mbathy(:,:) = INT( zbathy(:,:) ) 
    12181224        ENDIF  
    1219         ! remove pool of water stuck between ice shelf and bathymetry 
     1225        ! remove 1 cell pool of water stuck between ice shelf and bathymetry (need a 3D flood filling tools to do this properly) 
    12201226        DO jk = 1, jpk 
    12211227        WHERE (micedep==0) micedep=jpk 
    … …  
    13441350                e3w_0  (ji,jj,ik+1) = gdept_1d(ik+1) - gdept_0(ji,jj,ik) 
    13451351 
    1346                 IF( ik + 1 == mbathy(ji,jj) ) THEN               ! ice shelf point only  
     1352                IF( ik + 1 == mbathy(ji,jj) ) THEN               ! ice shelf point only (2 cell water column)  
    13471353                   e3w_0  (ji,jj,ik+1) = gdept_0(ji,jj,ik+1) - gdept_0(ji,jj,ik)  
    13481354                ENDIF  
    … …  
    13981404         DO jj = 1, jpjm1  
    13991405            DO ji = 1, fs_jpim1   ! vector opt.  
    1400 ! (ISF)** NEEDS CHANGING TO SECOND OPTION FOR ICE SHELF BUT WILL CHANGE RESULTS WITHOUT ICE (DIFFER AT THE 1e-13 LEVEL)  
    14011406               e3uw_0(ji,jj,jk) = MIN( gdept_0(ji,jj,jk), gdept_0(ji+1,jj  ,jk) ) - MAX( gdept_0(ji,jj,jk-1), gdept_0(ji+1,jj  ,jk-1) ) 
    14021407               e3vw_0(ji,jj,jk) = MIN( gdept_0(ji,jj,jk), gdept_0(ji  ,jj+1,jk) ) - MAX( gdept_0(ji,jj,jk-1), gdept_0(ji  ,jj+1,jk-1) ) 

  • branches/2014/dev_r4650_UKMO2_ice_shelves/NEMOGCM/NEMO/OPA_SRC/DYN/dynhpg.F90

    r4666 r4724  
    2929   !!---------------------------------------------------------------------- 
    3030   USE oce             ! ocean dynamics and tracers 
     31   USE sbc_oce         ! surface variable (only for the flag with ice shelf) 
    3132   USE dom_oce         ! ocean space and time domain 
    3233   USE phycst          ! physical constants 
    … …  
    176177      IF( ln_hpg_prj )   ioptio = ioptio + 1 
    177178      IF( ioptio /= 1 )   CALL ctl_stop( 'NO or several hydrostatic pressure gradient options used' ) 
     179      IF( (ln_hpg_zco .OR. ln_hpg_zps .OR. ln_hpg_djc .OR. ln_hpg_prj ) .AND. nn_isf .NE. 0 )   CALL ctl_stop( 'Only hpg_sco has been corrected to work with ice shelf cavity. Comparison in a GYRE simulation with bump in the middle show similar result than hpg_zps' ) 
    178180      ! 
    179181   END SUBROUTINE dyn_hpg_init 
    … …  
    395397      ! iniitialised to 0. zhpi zhpi  
    396398      zhpi(:,:,:)=0._wp ; zhpj(:,:,:)=0._wp 
    397  
     399       
     400      ! assume density of water displaced by the ice shelf is at T=-1.9 and S=34.4 (rude) 
    398401      ztstop(:,:,1)=-1.9_wp ; ztstop(:,:,2)=34.4_wp 
    399       zrhd = rhd 
     402      zrhd = rhd ! save rhd 
    400403      DO jk = 1, jpk 
    401404           zdept(:,:)=gdept_1d(jk) 
    … …  
    403406      END DO 
    404407      WHERE ( tmask(:,:,:) == 1._wp) 
    405         rhd(:,:,:) = zrhd(:,:,:) 
     408        rhd(:,:,:) = zrhd(:,:,:) ! replace wet cell by the saved rhd 
    406409      END WHERE 
    407410       
    408  
     411      ! compute rhd at the ice/oce interface (ice shelf side) 
    409412      CALL eos(ztstop,icedep,zrhdtop_isf) 
    410413 
     414      ! compute rhd at the ice/oce interface (ocean side) 
    411415      DO ji=1,jpi 
    412416        DO jj=1,jpj 
    … …  
    419423      ! 
    420424      ! Surface value + ice shelf gradient 
    421       ! compute pressure (used to compute hpgi/j for all the level from 1 to miku/v) 
     425      ! compute pressure due to ice shelf load (used to compute hpgi/j for all the level from 1 to miku/v) 
    422426      ziceload = 0._wp 
    423427      DO jj = 1, jpj 
    … …  
    433437         END DO 
    434438      END DO 
    435       ! compute zp from first level to first wet cell (blue and purple area) 
     439      riceload(:,:) = 0.0_wp ; riceload(:,:)=ziceload(:,:)  ! need to be saved for diaar5 
     440      ! compute zp from z=0 to first T wet point (correction due to zps not yet applied) 
    436441      DO jj = 2, jpjm1 
    437442         DO ji = fs_2, fs_jpim1   ! vector opt. 
    438443            ikt=mikt(ji,jj) ; iktp1i=mikt(ji+1,jj); iktp1j=mikt(ji,jj+1) 
    439444            ! hydrostatic pressure gradient along s-surfaces and ice shelf pressure 
    440             ! we assume ISF is in isostatic equilibrium with a density equal to the reference density 
     445            ! we assume ISF is in isostatic equilibrium 
    441446            zhpi(ji,jj,1) = zcoef0 / e1u(ji,jj) * ( 0.5_wp * fse3w(ji+1,jj  ,iktp1i)                                    & 
    442447               &                                  * ( 2._wp * znad + rhd(ji+1,jj  ,iktp1i) + zrhdtop_oce(ji+1,jj  ) )   & 
    … …  
    449454               &                                  * ( 2._wp * znad + rhd(ji  ,jj  ,ikt   ) + zrhdtop_oce(ji  ,jj  ) )   & 
    450455               &                                  + ( ziceload(ji,jj+1) - ziceload(ji,jj) )                             )  
    451             ! s-coordinate pressure gradient correction 
     456            ! s-coordinate pressure gradient correction (=0 if z coordinate) 
    452457            zuap = -zcoef0 * ( rhd   (ji+1,jj,1) + rhd   (ji,jj,1) + 2._wp * znad )   & 
    453458               &           * ( fsde3w(ji+1,jj,1) - fsde3w(ji,jj,1) ) / e1u(ji,jj) 
    … …  
    471476                      &                                 * ( rhd    (ji+1,jj,iku) + rhd   (ji,jj,iku)       & 
    472477                      &                                   + SIGN(1._wp,ze3wu) * grui(ji,jj) + 2._wp * znad ) 
    473   
    474                IF (mbku(ji,jj) == iku + 1) zpshpi(ji,jj)  = zhpi(ji,jj,iku) 
     478               ! corrective term ( = 0 if z coordinate ) 
    475479               zuap              = -zcoef0 * ( arui(ji,jj) + 2._wp * znad ) * gzui(ji,jj) / e1u(ji,jj) 
    476480               ! zhpi will be added in interior loop and zuapint will be removed in the interior loop 
    477481               ua(ji,jj,iku)     = ua(ji,jj,iku) + zuap 
     482               ! in case of 2 cell water column, need to save the pressure gradient to compute the bottom pressure   
     483               IF (mbku(ji,jj) == iku + 1) zpshpi(ji,jj)  = zhpi(ji,jj,iku) 
    478484 
    479485               ! zhpi will be added in interior loop and zuapint will be removed in the interior loop 
    … …  
    488494            ! v direction 
    489495            ikv = mikv(ji,jj) 
    490             !ze3wv  = 0.5 * (e3w_0(ikv+1) / e3t_0(ikv)) * ( fse3t(ji,jj+1,ikv)-fse3t(ji,jj,ikv) ) 
    491496            ze3wv  = fse3w(ji,jj+1,ikv+1)-fse3w(ji,jj,ikv+1)  
    492497            IF ( ikv .GT. 1 ) THEN 
    493498               ! case ikv 
    494          !      ze3wv             =  - (fsde3w(ji,jj+1,ikv) - fsde3w(ji,jj,ikv) ) 
    495499               zhpj(ji,jj,ikv)   =  zcoef0 / e2v(ji,jj) * ze3wv                                            & 
    496500                     &                                  * ( rhd(ji,jj+1,ikv) + rhd   (ji,jj,ikv)           & 
    497501                     &                                    + SIGN(1._wp,ze3wv) * grvi(ji,jj) + 2._wp * znad ) 
    498  
     502               ! corrective term ( = 0 if z coordinate ) 
    499503               zvap              = -zcoef0 * ( arvi(ji,jj) + 2._wp * znad ) * gzvi(ji,jj) / e2v(ji,jj) 
    500504               ! zhpi will be added in interior loop and zvapint will be removed in the interior loop 
    501505               va(ji,jj,ikv)      = va(ji,jj,ikv) + zvap 
     506               ! in case of 2 cell water column, need to save the pressure gradient to compute the bottom pressure   
    502507               IF (mbkv(ji,jj) == ikv + 1)  zpshpj(ji,jj)  =  zhpj(ji,jj,ikv)  
    503508               ! zhpi will be added in interior loop and zvapint will be removed in the interior loop 
    … …  
    526531                  &                                       - fse3w(ji  ,jj  ,jk) * ( (rhd(ji  ,jj,jk  ) + znad)                 & 
    527532                  &                                                     + (rhd(ji  ,jj,jk-1) + znad) ) * tmask(ji  ,jj,jk-1)   )  
    528                ! corrective term               
     533             !corrective term, we mask this term for the first wet level beneath the ice shelf (contribution done in the upper bloc)  
    529534               zuap = - zcoef0 * ( rhd   (ji+1,jj  ,jk) + rhd   (ji,jj,jk) + 2._wp * znad )   & 
    530535                  &            * ( fsde3w(ji+1,jj  ,jk) - fsde3w(ji,jj,jk) ) / e1u(ji,jj) * umask(ji,jj,jk-1) 
    … …  
    537542                  &                                       - fse3w(ji  ,jj  ,jk) * ( (rhd(ji,jj  ,jk  ) + znad)                 & 
    538543                  &                                                     + (rhd(ji,jj  ,jk-1) + znad) ) * tmask(ji,jj  ,jk-1)   ) 
    539                ! 
     544             !corrective term, we mask this term for the first wet level beneath the ice shelf (contribution done in the upper bloc) 
    540545               zvap = - zcoef0 * ( rhd   (ji  ,jj+1,jk) + rhd   (ji,jj,jk) + 2._wp * znad )   & 
    541546                  &            * ( fsde3w(ji  ,jj+1,jk) - fsde3w(ji,jj,jk) ) / e2v(ji,jj) * vmask(ji,jj,jk-1) 
    … …  
    563568               ua(ji,jj,iku)   = ua(ji,jj,iku) - zhpi(ji,jj,iku) - zuap 
    564569               ! put new value 
     570               ! -zpshpi to avoid double contribution of the partial step in the top layer  
    565571               zuap            = -zcoef0 * ( aru(ji,jj) + 2._wp * znad ) * gzu(ji,jj)  / e1u(ji,jj) 
    566572               zhpi(ji,jj,iku) =  zhpi(ji,jj,iku-1) + zcoef0 / e1u(ji,jj) * ge3ru(ji,jj) - zpshpi(ji,jj)  
    … …  
    574580               va(ji,jj,ikv)   = va(ji,jj,ikv) - zhpj(ji,jj,ikv) - zvap 
    575581               ! put new value 
     582               ! -zpshpj to avoid double contribution of the partial step in the top layer  
    576583               zvap            = -zcoef0 * ( arv(ji,jj) + 2._wp * znad ) * gzv(ji,jj)     / e2v(ji,jj) 
    577584               zhpj(ji,jj,ikv) =  zhpj(ji,jj,ikv-1) + zcoef0 / e2v(ji,jj) * ge3rv(ji,jj) - zpshpj(ji,jj)    

  • branches/2014/dev_r4650_UKMO2_ice_shelves/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg.F90

    r4496 r4724  
    246246      IF( ( ioptio > 1 .AND. .NOT. lk_esopa ) .OR. ( ioptio == 0 .AND. .NOT. lk_c1d ) )   & 
    247247           &   CALL ctl_stop( ' Choose only one surface pressure gradient scheme with a key cpp' ) 
     248      IF( ( lk_dynspg_ts .OR. lk_dynspg_exp ) .AND. nn_isf .NE. 0 )   & 
     249           &   CALL ctl_stop( ' dynspg_ts and dynspg_exp not tested with ice shelf cavity ' ) 
    248250      ! 
    249251      IF( lk_esopa     )   nspg = -1 

  • branches/2014/dev_r4650_UKMO2_ice_shelves/NEMOGCM/NEMO/OPA_SRC/DYN/sshwzv.F90

    r4486 r4724  
    111111      !  
    112112      z1_rau0 = 0.5_wp * r1_rau0 
    113       ssha(:,:) = (  sshb(:,:) - z2dt * ( z1_rau0 * ( emp_b(:,:) + emp(:,:) ) + zhdiv(:,:) )  ) * tmask(:,:,1) 
     113      ssha(:,:) = (  sshb(:,:) - z2dt * ( z1_rau0 * ( emp_b(:,:) + emp(:,:) ) + zhdiv(:,:) )  ) * lmask(:,:) 
    114114 
    115115#if ! defined key_dynspg_ts 
    … …  
    291291      ELSE                                         !** Leap-Frog time-stepping: Asselin filter + swap 
    292292         sshb(:,:) = sshn(:,:) + atfp * ( sshb(:,:) - 2 * sshn(:,:) + ssha(:,:) )     ! before <-- now filtered 
    293          IF( lk_vvl ) sshb(:,:) = sshb(:,:) - atfp * rdt / rau0 * ( emp_b(:,:) - emp(:,:) ) * tmask(:,:,1) 
     293         IF( lk_vvl ) sshb(:,:) = sshb(:,:) - atfp * rdt / rau0 * ( emp_b(:,:) - emp(:,:) ) * lmask(:,:) 
    294294         sshn(:,:) = ssha(:,:)                           ! now <-- after 
    295295      ENDIF 

  • branches/2014/dev_r4650_UKMO2_ice_shelves/NEMOGCM/NEMO/OPA_SRC/ICB/icbutl.F90

    r3821 r4724  
    7070      ! and ssh which is used to calculate gradients 
    7171 
    72       uo_e(:,:) = 0._wp ;   uo_e(1:jpi, 1:jpj) = ssu_m(:,:) 
    73       vo_e(:,:) = 0._wp ;   vo_e(1:jpi, 1:jpj) = ssv_m(:,:) 
     72      uo_e(:,:) = 0._wp ;   uo_e(1:jpi, 1:jpj) = ssu_m(:,:) * umask(:,:,1) 
     73      vo_e(:,:) = 0._wp ;   vo_e(1:jpi, 1:jpj) = ssv_m(:,:) * vmask(:,:,1) 
    7474      ff_e(:,:) = 0._wp ;   ff_e(1:jpi, 1:jpj) = ff   (:,:) 
    75       ua_e(:,:) = 0._wp ;   ua_e(1:jpi, 1:jpj) = utau (:,:) 
    76       va_e(:,:) = 0._wp ;   va_e(1:jpi, 1:jpj) = vtau (:,:) 
     75      ua_e(:,:) = 0._wp ;   ua_e(1:jpi, 1:jpj) = utau (:,:) * umask(:,:,1) ! maybe mask useless because mask applied in sbcblk 
     76      va_e(:,:) = 0._wp ;   va_e(1:jpi, 1:jpj) = vtau (:,:) * vmask(:,:,1) ! maybe mask useless because mask applied in sbcblk 
    7777 
    7878      CALL lbc_lnk_e( uo_e, 'U', -1._wp, 1, 1 ) 
    … …  
    9393      !! so fudge some numbers all the way around the boundary 
    9494 
    95       ssh_e(:,:) = 0._wp ;   ssh_e(1:jpi, 1:jpj) = ssh_m(:,:) 
     95      ssh_e(:,:) = 0._wp ;   ssh_e(1:jpi, 1:jpj) = ssh_m(:,:) * tmask(:,:,1) 
    9696      ssh_e(0    ,    :) = ssh_e(1  ,  :) 
    9797      ssh_e(jpi+1,    :) = ssh_e(jpi,  :) 

  • branches/2014/dev_r4650_UKMO2_ice_shelves/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_clio.F90

    r4624 r4724  
    257257         END DO 
    258258      END DO 
     259      utau(:,:) = utau(:,:) * umask(:,:,1) 
     260      vtau(:,:) = vtau(:,:) * vmask(:,:,1) 
     261      taum(:,:) = taum(:,:) * tmask(:,:,1) 
    259262      CALL lbc_lnk( taum, 'T', 1. ) 
    260263 
    … …  
    264267!CDIR COLLAPSE 
    265268      wndm(:,:) = sf(jp_wndm)%fnow(:,:,1) 
     269      wndm(:,:) = wndm(:,:) * tmask(:,:,1) 
    266270 
    267271      !------------------------------------------------! 
    … …  
    270274       
    271275      CALL blk_clio_qsr_oce( qsr ) 
    272  
     276      qsr(:,:) = qsr(:,:) * tmask(:,:,1) ! no shortwave radiation into the ocean beneath ice shelf 
    273277      !------------------------! 
    274278      !   Other ocean fluxes   ! 
    … …  
    376380         &     - zqla(:,:)             * pst(:,:) * zcevap                &   ! remove evap.   heat content at SST in Celcius 
    377381         &     + sf(jp_prec)%fnow(:,:,1) * sf(jp_tair)%fnow(:,:,1) * zcprec   ! add    precip. heat content at Tair in Celcius 
     382      qns(:,:) = qns(:,:) * tmask(:,:,1) 
    378383      ! NB: if sea-ice model, the snow precip are computed and the associated heat is added to qns (see blk_ice_clio) 
    379384 

  • branches/2014/dev_r4650_UKMO2_ice_shelves/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_mfs.F90

    r4624 r4724  
    233233         ! Interpolate utau, vtau into the grid_V and grid_V 
    234234         !------------------------------------------------- 
    235  
     235      !     Note the use of 0.5*(2-umask) in order to unmask the stress along coastlines 
     236      !     Note the use of MAX(tmask(i,j),tmask(i+1,j) is to mask tau over ice shelves 
    236237         DO jj = 1, jpjm1 
    237238            DO ji = 1, fs_jpim1 
    238239               utau(ji,jj) = 0.5 * ( 2. - umask(ji,jj,1) ) * ( utau(ji,jj) * tmask(ji,jj,1) & 
    239                &                                + utau(ji+1,jj) * tmask(ji+1,jj,1) ) 
     240               &                                + utau(ji+1,jj) * tmask(ji+1,jj,1) )        & 
     241               &                 * MAX(tmask(ji,jj,1),tmask(ji+1,jj  ,1)) 
    240242               vtau(ji,jj) = 0.5 * ( 2. - vmask(ji,jj,1) ) * ( vtau(ji,jj) * tmask(ji,jj,1) & 
    241                &                                + vtau(ji,jj+1) * tmask(ji,jj+1,1) ) 
     243               &                                + vtau(ji,jj+1) * tmask(ji,jj+1,1) )        & 
     244               &                 * MAX(tmask(ji,jj,1),tmask(ji  ,jj+1,1)) 
    242245            END DO 
    243246         END DO 

  • branches/2014/dev_r4650_UKMO2_ice_shelves/NEMOGCM/NEMO/OPA_SRC/SBC/sbcflx.F90

    r4624 r4724  
    156156            END DO 
    157157         END DO 
     158         taum(:,:) = taum(:,:) * tmask(:,:,1) ; wndm(:,:) = wndm(:,:) * tmask(:,:,1) 
    158159         CALL lbc_lnk( taum(:,:), 'T', 1. )   ;   CALL lbc_lnk( wndm(:,:), 'T', 1. ) 
    159160 

  • branches/2014/dev_r4650_UKMO2_ice_shelves/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim.F90

    r4333 r4724  
    179179         !                                           !----------------! 
    180180         ! 
    181          u_oce(:,:) = ssu_m(:,:)                     ! mean surface ocean current at ice velocity point 
    182          v_oce(:,:) = ssv_m(:,:)                     ! (C-grid dynamics :  U- & V-points as the ocean) 
     181         u_oce(:,:) = ssu_m(:,:) * umask(:,:,1)                     ! mean surface ocean current at ice velocity point 
     182         v_oce(:,:) = ssv_m(:,:) * vmask(:,:,1)                    ! (C-grid dynamics :  U- & V-points as the ocean) 
    183183         ! 
    184184         t_bo(:,:) = tfreez( sss_m ) +  rt0          ! masked sea surface freezing temperature [Kelvin] 

  • branches/2014/dev_r4650_UKMO2_ice_shelves/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim_2.F90

    r4621 r4724  
    126126            DO jj = 2, jpj 
    127127               DO ji = 2, jpi   ! NO vector opt. possible 
    128                   u_oce(ji,jj) = 0.5_wp * ( ssu_m(ji-1,jj  ) + ssu_m(ji-1,jj-1) ) * tmu(ji,jj) 
    129                   v_oce(ji,jj) = 0.5_wp * ( ssv_m(ji  ,jj-1) + ssv_m(ji-1,jj-1) ) * tmu(ji,jj) 
     128                  u_oce(ji,jj) = 0.5_wp * ( ssu_m(ji-1,jj  ) * umask(ji-1,jj  ,1) + ssu_m(ji-1,jj-1) * umask(ji-1,jj-1,1) ) * tmu(ji,jj) 
     129                  v_oce(ji,jj) = 0.5_wp * ( ssv_m(ji  ,jj-1) * vmask(ji  ,jj-1,1) + ssv_m(ji-1,jj-1) * vmask(ji-1,jj-1,1) ) * tmu(ji,jj) 
    130130               END DO 
    131131            END DO 
    … …  
    134134            ! 
    135135         CASE( 'C' )                  !== C-grid ice dynamics :   U & V-points (same as ocean) 
    136             u_oce(:,:) = ssu_m(:,:)                     ! mean surface ocean current at ice velocity point 
    137             v_oce(:,:) = ssv_m(:,:) 
     136            u_oce(:,:) = ssu_m(:,:) * umask(:,:,1)                     ! mean surface ocean current at ice velocity point 
     137            v_oce(:,:) = ssv_m(:,:) * vmask(:,:,1) 
    138138            ! 
    139139         END SELECT 

  • branches/2014/dev_r4650_UKMO2_ice_shelves/NEMOGCM/NEMO/OPA_SRC/TRA/traadv.F90

    r4624 r4724  
    3232   USE wrk_nemo        ! Memory Allocation 
    3333   USE timing          ! Timing 
     34   USE sbc_oce 
    3435 
    3536 
    … …  
    201202      IF( lk_esopa         )   ioptio =          1 
    202203 
     204      IF( ( ln_traadv_muscl .OR. ln_traadv_muscl2 .OR. ln_traadv_ubs .OR. ln_traadv_qck ) .AND. nn_isf .NE. 0 )  & 
     205      &   CALL ctl_stop( 'Only traadv_cen2 and traadv_tvd is compatible with ice shelf cavity') 
     206 
    203207      IF( ioptio /= 1 )   CALL ctl_stop( 'Choose ONE advection scheme in namelist namtra_adv' ) 
    204208 

  • branches/2014/dev_r4650_UKMO2_ice_shelves/NEMOGCM/NEMO/OPA_SRC/TRA/trasbc.F90

    r4666 r4724  
    233233                  zt_frz = -1.9 !tfreez1D( tsn(ji,jj,jk,jp_sal), zpress ) 
    234234               ! compute trend 
    235                   tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem)                                        & 
    236                      &           + zfact * (risf_tsc_b(ji,jj,jp_tem) + risf_tsc(ji,jj,jp_tem)        & 
     235                  tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem)                                          & 
     236                     &           + zfact * (risf_tsc_b(ji,jj,jp_tem) + risf_tsc(ji,jj,jp_tem)          & 
    237237                     &               - rdivisf * (fwfisf(ji,jj) + fwfisf_b(ji,jj)) * zt_frz * r1_rau0) & 
    238238                     &           * r1_hisf_tbl(ji,jj) 
    239                   tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal)                                        & 
     239                  tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal)                                          & 
    240240                     &           + zfact * (risf_tsc_b(ji,jj,jp_sal) + risf_tsc(ji,jj,jp_sal)) * r1_hisf_tbl(ji,jj) 
    241241               END DO 
    … …  
    243243               ! level partially include in ice shelf boundary layer  
    244244               zhk   = SUM( fse3t(ji, jj, ikt:ikb - 1)) * r1_hisf_tbl(ji,jj)  ! proportion of tbl cover by cell from ikt to ikb - 1 
    245                zalpha = rhisf_tbl(ji,jj) * ( 1 - zhk ) / fse3t(ji,jj,ikb)     ! proportion of bottom cell influenced by boundary layer 
     245               zalpha = rhisf_tbl(ji,jj) * ( 1._wp - zhk ) / fse3t(ji,jj,ikb)     ! proportion of bottom cell influenced by boundary layer 
    246246               ! compute tfreez for the temperature correction (we add water at freezing temperature) 
    247247               zpress = grav*rau0*fsdept(ji,jj,ikb)*1.e-04 
    248248               zt_frz = -1.9 !tfreez1D( tsn(ji,jj,ikb,jp_sal), zpress ) 
    249249               ! compute trend 
    250                tsa(ji,jj,ikb,jp_tem) = tsa(ji,jj,ikb,jp_tem)                                         & 
    251                   &              + zfact * (risf_tsc_b(ji,jj,jp_tem) + risf_tsc(ji,jj,jp_tem)        & 
     250               tsa(ji,jj,ikb,jp_tem) = tsa(ji,jj,ikb,jp_tem)                                           & 
     251                  &              + zfact * (risf_tsc_b(ji,jj,jp_tem) + risf_tsc(ji,jj,jp_tem)          & 
    252252                  &                  - rdivisf * (fwfisf(ji,jj) + fwfisf_b(ji,jj)) * zt_frz * r1_rau0) &  
    253253                  &              * r1_hisf_tbl(ji,jj) * zalpha 
    254                tsa(ji,jj,ikb,jp_sal) = tsa(ji,jj,ikb,jp_sal)                                         & 
     254               tsa(ji,jj,ikb,jp_sal) = tsa(ji,jj,ikb,jp_sal)                                           & 
    255255                  &              + zfact * (risf_tsc_b(ji,jj,jp_sal) + risf_tsc(ji,jj,jp_sal)) * r1_hisf_tbl(ji,jj) * zalpha  
    256256            END DO 
    … …  
    261261               &                    'at it= ', kt,' date= ', ndastp 
    262262            IF(lwp) WRITE(numout,*) '~~~~' 
    263             CALL iom_rstput( kt, nitrst, numrow, 'fwf_isf_b', fwfisf(:,:) ) 
    264             CALL iom_rstput( kt, nitrst, numrow, 'isf_hc_b', risf_tsc(:,:,jp_tem) ) 
    265             CALL iom_rstput( kt, nitrst, numrow, 'isf_sc_b', risf_tsc(:,:,jp_sal) ) 
     263            CALL iom_rstput( kt, nitrst, numrow, 'fwf_isf_b', fwfisf(:,:)          ) 
     264            CALL iom_rstput( kt, nitrst, numrow, 'isf_hc_b' , risf_tsc(:,:,jp_tem) ) 
     265            CALL iom_rstput( kt, nitrst, numrow, 'isf_sc_b' , risf_tsc(:,:,jp_sal) ) 
    266266         ENDIF 
    267267      END IF 

  • branches/2014/dev_r4650_UKMO2_ice_shelves/NEMOGCM/NEMO/OPA_SRC/TRA/zpshde.F90

    r4666 r4724  
    8181      !!          di(rho) = rd~ - rd(i,j,k) or rd(i+1,j,k) - rd~ 
    8282      !! 
    83       !! ** Action  : - pgtu, pgtv: horizontal gradient of tracer at u- & v-points 
    84       !!              - pgru, pgrv: horizontal gradient of rho (if present) at u- & v-points  
     83      !! ** Action  : compute for top and bottom interfaces 
     84      !!              - pgtu, pgtv, sgtu, sgtv: horizontal gradient of tracer at u- & v-points 
     85      !!              - pgru, pgrv, sgru, sgtv: horizontal gradient of rho (if present) at u- & v-points 
     86      !!              - pmru, pmrv, smru, smrv: horizontal sum of rho at u- & v- point (used in dynhpg with vvl) 
     87      !!              - pgzu, pgzv, sgzu, sgzv: horizontal gradient of z at u- and v- point (used in dynhpg with vvl) 
     88      !!              - pge3ru, pge3rv, sge3ru, sge3rv: horizontal gradient of rho weighted by local e3w at u- & v-points  
    8589      !!---------------------------------------------------------------------- 
    8690      ! 

  • branches/2014/dev_r4650_UKMO2_ice_shelves/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfbfr.F90

    r4704 r4724  
    188188         ! 
    189189         CALL lbc_lnk( bfrua, 'U', 1. )   ;   CALL lbc_lnk( bfrva, 'V', 1. )      ! Lateral boundary condition 
     190         CALL lbc_lnk( tfrua, 'U', 1. )   ;   CALL lbc_lnk( tfrva, 'V', 1. )      ! Lateral boundary condition 
    190191         ! 
    191192         IF(ln_ctl)   CALL prt_ctl( tab2d_1=bfrua, clinfo1=' bfr  - u: ', mask1=umask,        & 

  • branches/2014/dev_r4650_UKMO2_ice_shelves/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfini.F90

    r4624 r4724  
    1414   !!---------------------------------------------------------------------- 
    1515   USE par_oce         ! mesh and scale factors 
     16   USE sbc_oce         ! surface module (only for nn_isf in the option compatibility test) 
    1617   USE ldftra_oce      ! ocean active tracers: lateral physics 
    1718   USE ldfdyn_oce      ! ocean dynamics lateral physics 
    … …  
    117118      IF( ioptio == 0 .OR. ioptio > 1 .AND. .NOT. lk_esopa )   & 
    118119         &   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 )   & 
     121         &   CALL ctl_stop( ' only zdfcst and zdftke were tested with ice shelves cavities ' ) 
    119122      ! 
    120123      !                               ! ... Convection 

  • branches/2014/dev_r4650_UKMO2_ice_shelves/NEMOGCM/NEMO/OPA_SRC/oce.F90

    r4666 r4724  
    5757   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   gzui, gzvi   !: horizontal gradient of z           at top v-point   
    5858   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   ge3rui, ge3rvi   !: horizontal gradient of T, S and rd at top v-point   
     59   !! (ISF) ice load 
     60   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   riceload 
    5961 
    6062   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   rke          !: kinetic energy 
    … …  
    109111         &     gzui(jpi,jpj)     , gzvi(jpi,jpj)     ,                     & 
    110112         &     ge3rui(jpi,jpj)   , ge3rvi(jpi,jpj)   ,                     & 
    111          &     grui(jpi,jpj)     , grvi(jpi,jpj)     ,                           STAT=ierr(2) ) 
     113         &     grui(jpi,jpj)     , grvi(jpi,jpj)     ,                     & 
     114         &     riceload(jpi,jpj),                             STAT=ierr(2) ) 
    112115         ! 
    113116      ALLOCATE( snwice_mass(jpi,jpj) , snwice_mass_b(jpi,jpj), snwice_fmass(jpi,jpj) , STAT=ierr(3) ) 

  • branches/2014/dev_r4650_UKMO2_ice_shelves/NEMOGCM/NEMO/OPA_SRC/step.F90

    r4666 r4724  
    252252            &                               gtui, gtvi, grui, grvi, arui, arvi, gzui, gzvi, ge3rui, ge3rvi  )  ! of t, s, rd at the last ocean level 
    253253      ELSE                                                  ! centered hpg  (eos then time stepping) 
    254          IF ( .NOT. lk_dynspg_ts ) THEN                     ! eos already called in time-split case 
    255                                 CALL eos    ( tsn, rhd, rhop, fsdept_n(:,:,:) )  ! now in situ density for hpg computation 
    256             IF( ln_zps      )   CALL zps_hde( kstp, jpts, tsn, gtsu, gtsv,  &        ! Partial steps: before horizontal gradient 
    257             &                                      rhd, gru , grv , aru , arv , gzu , gzv , ge3ru , ge3rv ,   &             ! 
    258             &                               gtui, gtvi, grui, grvi, arui, arvi, gzui, gzvi, ge3rui, ge3rvi  )  ! of t, s, rd at the last ocean level 
    259          ENDIF 
    260254         IF( ln_zdfnpc   )   CALL tra_npc( kstp )                ! update after fields by non-penetrative convection 
    261255                             CALL tra_nxt( kstp )                ! tracer fields at next time step 
     256         IF ( .NOT. lk_dynspg_ts ) THEN                     ! eos already called in time-split case 
     257                             CALL eos    ( tsb, rhd, rhop, fsdept_n(:,:,:) )  ! now in situ density for hpg computation 
     258         IF( ln_zps      )   CALL zps_hde( kstp, jpts, tsb, gtsu, gtsv,  &        ! Partial steps: before horizontal gradient 
     259         &                                      rhd, gru , grv , aru , arv , gzu , gzv , ge3ru , ge3rv ,   &             ! 
     260         &                               gtui, gtvi, grui, grvi, arui, arvi, gzui, gzvi, ge3rui, ge3rvi  )  ! of t, s, rd at the last ocean level 
     261         ENDIF 
    262262      ENDIF 
    263263 
Note: See TracChangeset for help on using the changeset viewer.