Changeset 9271


Ignore:
Timestamp:
2018-01-19T18:56:15+01:00 (3 years ago)
Author:
clem
Message:

first steps for having more than 1 snow layers in the ice (in theory). There is still icethd_dh.F90 routine to change

Location:
branches/2017/dev_merge_2017/NEMOGCM/NEMO/LIM_SRC_3
Files:
7 edited

Legend:

Unmodified
Added
Removed
  • branches/2017/dev_merge_2017/NEMOGCM/NEMO/LIM_SRC_3/icectl.F90

    r9019 r9271  
    475475                  WRITE(numout,*) ' v_i           : ', v_i(ji,jj,jl) 
    476476                  WRITE(numout,*) ' v_s           : ', v_s(ji,jj,jl) 
    477                   WRITE(numout,*) ' e_s           : ', e_s(ji,jj,1,jl) 
     477                  WRITE(numout,*) ' e_s           : ', e_s(ji,jj,1:nlay_s,jl) 
    478478                  WRITE(numout,*) ' e_i           : ', e_i(ji,jj,1:nlay_i,jl) 
    479479                  WRITE(numout,*) ' t_su          : ', t_su(ji,jj,jl) 
    480                   WRITE(numout,*) ' t_snow        : ', t_s(ji,jj,1,jl) 
     480                  WRITE(numout,*) ' t_snow        : ', t_s(ji,jj,1:nlay_s,jl) 
    481481                  WRITE(numout,*) ' t_i           : ', t_i(ji,jj,1:nlay_i,jl) 
    482482                  WRITE(numout,*) ' s_i           : ', s_i(ji,jj,jl) 
     
    524524                  WRITE(numout,*) ' h_i        : ', h_i(ji,jj,jl)              , ' h_s        : ', h_s(ji,jj,jl) 
    525525                  WRITE(numout,*) ' t_i        : ', t_i(ji,jj,1:nlay_i,jl) 
    526                   WRITE(numout,*) ' t_su       : ', t_su(ji,jj,jl)             , ' t_s        : ', t_s(ji,jj,1,jl) 
     526                  WRITE(numout,*) ' t_su       : ', t_su(ji,jj,jl)             , ' t_s        : ', t_s(ji,jj,1:nlay_s,jl) 
    527527                  WRITE(numout,*) ' s_i        : ', s_i(ji,jj,jl)              , ' o_i        : ', o_i(ji,jj,jl) 
    528528                  WRITE(numout,*) ' a_i        : ', a_i(ji,jj,jl)              , ' a_i_b      : ', a_i_b(ji,jj,jl)    
  • branches/2017/dev_merge_2017/NEMOGCM/NEMO/LIM_SRC_3/icedyn_adv_pra.F90

    r9123 r9271  
    3838   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   sxsn , sysn , sxxsn , syysn , sxysn    ! snow thickness 
    3939   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   sxa  , sya  , sxxa  , syya  , sxya     ! lead fraction 
    40    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   sxc0 , syc0 , sxxc0 , syyc0 , sxyc0    ! snow thermal content 
    4140   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   sxsal, sysal, sxxsal, syysal, sxysal   ! ice salinity 
    4241   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   sxage, syage, sxxage, syyage, sxyage   ! ice age 
    4342   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)     ::   sxopw, syopw, sxxopw, syyopw, sxyopw   ! open water in sea ice 
     43   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) ::   sxc0 , syc0 , sxxc0 , syyc0 , sxyc0    ! snow layers heat content 
    4444   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) ::   sxe  , sye  , sxxe  , syye  , sxye     ! ice layers heat content 
    4545   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   sxap , syap , sxxap , syyap , sxyap    ! melt pond fraction 
     
    8787      REAL(wp), ALLOCATABLE, DIMENSION(:,:)     ::   zarea 
    8888      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:)   ::   z0opw 
    89       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:)   ::   z0ice, z0snw, z0ai, z0es , z0smi , z0oi 
     89      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:)   ::   z0ice, z0snw, z0ai, z0smi, z0oi 
    9090      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:)   ::   z0ap , z0vp 
     91      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) ::   z0es 
    9192      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) ::   z0ei 
    9293      !!---------------------------------------------------------------------- 
     
    9495      IF( kt == nit000 .AND. lwp )   WRITE(numout,*) '-- ice_dyn_adv_pra: Prather advection scheme' 
    9596      ! 
    96       ALLOCATE( zarea(jpi,jpj)     , z0opw(jpi,jpj, 1 ) ,                                           & 
    97          &      z0ice(jpi,jpj,jpl) , z0snw(jpi,jpj,jpl) , z0ai(jpi,jpj,jpl) , z0es(jpi,jpj,jpl) ,   & 
    98          &      z0smi(jpi,jpj,jpl) , z0oi (jpi,jpj,jpl) , z0ap(jpi,jpj,jpl) , z0vp(jpi,jpj,jpl) ,   & 
    99          &      z0ei (jpi,jpj,nlay_i,jpl)                                                           ) 
     97      ALLOCATE( zarea(jpi,jpj)    , z0opw(jpi,jpj, 1 ) , z0ice(jpi,jpj,jpl) , z0snw(jpi,jpj,jpl) ,                       & 
     98         &      z0ai(jpi,jpj,jpl) , z0smi(jpi,jpj,jpl) , z0oi (jpi,jpj,jpl) , z0ap (jpi,jpj,jpl) , z0vp(jpi,jpj,jpl) ,   & 
     99         &      z0es (jpi,jpj,nlay_s,jpl), z0ei(jpi,jpj,nlay_i,jpl) ) 
    100100      ! 
    101101      ! --- If ice drift field is too fast, use an appropriate time step for advection (CFL test for stability) --- !         
     
    119119         z0smi(:,:,jl) = psv_i(:,:,  jl) * e1e2t(:,:)     ! Salt content 
    120120         z0oi (:,:,jl) = poa_i(:,:,  jl) * e1e2t(:,:)     ! Age content 
    121          z0es (:,:,jl) = pe_s (:,:,1,jl) * e1e2t(:,:)     ! Snow heat content 
     121         DO jk = 1, nlay_s 
     122            z0es(:,:,jk,jl) = pe_s(:,:,jk,jl) * e1e2t(:,:) ! Snow heat content 
     123         END DO 
    122124         DO jk = 1, nlay_i 
    123125            z0ei(:,:,jk,jl) = pe_i(:,:,jk,jl) * e1e2t(:,:) ! Ice  heat content 
     
    158160               CALL adv_y( zusnit, pv_ice, 0._wp, zarea, z0ai  (:,:,jl), sxa  (:,:,jl),   &  
    159161                  &                                      sxxa  (:,:,jl), sya  (:,:,jl), syya  (:,:,jl), sxya  (:,:,jl)  ) 
    160                CALL adv_x( zusnit, pu_ice, 1._wp, zarea, z0es  (:,:,jl), sxc0 (:,:,jl),   &    !--- snow heat contents --- 
    161                   &                                      sxxc0 (:,:,jl), syc0 (:,:,jl), syyc0 (:,:,jl), sxyc0 (:,:,jl)  ) 
    162                CALL adv_y( zusnit, pv_ice, 0._wp, zarea, z0es  (:,:,jl), sxc0 (:,:,jl),   & 
    163                   &                                      sxxc0 (:,:,jl), syc0 (:,:,jl), syyc0 (:,:,jl), sxyc0 (:,:,jl)  ) 
     162               DO jk = 1, nlay_s                                                               !--- snow heat contents --- 
     163                  CALL adv_x( zusnit, pu_ice, 1._wp, zarea, z0es (:,:,jk,jl), sxc0(:,:,jk,jl),   & 
     164                     &                                      sxxc0(:,:,jk,jl), syc0(:,:,jk,jl), syyc0(:,:,jk,jl), sxyc0(:,:,jk,jl) ) 
     165                  CALL adv_y( zusnit, pv_ice, 0._wp, zarea, z0es (:,:,jk,jl), sxc0(:,:,jk,jl),   & 
     166                     &                                      sxxc0(:,:,jk,jl), syc0(:,:,jk,jl), syyc0(:,:,jk,jl), sxyc0(:,:,jk,jl) ) 
     167               END DO 
    164168               DO jk = 1, nlay_i                                                               !--- ice heat contents --- 
    165                   CALL adv_x( zusnit, pu_ice, 1._wp, zarea, z0ei(:,:,jk,jl), sxe (:,:,jk,jl),   &  
    166                      &                                      sxxe(:,:,jk,jl), sye (:,:,jk,jl),   & 
    167                      &                                      syye(:,:,jk,jl), sxye(:,:,jk,jl) ) 
    168                   CALL adv_y( zusnit, pv_ice, 0._wp, zarea, z0ei(:,:,jk,jl), sxe (:,:,jk,jl),   &  
    169                      &                                      sxxe(:,:,jk,jl), sye (:,:,jk,jl),   & 
    170                      &                                      syye(:,:,jk,jl), sxye(:,:,jk,jl) ) 
     169                  CALL adv_x( zusnit, pu_ice, 1._wp, zarea, z0ei(:,:,jk,jl), sxe(:,:,jk,jl),   &  
     170                     &                                      sxxe(:,:,jk,jl), sye(:,:,jk,jl), syye(:,:,jk,jl), sxye(:,:,jk,jl) ) 
     171                  CALL adv_y( zusnit, pv_ice, 0._wp, zarea, z0ei(:,:,jk,jl), sxe(:,:,jk,jl),   &  
     172                     &                                      sxxe(:,:,jk,jl), sye(:,:,jk,jl), syye(:,:,jk,jl), sxye(:,:,jk,jl) ) 
    171173               END DO 
    172174               IF ( ln_pnd_H12 ) THEN 
     
    211213               CALL adv_x( zusnit, pu_ice, 0._wp, zarea, z0ai  (:,:,jl), sxa  (:,:,jl),   & 
    212214                  &                                      sxxa  (:,:,jl), sya  (:,:,jl), syya  (:,:,jl), sxya  (:,:,jl)  ) 
    213                CALL adv_y( zusnit, pv_ice, 1._wp, zarea, z0es  (:,:,jl), sxc0 (:,:,jl),   &  !--- snow heat contents --- 
    214                   &                                      sxxc0 (:,:,jl), syc0 (:,:,jl), syyc0 (:,:,jl), sxyc0 (:,:,jl)  ) 
    215                CALL adv_x( zusnit, pu_ice, 0._wp, zarea, z0es  (:,:,jl), sxc0 (:,:,jl),   & 
    216                   &                                      sxxc0 (:,:,jl), syc0 (:,:,jl), syyc0 (:,:,jl), sxyc0 (:,:,jl)  ) 
     215               DO jk = 1, nlay_s                                                             !--- snow heat contents --- 
     216                  CALL adv_y( zusnit, pv_ice, 1._wp, zarea, z0es (:,:,jk,jl), sxc0(:,:,jk,jl),   & 
     217                     &                                      sxxc0(:,:,jk,jl), syc0(:,:,jk,jl), syyc0(:,:,jk,jl), sxyc0(:,:,jk,jl) ) 
     218                  CALL adv_x( zusnit, pu_ice, 0._wp, zarea, z0es (:,:,jk,jl), sxc0(:,:,jk,jl),   & 
     219                     &                                      sxxc0(:,:,jk,jl), syc0(:,:,jk,jl), syyc0(:,:,jk,jl), sxyc0(:,:,jk,jl) ) 
     220               END DO 
    217221               DO jk = 1, nlay_i                                                             !--- ice heat contents --- 
    218                   CALL adv_y( zusnit, pv_ice, 1._wp, zarea, z0ei(:,:,jk,jl), sxe (:,:,jk,jl),   &  
    219                      &                                      sxxe(:,:,jk,jl), sye (:,:,jk,jl),   & 
    220                      &                                      syye(:,:,jk,jl), sxye(:,:,jk,jl) ) 
    221                   CALL adv_x( zusnit, pu_ice, 0._wp, zarea, z0ei(:,:,jk,jl), sxe (:,:,jk,jl),   &  
    222                      &                                      sxxe(:,:,jk,jl), sye (:,:,jk,jl),   & 
    223                      &                                      syye(:,:,jk,jl), sxye(:,:,jk,jl) ) 
     222                  CALL adv_y( zusnit, pv_ice, 1._wp, zarea, z0ei(:,:,jk,jl), sxe(:,:,jk,jl),   &  
     223                     &                                      sxxe(:,:,jk,jl), sye(:,:,jk,jl), syye(:,:,jk,jl), sxye(:,:,jk,jl) ) 
     224                  CALL adv_x( zusnit, pu_ice, 0._wp, zarea, z0ei(:,:,jk,jl), sxe(:,:,jk,jl),   &  
     225                     &                                      sxxe(:,:,jk,jl), sye(:,:,jk,jl), syye(:,:,jk,jl), sxye(:,:,jk,jl) ) 
    224226               END DO 
    225227               IF ( ln_pnd_H12 ) THEN 
     
    247249         poa_i(:,:,  jl) = z0oi (:,:,jl) * r1_e1e2t(:,:) 
    248250         pa_i (:,:,  jl) = z0ai (:,:,jl) * r1_e1e2t(:,:) 
    249          pe_s (:,:,1,jl) = z0es (:,:,jl) * r1_e1e2t(:,:) 
     251         DO jk = 1, nlay_s 
     252            pe_s(:,:,jk,jl) = z0es(:,:,jk,jl) * r1_e1e2t(:,:) 
     253         END DO 
    250254         DO jk = 1, nlay_i 
    251255            pe_i(:,:,jk,jl) = z0ei(:,:,jk,jl) * r1_e1e2t(:,:) 
     
    257261      END DO 
    258262      ! 
    259       DEALLOCATE( zarea , z0opw , z0ice, z0snw , z0ai , z0es , z0smi , z0oi , z0ap , z0vp , z0ei ) 
     263      DEALLOCATE( zarea , z0opw , z0ice, z0snw , z0ai , z0smi , z0oi , z0ap , z0vp , z0es, z0ei ) 
    260264      ! 
    261265      IF( lrst_ice )   CALL adv_pra_rst( 'WRITE', kt )   !* write Prather fields in the restart file 
     
    625629         &      sxsn (jpi,jpj,jpl) , sysn (jpi,jpj,jpl) , sxxsn (jpi,jpj,jpl) , syysn (jpi,jpj,jpl) , sxysn (jpi,jpj,jpl) ,   & 
    626630         &      sxa  (jpi,jpj,jpl) , sya  (jpi,jpj,jpl) , sxxa  (jpi,jpj,jpl) , syya  (jpi,jpj,jpl) , sxya  (jpi,jpj,jpl) ,   & 
    627          &      sxc0 (jpi,jpj,jpl) , syc0 (jpi,jpj,jpl) , sxxc0 (jpi,jpj,jpl) , syyc0 (jpi,jpj,jpl) , sxyc0 (jpi,jpj,jpl) ,   & 
    628631         &      sxsal(jpi,jpj,jpl) , sysal(jpi,jpj,jpl) , sxxsal(jpi,jpj,jpl) , syysal(jpi,jpj,jpl) , sxysal(jpi,jpj,jpl) ,   & 
    629632         &      sxage(jpi,jpj,jpl) , syage(jpi,jpj,jpl) , sxxage(jpi,jpj,jpl) , syyage(jpi,jpj,jpl) , sxyage(jpi,jpj,jpl) ,   & 
    630633         &      sxap(jpi,jpj,jpl)  , syap (jpi,jpj,jpl) , sxxap (jpi,jpj,jpl) , syyap (jpi,jpj,jpl) , sxyap (jpi,jpj,jpl) ,   & 
    631634         &      sxvp(jpi,jpj,jpl)  , syvp (jpi,jpj,jpl) , sxxvp (jpi,jpj,jpl) , syyvp (jpi,jpj,jpl) , sxyvp (jpi,jpj,jpl) ,   & 
    632          &      sxe (jpi,jpj,nlay_i,jpl) , sye (jpi,jpj,nlay_i,jpl) , sxxe(jpi,jpj,nlay_i,jpl) , & 
    633          &      syye(jpi,jpj,nlay_i,jpl) , sxye(jpi,jpj,nlay_i,jpl)                            , & 
     635         ! 
     636         &      sxc0 (jpi,jpj,nlay_s,jpl) , syc0 (jpi,jpj,nlay_s,jpl) , sxxc0(jpi,jpj,nlay_s,jpl) , & 
     637         &      syyc0(jpi,jpj,nlay_s,jpl) , sxyc0(jpi,jpj,nlay_s,jpl)                             , & 
     638         ! 
     639         &      sxe  (jpi,jpj,nlay_i,jpl) , sye  (jpi,jpj,nlay_i,jpl) , sxxe (jpi,jpj,nlay_i,jpl) , & 
     640         &      syye (jpi,jpj,nlay_i,jpl) , sxye (jpi,jpj,nlay_i,jpl)                             , & 
    634641         &      STAT = ierr ) 
    635642      ! 
     
    689696            CALL iom_get( numrir, jpdom_autoglo, 'syya'  , syya   ) 
    690697            CALL iom_get( numrir, jpdom_autoglo, 'sxya'  , sxya   ) 
    691             !                                                        ! snow thermal content 
    692             CALL iom_get( numrir, jpdom_autoglo, 'sxc0'  , sxc0   ) 
    693             CALL iom_get( numrir, jpdom_autoglo, 'syc0'  , syc0   ) 
    694             CALL iom_get( numrir, jpdom_autoglo, 'sxxc0' , sxxc0  ) 
    695             CALL iom_get( numrir, jpdom_autoglo, 'syyc0' , syyc0  ) 
    696             CALL iom_get( numrir, jpdom_autoglo, 'sxyc0' , sxyc0  ) 
    697698            !                                                        ! ice salinity 
    698699            CALL iom_get( numrir, jpdom_autoglo, 'sxsal' , sxsal  ) 
     
    713714            CALL iom_get( numrir, jpdom_autoglo, 'syyopw', syyopw ) 
    714715            CALL iom_get( numrir, jpdom_autoglo, 'sxyopw', sxyopw ) 
     716            !                                                        ! snow layers heat content 
     717            DO jk = 1, nlay_s 
     718               WRITE(zchar1,'(I2.2)') jk 
     719               znam = 'sxc0'//'_l'//zchar1  ;   CALL iom_get( numrir, jpdom_autoglo, znam , z3d )   ;   sxc0 (:,:,jk,:) = z3d(:,:,:) 
     720               znam = 'syc0'//'_l'//zchar1  ;   CALL iom_get( numrir, jpdom_autoglo, znam , z3d )   ;   syc0 (:,:,jk,:) = z3d(:,:,:) 
     721               znam = 'sxxc0'//'_l'//zchar1 ;   CALL iom_get( numrir, jpdom_autoglo, znam , z3d )   ;   sxxc0(:,:,jk,:) = z3d(:,:,:) 
     722               znam = 'syyc0'//'_l'//zchar1 ;   CALL iom_get( numrir, jpdom_autoglo, znam , z3d )   ;   syyc0(:,:,jk,:) = z3d(:,:,:) 
     723               znam = 'sxyc0'//'_l'//zchar1 ;   CALL iom_get( numrir, jpdom_autoglo, znam , z3d )   ;   sxyc0(:,:,jk,:) = z3d(:,:,:) 
     724            END DO 
    715725            !                                                        ! ice layers heat content 
    716             DO jk = 1, nlay_i  
     726            DO jk = 1, nlay_i 
    717727               WRITE(zchar1,'(I2.2)') jk 
    718728               znam = 'sxe'//'_l'//zchar1   ;   CALL iom_get( numrir, jpdom_autoglo, znam , z3d )   ;   sxe (:,:,jk,:) = z3d(:,:,:) 
     
    744754            sxsn  = 0._wp   ;   sysn  = 0._wp   ;   sxxsn  = 0._wp   ;   syysn  = 0._wp   ;   sxysn  = 0._wp      ! snow thickness 
    745755            sxa   = 0._wp   ;   sya   = 0._wp   ;   sxxa   = 0._wp   ;   syya   = 0._wp   ;   sxya   = 0._wp      ! lead fraction 
    746             sxc0  = 0._wp   ;   syc0  = 0._wp   ;   sxxc0  = 0._wp   ;   syyc0  = 0._wp   ;   sxyc0  = 0._wp      ! snow thermal content 
    747756            sxsal = 0._wp   ;   sysal = 0._wp   ;   sxxsal = 0._wp   ;   syysal = 0._wp   ;   sxysal = 0._wp      ! ice salinity 
    748757            sxage = 0._wp   ;   syage = 0._wp   ;   sxxage = 0._wp   ;   syyage = 0._wp   ;   sxyage = 0._wp      ! ice age 
    749758            sxopw = 0._wp   ;   syopw = 0._wp   ;   sxxopw = 0._wp   ;   syyopw = 0._wp   ;   sxyopw = 0._wp      ! open water in sea ice 
     759            sxc0  = 0._wp   ;   syc0  = 0._wp   ;   sxxc0  = 0._wp   ;   syyc0  = 0._wp   ;   sxyc0  = 0._wp      ! snow layers heat content 
    750760            sxe   = 0._wp   ;   sye   = 0._wp   ;   sxxe   = 0._wp   ;   syye   = 0._wp   ;   sxye   = 0._wp      ! ice layers heat content 
    751761            IF( ln_pnd_H12 ) THEN 
     
    783793         CALL iom_rstput( iter, nitrst, numriw, 'syya'  , syya   ) 
    784794         CALL iom_rstput( iter, nitrst, numriw, 'sxya'  , sxya   ) 
    785          !                                                           ! snow thermal content 
    786          CALL iom_rstput( iter, nitrst, numriw, 'sxc0'  , sxc0   ) 
    787          CALL iom_rstput( iter, nitrst, numriw, 'syc0'  , syc0   ) 
    788          CALL iom_rstput( iter, nitrst, numriw, 'sxxc0' , sxxc0  ) 
    789          CALL iom_rstput( iter, nitrst, numriw, 'syyc0' , syyc0  ) 
    790          CALL iom_rstput( iter, nitrst, numriw, 'sxyc0' , sxyc0  ) 
    791795         !                                                           ! ice salinity 
    792796         CALL iom_rstput( iter, nitrst, numriw, 'sxsal' , sxsal  ) 
     
    807811         CALL iom_rstput( iter, nitrst, numriw, 'syyopw', syyopw ) 
    808812         CALL iom_rstput( iter, nitrst, numriw, 'sxyopw', sxyopw ) 
     813         !                                                           ! snow layers heat content 
     814         DO jk = 1, nlay_s 
     815            WRITE(zchar1,'(I2.2)') jk 
     816            znam = 'sxc0'//'_l'//zchar1  ;   z3d(:,:,:) = sxc0 (:,:,jk,:)  ;   CALL iom_rstput( iter, nitrst, numriw, znam , z3d ) 
     817            znam = 'syc0'//'_l'//zchar1  ;   z3d(:,:,:) = syc0 (:,:,jk,:)  ;   CALL iom_rstput( iter, nitrst, numriw, znam , z3d ) 
     818            znam = 'sxxc0'//'_l'//zchar1 ;   z3d(:,:,:) = sxxc0(:,:,jk,:)  ;   CALL iom_rstput( iter, nitrst, numriw, znam , z3d ) 
     819            znam = 'syyc0'//'_l'//zchar1 ;   z3d(:,:,:) = syyc0(:,:,jk,:)  ;   CALL iom_rstput( iter, nitrst, numriw, znam , z3d ) 
     820            znam = 'sxyc0'//'_l'//zchar1 ;   z3d(:,:,:) = sxyc0(:,:,jk,:)  ;   CALL iom_rstput( iter, nitrst, numriw, znam , z3d ) 
     821         END DO 
    809822         !                                                           ! ice layers heat content 
    810          DO jk = 1, nlay_i  
     823         DO jk = 1, nlay_i 
    811824            WRITE(zchar1,'(I2.2)') jk 
    812825            znam = 'sxe'//'_l'//zchar1   ;   z3d(:,:,:) = sxe (:,:,jk,:)   ;   CALL iom_rstput( iter, nitrst, numriw, znam , z3d ) 
  • branches/2017/dev_merge_2017/NEMOGCM/NEMO/LIM_SRC_3/icedyn_adv_umx.F90

    r9019 r9271  
    122122            END DO 
    123123            CALL adv_umx( k_order, kt, zdt, zudy, zvdx, zcu_box, zcv_box, pv_s(:,:,jl) )         ! Snow volume 
    124             CALL adv_umx( k_order, kt, zdt, zudy, zvdx, zcu_box, zcv_box, pe_s(:,:,1,jl) )       ! Snow heat content 
     124            DO jk = 1, nlay_s 
     125               CALL adv_umx( k_order, kt, zdt, zudy, zvdx, zcu_box, zcv_box, pe_s(:,:,jk,jl) )   ! Snow heat content 
     126            END DO 
    125127            IF ( ln_pnd_H12 ) THEN 
    126128               CALL adv_umx( k_order, kt, zdt, zudy, zvdx, zcu_box, zcv_box, pa_ip(:,:,jl) )     ! Melt pond fraction 
  • branches/2017/dev_merge_2017/NEMOGCM/NEMO/LIM_SRC_3/icedyn_rdgrft.F90

    r9169 r9271  
    570570      REAL(wp)                  ::   airdg1, oirdg1, aprdg1, virdg1, sirdg1 
    571571      REAL(wp)                  ::   airft1, oirft1, aprft1 
    572       REAL(wp), DIMENSION(jpij) ::   airdg2, oirdg2, aprdg2, virdg2, sirdg2, vsrdg, vprdg, esrdg  ! area etc of new ridges 
    573       REAL(wp), DIMENSION(jpij) ::   airft2, oirft2, aprft2, virft , sirft , vsrft, vprft, esrft  ! area etc of rafted ice 
     572      REAL(wp), DIMENSION(jpij) ::   airdg2, oirdg2, aprdg2, virdg2, sirdg2, vsrdg, vprdg  ! area etc of new ridges 
     573      REAL(wp), DIMENSION(jpij) ::   airft2, oirft2, aprft2, virft , sirft , vsrft, vprft  ! area etc of rafted ice 
    574574      ! 
    575575      REAL(wp), DIMENSION(jpij) ::   ersw             ! enth of water trapped into ridges 
     
    577577      REAL(wp), DIMENSION(jpij) ::   z1_ai            ! 1 / a 
    578578      ! 
    579       REAL(wp), DIMENSION(jpij,nlay_i) ::   eirft     ! ice energy of rafting ice 
     579      REAL(wp), DIMENSION(jpij,nlay_s) ::   esrft     ! snow energy of rafting ice 
     580      REAL(wp), DIMENSION(jpij,nlay_i) ::   eirft     ! ice  energy of rafting ice 
     581      REAL(wp), DIMENSION(jpij,nlay_s) ::   esrdg     ! enth*volume of new ridges       
    580582      REAL(wp), DIMENSION(jpij,nlay_i) ::   eirdg     ! enth*volume of new ridges 
    581583      !!------------------------------------------------------------------- 
     
    633635               oirdg1     = oa_i_2d(ji,jl1)   * afrdg 
    634636               oirdg2(ji) = oa_i_2d(ji,jl1)   * afrdg * hi_hrdg(ji,jl1)  
    635                esrdg(ji)  = ze_s_2d(ji,1,jl1) * afrdg 
    636637 
    637638               virft(ji)  = v_i_2d (ji,jl1)   * afrft 
     
    640641               oirft1     = oa_i_2d(ji,jl1)   * afrft  
    641642               oirft2(ji) = oa_i_2d(ji,jl1)   * afrft * hi_hrft  
    642                esrft(ji)  = ze_s_2d(ji,1,jl1) * afrft 
    643643 
    644644               IF ( ln_pnd_H12 ) THEN 
     
    663663               wfx_snw_dyn_1d(ji) = wfx_snw_dyn_1d(ji) + ( rhosn * vsrdg(ji) * ( 1._wp - rn_fsnwrdg )   &   ! fresh water source for ocean 
    664664                  &                                      + rhosn * vsrft(ji) * ( 1._wp - rn_fsnwrft ) ) * r1_rdtice 
    665  
    666                hfx_dyn_1d(ji) = hfx_dyn_1d(ji) + ( - esrdg(ji) * ( 1._wp - rn_fsnwrdg )   &                 ! heat sink for ocean (<0, W.m-2) 
    667                   &                                - esrft(ji) * ( 1._wp - rn_fsnwrft ) ) * r1_rdtice 
    668665 
    669666               ! Put the melt pond water into the ocean 
     
    697694                  v_ip_2d(ji,jl1) = v_ip_2d(ji,jl1) - vprdg(ji) - vprft(ji) 
    698695               ENDIF 
    699                ze_s_2d(ji,1,jl1) = ze_s_2d(ji,1,jl1) - esrdg (ji) - esrft (ji) 
    700696            ENDIF 
    701697 
    702698         END DO ! ji 
    703699 
     700         ! special loop for e_s because of layers jk 
     701         DO jk = 1, nlay_s 
     702            DO ji = 1, npti 
     703               IF( apartf(ji,jl1) > 0._wp .AND. closing_gross(ji) > 0._wp ) THEN 
     704                  ! Compute ridging /rafting fractions 
     705                  afrdg = aridge(ji,jl1) * closing_gross(ji) * rdt_ice * z1_ai(ji) 
     706                  afrft = araft (ji,jl1) * closing_gross(ji) * rdt_ice * z1_ai(ji) 
     707                  ! Compute ridging ice and new ridges for es 
     708                  esrdg(ji,jk) = ze_s_2d (ji,jk,jl1) * afrdg 
     709                  esrft(ji,jk) = ze_s_2d (ji,jk,jl1) * afrft 
     710                  ! Put the snow lost by ridging into the ocean 
     711                  hfx_dyn_1d(ji) = hfx_dyn_1d(ji) + ( - esrdg(ji,jk) * ( 1._wp - rn_fsnwrdg )   &                 ! heat sink for ocean (<0, W.m-2) 
     712                     &                                - esrft(ji,jk) * ( 1._wp - rn_fsnwrft ) ) * r1_rdtice 
     713                  ! Update jl1 
     714                  ze_s_2d(ji,jk,jl1) = ze_s_2d(ji,jk,jl1) * ( 1._wp - afrdg - afrft )  
     715               ENDIF 
     716            END DO 
     717         END DO 
     718                   
    704719         ! special loop for e_i because of layers jk 
    705720         DO jk = 1, nlay_i 
     
    754769                  sv_i_2d(ji,jl2) = sv_i_2d(ji,jl2) + ( sirdg2(ji) * fvol(ji) + sirft (ji) * zswitch(ji) ) 
    755770                  v_s_2d (ji,jl2) = v_s_2d (ji,jl2) + ( vsrdg (ji) * rn_fsnwrdg * fvol(ji)  +  & 
    756                      &                                    vsrft (ji) * rn_fsnwrft * zswitch(ji) ) 
     771                     &                                  vsrft (ji) * rn_fsnwrft * zswitch(ji) ) 
    757772                  IF ( ln_pnd_H12 ) THEN 
    758773                     v_ip_2d (ji,jl2) = v_ip_2d(ji,jl2) + (   vprdg (ji) * rn_fpndrdg * fvol   (ji)   & 
     
    761776                        &                                   + aprft2(ji) * rn_fpndrft * zswitch(ji)   ) 
    762777                  ENDIF 
    763                   ze_s_2d(ji,1,jl2) = ze_s_2d(ji,1,jl2) + ( esrdg (ji) * rn_fsnwrdg * fvol(ji)  +  & 
    764                      &                                      esrft (ji) * rn_fsnwrft * zswitch(ji) ) 
    765778                   
    766779               ENDIF 
    767780 
    768781            END DO 
    769  
     782            ! for snow enthalpy 
     783            DO jk = 1, nlay_s 
     784               DO ji = 1, npti 
     785                  IF( apartf(ji,jl1) > 0._wp .AND. closing_gross(ji) > 0._wp )   & 
     786                     &   ze_s_2d(ji,jk,jl2) = ze_s_2d(ji,jk,jl2) + ( esrdg(ji,jk) * rn_fsnwrdg * fvol(ji)  +  & 
     787                     &                                               esrft(ji,jk) * rn_fsnwrft * zswitch(ji) ) 
     788               END DO 
     789            END DO 
     790            ! for ice enthalpy 
    770791            DO jk = 1, nlay_i 
    771792               DO ji = 1, npti 
  • branches/2017/dev_merge_2017/NEMOGCM/NEMO/LIM_SRC_3/icerst.F90

    r9123 r9271  
    132132      CALL iom_rstput( iter, nitrst, numriw, 'a_ip', a_ip ) 
    133133      CALL iom_rstput( iter, nitrst, numriw, 'v_ip', v_ip ) 
    134 !!gm dangerous !!!!!  ===>>>> better reading writing all snow layers ! 
    135       ! Snow enthalpy (1st snow layer only) 
    136       z3d = e_s(:,:,1,:) 
    137       CALL iom_rstput( iter, nitrst, numriw, 'e_s_l01' , z3d ) 
    138       ! Ice enthalpy (all ice layers) 
     134      ! Snow enthalpy 
     135      DO jk = 1, nlay_s  
     136         WRITE(zchar1,'(I2.2)') jk 
     137         znam = 'e_s'//'_l'//zchar1 
     138         z3d(:,:,:) = e_s(:,:,jk,:) 
     139         CALL iom_rstput( iter, nitrst, numriw, znam , z3d ) 
     140      END DO 
     141      ! Ice enthalpy 
    139142      DO jk = 1, nlay_i  
    140143         WRITE(zchar1,'(I2.2)') jk 
     
    219222         v_ip(:,:,:) = 0._wp 
    220223      ENDIF 
    221 !!gm dangerous !!!!!  ===>>>> better reading writing all snow layers ! 
    222       ! Snow enthalpy (1st snow layer only) 
    223       CALL iom_get( numrir, jpdom_autoglo, 'e_s_l01' , z3d ) 
    224       e_s(:,:,1,:) = z3d 
    225       ! Ice enthalpy (all ice layers) 
    226       DO jk = 1, nlay_i  
     224      ! Snow enthalpy 
     225      DO jk = 1, nlay_s 
     226         WRITE(zchar1,'(I2.2)') jk 
     227         znam = 'e_s'//'_l'//zchar1 
     228         CALL iom_get( numrir, jpdom_autoglo, znam , z3d ) 
     229         e_s(:,:,jk,:) = z3d(:,:,:) 
     230      END DO 
     231      ! Ice enthalpy 
     232      DO jk = 1, nlay_i 
    227233         WRITE(zchar1,'(I2.2)') jk 
    228234         znam = 'e_i'//'_l'//zchar1 
  • branches/2017/dev_merge_2017/NEMOGCM/NEMO/LIM_SRC_3/icethd_dh.F90

    r9019 r9271  
    9898      REAL(wp), DIMENSION(jpij) ::   zdh_s_sub   ! snow sublimation 
    9999 
     100      REAL(wp), DIMENSION(jpij,nlay_s) ::   zh_s      ! snw layer thickness 
     101      REAL(wp), DIMENSION(jpij,nlay_i) ::   zh_i      ! ice layer thickness 
    100102      REAL(wp), DIMENSION(jpij,nlay_i) ::   zdeltah 
    101       REAL(wp), DIMENSION(jpij,nlay_i) ::   zh_i      ! ice layer thickness 
    102103      INTEGER , DIMENSION(jpij,nlay_i) ::   icount    ! number of layers vanished by melting  
    103104 
    104       REAL(wp), DIMENSION(jpij) ::   zeh_i       ! total ice heat content  (J.m-2) 
    105105      REAL(wp), DIMENSION(jpij) ::   zsnw        ! distribution of snow after wind blowing 
    106106 
     
    161161      END DO 
    162162 
     163      !--------------------------------------------! 
     164      !  2) Computing layer thicknesses            ! 
     165      !--------------------------------------------! 
     166      DO jk = 1, nlay_i 
     167         DO ji = 1, npti 
     168            zh_i(ji,jk) = h_i_1d(ji) * r1_nlay_i 
     169         END DO 
     170      END DO 
     171 
     172      DO jk = 1, nlay_s 
     173         DO ji = 1, npti 
     174            zh_s(ji,jk) = h_s_1d(ji) * r1_nlay_s 
     175         END DO 
     176      END DO 
     177 
    163178      !------------------------------------------------------------------------------! 
    164179      ! If snow temperature is above freezing point, then snow melts  
    165180      ! (should not happen but sometimes it does) 
    166181      !------------------------------------------------------------------------------! 
    167       DO ji = 1, npti 
    168          IF( t_s_1d(ji,1) > rt0 ) THEN !!! Internal melting 
    169             ! Contribution to heat flux to the ocean [W.m-2], < 0   
    170             hfx_res_1d(ji) = hfx_res_1d(ji) + e_s_1d(ji,1) * h_s_1d(ji) * a_i_1d(ji) * r1_rdtice 
    171             ! Contribution to mass flux 
    172             wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) + rhosn * h_s_1d(ji) * a_i_1d(ji) * r1_rdtice 
    173             dh_s_mlt(ji)       = dh_s_mlt(ji) - h_s_1d(ji) 
    174             ! updates 
    175             h_s_1d(ji)    = 0._wp 
    176             e_s_1d (ji,1) = 0._wp 
    177             t_s_1d (ji,1) = rt0 
    178          END IF 
    179       END DO 
    180  
    181       !------------------------------------------------------------! 
    182       !  2) Computing layer thicknesses and enthalpies.            ! 
    183       !------------------------------------------------------------! 
    184       DO jk = 1, nlay_i 
    185          DO ji = 1, npti 
    186             zh_i(ji,jk) = h_i_1d(ji) * r1_nlay_i 
    187             zeh_i(ji)   = zeh_i(ji) + e_i_1d(ji,jk) * zh_i(ji,jk) 
    188          END DO 
    189       END DO 
     182      DO jk = 1, nlay_s 
     183         DO ji = 1, npti 
     184            IF( t_s_1d(ji,jk) > rt0 ) THEN !!! Internal melting 
     185               ! Contribution to heat flux to the ocean [W.m-2], < 0   
     186               hfx_res_1d(ji) = hfx_res_1d(ji) + e_s_1d(ji,jk) * zh_s(ji,jk) * a_i_1d(ji) * r1_rdtice 
     187               ! Contribution to mass flux 
     188               wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) + rhosn * zh_s(ji,jk) * a_i_1d(ji) * r1_rdtice 
     189               dh_s_mlt(ji)       = dh_s_mlt(ji) - zh_s(ji,jk) 
     190               ! updates 
     191               h_s_1d(ji)    = h_s_1d(ji) - zh_s(ji,jk) 
     192               zh_s  (ji,jk) = 0._wp 
     193               e_s_1d(ji,jk) = 0._wp 
     194               t_s_1d(ji,jk) = rt0 
     195            END IF 
     196         END DO 
     197      END DO 
     198          
    190199 
    191200      !------------------------------------------------------------------------------| 
     
    244253      END DO 
    245254 
     255      ! recalculate snow layers 
     256      DO jk = 1, nlay_s 
     257         DO ji = 1, npti 
     258            zh_s(ji,jk) = h_s_1d(ji) * r1_nlay_s 
     259         END DO 
     260      END DO 
     261 
    246262      ! If heat still available (zq_su > 0), then melt more snow 
    247263      zdeltah(1:npti,:) = 0._wp 
     
    250266         DO ji = 1, npti 
    251267            ! thickness change 
    252             rswitch          = 1._wp - MAX( 0._wp, SIGN( 1._wp, - h_s_1d(ji) ) )  
     268            rswitch          = 1._wp - MAX( 0._wp, SIGN( 1._wp, - zh_s(ji,jk) ) )  
    253269            rswitch          = rswitch * ( MAX( 0._wp, SIGN( 1._wp, e_s_1d(ji,jk) - epsi20 ) ) )  
    254270            zdeltah  (ji,jk) = - rswitch * zq_su(ji) / MAX( e_s_1d(ji,jk), epsi20 ) 
    255             zdeltah  (ji,jk) = MAX( zdeltah(ji,jk) , - h_s_1d(ji) ) ! bound melting 
     271            zdeltah  (ji,jk) = MAX( zdeltah(ji,jk) , - zh_s(ji,jk) ) ! bound melting 
    256272            zdh_s_mel(ji)    = zdh_s_mel(ji) + zdeltah(ji,jk)     
    257273            ! heat used to melt snow(W.m-2, >0) 
     
    263279            zq_su (ji) = MAX( 0._wp , zq_su (ji) + zdeltah(ji,jk) * e_s_1d(ji,jk) ) 
    264280            h_s_1d(ji) = MAX( 0._wp , h_s_1d(ji) + zdeltah(ji,jk) ) 
     281            zh_s(ji,jk)= MAX( 0._wp , zh_s(ji,jk) + zdeltah(ji,jk) ) 
    265282         END DO 
    266283      END DO 
  • branches/2017/dev_merge_2017/NEMOGCM/NEMO/LIM_SRC_3/icevar.F90

    r9121 r9271  
    474474         END DO 
    475475         ! 
     476         DO jk = 1, nlay_s 
     477            DO jj = 1 , jpj 
     478               DO ji = 1 , jpi 
     479                  ! update exchanges with ocean 
     480                  hfx_res(ji,jj)   = hfx_res(ji,jj) - (1._wp - zswitch(ji,jj) ) * e_s(ji,jj,jk,jl) * r1_rdtice ! W.m-2 <0 
     481                  e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) * zswitch(ji,jj) 
     482                  t_s(ji,jj,jk,jl) = t_s(ji,jj,jk,jl) * zswitch(ji,jj) + rt0 * ( 1._wp - zswitch(ji,jj) ) 
     483               END DO 
     484            END DO 
     485         END DO 
     486         ! 
    476487         DO jj = 1 , jpj 
    477488            DO ji = 1 , jpi 
     
    480491               wfx_res(ji,jj)  = wfx_res(ji,jj) + (1._wp - zswitch(ji,jj) ) * v_i (ji,jj,jl)   * rhoic * r1_rdtice 
    481492               wfx_res(ji,jj)  = wfx_res(ji,jj) + (1._wp - zswitch(ji,jj) ) * v_s (ji,jj,jl)   * rhosn * r1_rdtice 
    482                hfx_res(ji,jj)  = hfx_res(ji,jj) - (1._wp - zswitch(ji,jj) ) * e_s (ji,jj,1,jl) * r1_rdtice ! W.m-2 <0 
    483                !----------------------------------------------------------------- 
    484                ! Zap snow energy  
    485                !----------------------------------------------------------------- 
    486                t_s(ji,jj,1,jl) = t_s(ji,jj,1,jl) * zswitch(ji,jj) + rt0 * ( 1._wp - zswitch(ji,jj) ) 
    487                e_s(ji,jj,1,jl) = e_s(ji,jj,1,jl) * zswitch(ji,jj) 
    488493               ! 
    489494               !----------------------------------------------------------------- 
Note: See TracChangeset for help on using the changeset viewer.