Changeset 10014 for NEMO/branches/UKMO/dev_r9950_old_tidal_mixing/src

Timestamp:

2018-07-30T12:27:40+02:00 (6 years ago)

Author:

davestorkey

Message:

UKMO dev_r9950_old_tidal_mixing branch: update to be relative to rev 10011 of NEMO4_beta_mirror branch.

Location:

NEMO/branches/UKMO/dev_r9950_old_tidal_mixing/src

Files:

• ICE/icealb.F90 (modified) (1 diff)
• OCE/DYN/dynzad.F90 (modified) (1 diff)
• OCE/ICB/icbdia.F90 (modified) (8 diffs)
• OCE/ICB/icbthm.F90 (modified) (4 diffs)
• OCE/IOM/iom.F90 (modified) (10 diffs)
• OCE/LBC/mpp_nfd_generic.h90 (modified) (2 diffs)
• OCE/SBC/sbcwave.F90 (modified) (1 diff)

NEMO/branches/UKMO/dev_r9950_old_tidal_mixing/src/ICE/icealb.F90

@@ -110 +110 @@
       IF( ln_timing )   CALL timing_start('icealb')
       !
-      z1_href_pnd = 0.05
+      z1_href_pnd = 1. / 0.05
       z1_c1 = 1. / ( LOG(1.5) - LOG(0.05) ) 
       z1_c2 = 1. / 0.05

NEMO/branches/UKMO/dev_r9950_old_tidal_mixing/src/OCE/DYN/dynzad.F90

@@ -89 +89 @@
       !
       ! Surface and bottom advective fluxes set to zero
-      IF( ln_isfcav ) THEN
-         DO jj = 2, jpjm1
-            DO ji = fs_2, fs_jpim1           ! vector opt.
-               zwuw(ji,jj, 1:miku(ji,jj) ) = 0._wp
-               zwvw(ji,jj, 1:mikv(ji,jj) ) = 0._wp
-               zwuw(ji,jj,jpk) = 0._wp
-               zwvw(ji,jj,jpk) = 0._wp
-            END DO
-         END DO
-      ELSE
-         DO jj = 2, jpjm1        
-            DO ji = fs_2, fs_jpim1           ! vector opt.
-               zwuw(ji,jj, 1 ) = 0._wp
-               zwvw(ji,jj, 1 ) = 0._wp
-               zwuw(ji,jj,jpk) = 0._wp
-               zwvw(ji,jj,jpk) = 0._wp
-            END DO  
-         END DO
-      END IF
-
+      DO jj = 2, jpjm1        
+         DO ji = fs_2, fs_jpim1           ! vector opt.
+            zwuw(ji,jj, 1 ) = 0._wp
+            zwvw(ji,jj, 1 ) = 0._wp
+            zwuw(ji,jj,jpk) = 0._wp
+            zwvw(ji,jj,jpk) = 0._wp
+         END DO  
+      END DO
+      !
       DO jk = 1, jpkm1              ! Vertical momentum advection at u- and v-points
          DO jj = 2, jpjm1

NEMO/branches/UKMO/dev_r9950_old_tidal_mixing/src/OCE/ICB/icbdia.F90

@@ -51 +51 @@
 
    REAL(wp), DIMENSION(:,:)  , ALLOCATABLE, PUBLIC  ::   berg_melt       ! Melting+erosion rate of icebergs     [kg/s/m2]
+   REAL(wp), DIMENSION(:,:)  , ALLOCATABLE, PUBLIC  ::   berg_melt_hcflx ! Heat flux to ocean due to heat content of melting icebergs [J/s/m2]
+   REAL(wp), DIMENSION(:,:)  , ALLOCATABLE, PUBLIC  ::   berg_melt_qlat  ! Heat flux to ocean due to latent heat of melting icebergs [J/s/m2]
    REAL(wp), DIMENSION(:,:)  , ALLOCATABLE, PUBLIC  ::   buoy_melt       ! Buoyancy component of melting rate   [kg/s/m2]
    REAL(wp), DIMENSION(:,:)  , ALLOCATABLE, PUBLIC  ::   eros_melt       ! Erosion component of melting rate    [kg/s/m2]
@@ -101 +103 @@
 
       ALLOCATE( berg_melt    (jpi,jpj)   )           ;   berg_melt   (:,:)   = 0._wp
+      ALLOCATE( berg_melt_hcflx(jpi,jpj) )           ;   berg_melt_hcflx(:,:)   = 0._wp
+      ALLOCATE( berg_melt_qlat(jpi,jpj)  )           ;   berg_melt_qlat(:,:)   = 0._wp
       ALLOCATE( buoy_melt    (jpi,jpj)   )           ;   buoy_melt   (:,:)   = 0._wp
       ALLOCATE( eros_melt    (jpi,jpj)   )           ;   eros_melt   (:,:)   = 0._wp
@@ -364 +368 @@
       IF( .NOT.ln_bergdia )   RETURN
       berg_melt   (:,:)   = 0._wp
+      berg_melt_hcflx(:,:)   = 0._wp
+      berg_melt_qlat(:,:)   = 0._wp
       buoy_melt   (:,:)   = 0._wp
       eros_melt   (:,:)   = 0._wp
@@ -384 +390 @@
       !
       CALL iom_put( "berg_melt"        , berg_melt   (:,:)   )   ! Melt rate of icebergs                     [kg/m2/s]
+      !! NB. The berg_melt_hcflx field is currently always zero - see comment in icbthm.F90
+      CALL iom_put( "berg_melt_hcflx"  , berg_melt_hcflx(:,:))   ! Heat flux to ocean due to heat content of melting icebergs [J/m2/s]
+      CALL iom_put( "berg_melt_qlat"   , berg_melt_qlat(:,:) )   ! Heat flux to ocean due to latent heat of melting icebergs [J/m2/s]
       CALL iom_put( "berg_buoy_melt"   , buoy_melt   (:,:)   )   ! Buoyancy component of iceberg melt rate   [kg/m2/s]
       CALL iom_put( "berg_eros_melt"   , eros_melt   (:,:)   )   ! Erosion component of iceberg melt rate    [kg/m2/s]
@@ -471 +480 @@
 
 
-   SUBROUTINE icb_dia_melt(ki, kj, pmnew, pheat, pmass_scale,     &
+   SUBROUTINE icb_dia_melt(ki, kj, pmnew, pheat_hcflux, pheat_latent, pmass_scale,     &
       &                    pdM, pdMbitsE, pdMbitsM, pdMb, pdMe,   &
       &                    pdMv, pz1_dt_e1e2 )
@@ -477 +486 @@
       !!----------------------------------------------------------------------
       INTEGER , INTENT(in) ::   ki, kj
-      REAL(wp), INTENT(in) ::   pmnew, pheat, pmass_scale
+      REAL(wp), INTENT(in) ::   pmnew, pheat_hcflux, pheat_latent, pmass_scale
       REAL(wp), INTENT(in) ::   pdM, pdMbitsE, pdMbitsM, pdMb, pdMe, pdMv, pz1_dt_e1e2
       !!----------------------------------------------------------------------
@@ -484 +493 @@
       !
       berg_melt (ki,kj) = berg_melt (ki,kj) + pdM      * pz1_dt_e1e2   ! kg/m2/s
+      berg_melt_hcflx (ki,kj) = berg_melt_hcflx (ki,kj) + pheat_hcflux * pz1_dt_e1e2   ! J/m2/s
+      berg_melt_qlat (ki,kj) = berg_melt_qlat (ki,kj) + pheat_latent * pz1_dt_e1e2   ! J/m2/s
       bits_src  (ki,kj) = bits_src  (ki,kj) + pdMbitsE * pz1_dt_e1e2   ! mass flux into bergy bitskg/m2/s
       bits_melt (ki,kj) = bits_melt (ki,kj) + pdMbitsM * pz1_dt_e1e2   ! melt rate of bergy bits kg/m2/s
@@ -489 +500 @@
       eros_melt (ki,kj) = eros_melt (ki,kj) + pdMe     * pz1_dt_e1e2   ! erosion rate kg/m2/s
       conv_melt (ki,kj) = conv_melt (ki,kj) + pdMv     * pz1_dt_e1e2   ! kg/m2/s
-      heat_to_ocean_net = heat_to_ocean_net + pheat * pmass_scale * berg_dt         ! J
+      heat_to_ocean_net = heat_to_ocean_net + (pheat_hcflux + pheat_latent) * pmass_scale * berg_dt         ! J
       IF( pmnew <= 0._wp ) nbergs_melted = nbergs_melted + 1                        ! Delete the berg if completely melted
       !

NEMO/branches/UKMO/dev_r9950_old_tidal_mixing/src/OCE/ICB/icbthm.F90

@@ -50 +50 @@
       REAL(wp) ::   zM, zT, zW, zL, zSST, zVol, zLn, zWn, zTn, znVol, zIC, zDn
       REAL(wp) ::   zMv, zMe, zMb, zmelt, zdvo, zdva, zdM, zSs, zdMe, zdMb, zdMv
-      REAL(wp) ::   zMnew, zMnew1, zMnew2, zheat, z1_12
+      REAL(wp) ::   zMnew, zMnew1, zMnew2, zheat_hcflux, zheat_latent, z1_12
       REAL(wp) ::   zMbits, znMbits, zdMbitsE, zdMbitsM, zLbits, zAbits, zMbb
       REAL(wp) ::   zxi, zyj, zff, z1_rday, z1_e1e2, zdt, z1_dt, z1_dt_e1e2
@@ -68 +68 @@
       !
       berg_grid%floating_melt(:,:) = 0._wp
+      ! calving_hflx re-used here as temporary workspace for the heat flux associated with melting
       berg_grid%calving_hflx(:,:)  = 0._wp
       !
@@ -166 +167 @@
             zmelt    = ( zdM - ( zdMbitsE - zdMbitsM ) ) * z1_dt   ! kg/s
             berg_grid%floating_melt(ii,ij) = berg_grid%floating_melt(ii,ij) + zmelt    * z1_e1e2    ! kg/m2/s
-            zheat = zmelt * pt%heat_density              ! kg/s x J/kg = J/s
-            berg_grid%calving_hflx (ii,ij) = berg_grid%calving_hflx (ii,ij) + zheat    * z1_e1e2    ! W/m2
-            CALL icb_dia_melt( ii, ij, zMnew, zheat, this%mass_scaling,       &
+            !! NB. The src_calving_hflx field is currently hardwired to zero in icb_stp, which means that the
+            !!     heat density of the icebergs is zero and the heat content flux to the ocean from iceberg
+            !!     melting is always zero. Leaving the term in the code until such a time as this is fixed. DS.
+            zheat_hcflux = zmelt * pt%heat_density       ! heat content flux : kg/s x J/kg = J/s
+            zheat_latent = - zmelt * rLfus               ! latent heat flux:  kg/s x J/kg = J/s
+            berg_grid%calving_hflx (ii,ij) = berg_grid%calving_hflx (ii,ij) + ( zheat_hcflux + zheat_latent ) * z1_e1e2    ! W/m2
+            CALL icb_dia_melt( ii, ij, zMnew, zheat_hcflux, zheat_latent, this%mass_scaling,       &
                &                       zdM, zdMbitsE, zdMbitsM, zdMb, zdMe,   &
                &                       zdMv, z1_dt_e1e2 )
@@ -214 +219 @@
       IF(.NOT. ln_passive_mode ) THEN
          emp (:,:) = emp (:,:) - berg_grid%floating_melt(:,:)
-!!       qns (:,:) = qns (:,:) + berg_grid%calving_hflx (:,:)  !!gm heat flux not yet properly coded ==>> need it, SOLVE that!
+         qns (:,:) = qns (:,:) + berg_grid%calving_hflx (:,:)  
       ENDIF
       !

NEMO/branches/UKMO/dev_r9950_old_tidal_mixing/src/OCE/IOM/iom.F90

@@ -150 +150 @@
       !
       IF( TRIM(cdname) == TRIM(cxios_context) ) THEN  
-         CALL set_grid( "T", glamt, gphit, .FALSE. ) 
-         CALL set_grid( "U", glamu, gphiu, .FALSE. )
-         CALL set_grid( "V", glamv, gphiv, .FALSE. )
-         CALL set_grid( "W", glamt, gphit, .FALSE. )
+         CALL set_grid( "T", glamt, gphit, .FALSE., .FALSE. ) 
+         CALL set_grid( "U", glamu, gphiu, .FALSE., .FALSE. )
+         CALL set_grid( "V", glamv, gphiv, .FALSE., .FALSE. )
+         CALL set_grid( "W", glamt, gphit, .FALSE., .FALSE. )
          CALL set_grid_znl( gphit )
          !
@@ -171 +171 @@
          CALL dom_grid_crs   ! Save the parent grid information  & Switch to coarse grid domain
          !
-         CALL set_grid( "T", glamt_crs, gphit_crs, .FALSE. ) 
-         CALL set_grid( "U", glamu_crs, gphiu_crs, .FALSE. ) 
-         CALL set_grid( "V", glamv_crs, gphiv_crs, .FALSE. ) 
-         CALL set_grid( "W", glamt_crs, gphit_crs, .FALSE. ) 
+         CALL set_grid( "T", glamt_crs, gphit_crs, .FALSE., .FALSE. ) 
+         CALL set_grid( "U", glamu_crs, gphiu_crs, .FALSE., .FALSE. ) 
+         CALL set_grid( "V", glamv_crs, gphiv_crs, .FALSE., .FALSE. ) 
+         CALL set_grid( "W", glamt_crs, gphit_crs, .FALSE., .FALSE. ) 
          CALL set_grid_znl( gphit_crs )
           !
@@ -227 +227 @@
       IF( TRIM(cdname) == TRIM(crxios_context) ) THEN
 !set names of the fields in restart file IF using XIOS to read data
-          CALL iom_set_rst_context()
+          CALL iom_set_rst_context(.TRUE.)
           CALL iom_set_rst_vars(rst_rfields)
 !set which fields are to be read from restart file
@@ -233 +233 @@
       ELSE IF( TRIM(cdname) == TRIM(cwxios_context) ) THEN
 !set names of the fields in restart file IF using XIOS to write data
-          CALL iom_set_rst_context()
+          CALL iom_set_rst_context(.FALSE.)
           CALL iom_set_rst_vars(rst_wfields)
 !set which fields are to be written to a restart file
@@ -571 +571 @@
    END SUBROUTINE iom_set_rstw_active
 
-   SUBROUTINE iom_set_rst_context( ) 
+   SUBROUTINE iom_set_rst_context(ld_rstr) 
      !!---------------------------------------------------------------------
       !!                   ***  SUBROUTINE  iom_set_rst_context  ***
@@ -579 +579 @@
       !!               
       !!---------------------------------------------------------------------
+   LOGICAL, INTENT(IN)               :: ld_rstr
+!ld_rstr is true for restart context. There is no need to define grid for 
+!restart read, because it's read from file
 #if defined key_iomput
    TYPE(xios_domaingroup)            :: domaingroup_hdl 
@@ -589 +592 @@
      CALL xios_get_handle("domain_definition",domaingroup_hdl) 
      CALL xios_add_child(domaingroup_hdl, domain_hdl, "grid_N") 
-     CALL set_grid("N", glamt, gphit, .TRUE.) 
+     CALL set_grid("N", glamt, gphit, .TRUE., ld_rstr) 
  
      CALL xios_get_handle("axis_definition",axisgroup_hdl) 
@@ -1935 +1938 @@
 
 
-   SUBROUTINE set_grid( cdgrd, plon, plat, ldxios )
+   SUBROUTINE set_grid( cdgrd, plon, plat, ldxios, ldrxios )
       !!----------------------------------------------------------------------
       !!                     ***  ROUTINE set_grid  ***
@@ -1947 +1950 @@
       INTEGER  :: ni, nj
       REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zmask
-      LOGICAL, INTENT(IN) :: ldxios
+      LOGICAL, INTENT(IN) :: ldxios, ldrxios
       !!----------------------------------------------------------------------
       !
@@ -1955 +1958 @@
       CALL iom_set_domain_attr("grid_"//cdgrd, ni_glo=jpiglo, nj_glo=jpjglo, ibegin=nimpp+nldi-2, jbegin=njmpp+nldj-2, ni=ni, nj=nj)
       CALL iom_set_domain_attr("grid_"//cdgrd, data_dim=2, data_ibegin = 1-nldi, data_ni = jpi, data_jbegin = 1-nldj, data_nj = jpj)
-      CALL iom_set_domain_attr("grid_"//cdgrd, lonvalue = RESHAPE(plon(nldi:nlei, nldj:nlej),(/ ni*nj /)),   &
+!don't define lon and lat for restart reading context. 
+      IF ( .NOT.ldrxios ) &
+         CALL iom_set_domain_attr("grid_"//cdgrd, lonvalue = RESHAPE(plon(nldi:nlei, nldj:nlej),(/ ni*nj /)),   &
          &                                     latvalue = RESHAPE(plat(nldi:nlei, nldj:nlej),(/ ni*nj /)))  
       !

NEMO/branches/UKMO/dev_r9950_old_tidal_mixing/src/OCE/LBC/mpp_nfd_generic.h90

@@ -95 +95 @@
          ALLOCATE( ztabl(jpimax   ,4,ipk,ipl,ipf) , ztabr(jpimax*jpmaxngh,4,ipk,ipl,ipf) ) 
          !
-         ztabr(:,:,:,:,:) = 0._wp
-         ztabl(:,:,:,:,:) = 0._wp
+         ! when some processors of the north fold are suppressed, 
+         ! values of ztab* arrays corresponding to these suppressed domain won't be defined 
+         ! and we need a default definition to 0.
+         ! a better test should be: a testing if "suppressed land-processors" belongs to the north-pole folding
+         IF ( jpni*jpnj /= jpnij ) THEN
+            ztabr(:,:,:,:,:) = 0._wp
+            ztabl(:,:,:,:,:) = 0._wp
+         END IF
          !
          DO jf = 1, ipf
@@ -183 +189 @@
          ALLOCATE( znorthgloio(jpimax,4,ipk,ipl,ipf,jpni) )
          !
-         ztab(:,:,:,:,:)=0._wp
+         ! when some processors of the north fold are suppressed,
+         ! values of ztab* arrays corresponding to these suppressed domain won't be defined
+         ! and we need a default definition to 0.
+         ! a better test should be: a testing if "suppressed land-processors" belongs to the north-pole folding
+         IF ( jpni*jpnj /= jpnij ) ztab(:,:,:,:,:)=0._wp
          !
          CALL MPI_ALLGATHER( znorthloc  , itaille, MPI_DOUBLE_PRECISION,                &

NEMO/branches/UKMO/dev_r9950_old_tidal_mixing/src/OCE/SBC/sbcwave.F90

@@ -131 +131 @@
          END DO
       ELSE IF( ll_st_peakfr ) THEN    ! peak wave number calculated from the peak frequency received by the wave model
+         DO jj = 1, jpj
+            DO ji = 1, jpi
+               zk_t(ji,jj) = ( 2.0_wp * rpi * wfreq(ji,jj) ) * ( 2.0_wp * rpi * wfreq(ji,jj) ) / grav
+            END DO
+         END DO
          DO jj = 1, jpjm1
             DO ji = 1, jpim1
-               zk_u(ji,jj) = 0.5_wp * ( wfreq(ji,jj)*wfreq(ji,jj) + wfreq(ji+1,jj)*wfreq(ji+1,jj) ) / grav
-               zk_v(ji,jj) = 0.5_wp * ( wfreq(ji,jj)*wfreq(ji,jj) + wfreq(ji,jj+1)*wfreq(ji,jj+1) ) / grav
+               zk_u(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji+1,jj) )
+               zk_v(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji,jj+1) )
                !
                zu0_sd(ji,jj) = 0.5_wp * ( ut0sd(ji,jj) + ut0sd(ji+1,jj) )