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 9495 – NEMO

Changeset 9495


Ignore:
Timestamp:
2018-04-23T16:34:55+02:00 (6 years ago)
Author:
davestorkey
Message:

branches/UKMO/dev_merge_2017_restart_datestamp_GO6_mixing : revert previous change so I can remove SVN keywords first.

Location:
branches/UKMO/dev_merge_2017_restart_datestamp_GO6_mixing/NEMOGCM
Files:
8 edited

Legend:

Unmodified
Added
Removed
  • branches/UKMO/dev_merge_2017_restart_datestamp_GO6_mixing/NEMOGCM/CONFIG/SHARED/field_def_nemo-opa.xml

    r9493 r9495  
    6262         <field id="mldr10_1max"  long_name="Max of Mixed Layer Depth (dsigma = 0.01 wrt 10m)"   field_ref="mldr10_1"   operation="maximum"                                                                          /> 
    6363         <field id="mldr10_1min"  long_name="Min of Mixed Layer Depth (dsigma = 0.01 wrt 10m)"   field_ref="mldr10_1"   operation="minimum"                                                                          /> 
    64          <field id="mldzint_1"    long_name="Mixed Layer Depth interpolated"                     standard_name="ocean_mixed_layer_thickness"                                                       unit="m"          /> 
    65          <field id="mldzint_2"    long_name="Mixed Layer Depth interpolated"                     standard_name="ocean_mixed_layer_thickness"                                                       unit="m"          /> 
    66          <field id="mldzint_3"    long_name="Mixed Layer Depth interpolated"                     standard_name="ocean_mixed_layer_thickness"                                                       unit="m"          /> 
    67          <field id="mldzint_4"    long_name="Mixed Layer Depth interpolated"                     standard_name="ocean_mixed_layer_thickness"                                                       unit="m"          /> 
    68          <field id="mldzint_5"    long_name="Mixed Layer Depth interpolated"                     standard_name="ocean_mixed_layer_thickness"                                                       unit="m"          /> 
    69          <field id="mldhtc_1"     long_name="Mixed Layer Depth integrated heat content"          standard_name="integral_of_sea_water_potential_temperature_wrt_depth_expressed_as_heat_content"   unit="J/m2"       /> 
    70          <field id="mldhtc_2"     long_name="Mixed Layer Depth integrated heat content"          standard_name="integral_of_sea_water_potential_temperature_wrt_depth_expressed_as_heat_content"   unit="J/m2"       /> 
    71          <field id="mldhtc_3"     long_name="Mixed Layer Depth integrated heat content"          standard_name="integral_of_sea_water_potential_temperature_wrt_depth_expressed_as_heat_content"   unit="J/m2"       /> 
    72          <field id="mldhtc_4"     long_name="Mixed Layer Depth integrated heat content"          standard_name="integral_of_sea_water_potential_temperature_wrt_depth_expressed_as_heat_content"   unit="J/m2"       /> 
    73          <field id="mldhtc_5"     long_name="Mixed Layer Depth integrated heat content"          standard_name="integral_of_sea_water_potential_temperature_wrt_depth_expressed_as_heat_content"   unit="J/m2"       /> 
    7464         <field id="heatc"        long_name="Heat content vertically integrated"                 standard_name="integral_of_sea_water_potential_temperature_wrt_depth_expressed_as_heat_content"   unit="J/m2"       /> 
    7565         <field id="saltc"        long_name="Salt content vertically integrated"                                                                                                                   unit="1e-3*kg/m2" /> 
  • branches/UKMO/dev_merge_2017_restart_datestamp_GO6_mixing/NEMOGCM/NEMO/LIM_SRC_3/icerst.F90

    r9493 r9495  
    2424   USE in_out_manager ! I/O manager 
    2525   USE iom            ! I/O manager library 
    26    USE ioipsl, ONLY : ju2ymds    ! for calendar 
    2726   USE lib_mpp        ! MPP library 
    2827   USE lib_fortran    ! fortran utilities (glob_sum + no signed zero) 
     
    5049      INTEGER, INTENT(in) ::   kt       ! number of iteration 
    5150      ! 
    52       INTEGER             ::   iyear, imonth, iday 
    53       REAL (wp)           ::   zsec 
    54       REAL (wp)           ::   zfjulday 
    5551      CHARACTER(len=20)   ::   clkt     ! ocean time-step define as a character 
    5652      CHARACTER(len=50)   ::   clname   ! ice output restart file name 
     
    6763         IF( nitrst <= nitend .AND. nitrst > 0 ) THEN 
    6864            ! beware of the format used to write kt (default is i8.8, that should be large enough...) 
    69             IF ( ln_rstdate ) THEN 
    70                zfjulday = fjulday + (2*nn_fsbc+1)*rdt / rday 
    71                IF( ABS(zfjulday - REAL(NINT(zfjulday),wp)) < 0.1 / rday )   zfjulday = REAL(NINT(zfjulday),wp)   ! avoid truncation error 
    72                CALL ju2ymds( zfjulday, iyear, imonth, iday, zsec )            
    73                WRITE(clkt, '(i4.4,2i2.2)') iyear, imonth, iday 
    74             ELSE 
    75                IF( nitrst > 99999999 ) THEN   ;   WRITE(clkt, *       ) nitrst 
    76                ELSE                           ;   WRITE(clkt, '(i8.8)') nitrst 
    77                ENDIF 
     65            IF( nitrst > 99999999 ) THEN   ;   WRITE(clkt, *       ) nitrst 
     66            ELSE                           ;   WRITE(clkt, '(i8.8)') nitrst 
    7867            ENDIF 
    7968            ! create the file 
  • branches/UKMO/dev_merge_2017_restart_datestamp_GO6_mixing/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90

    r9493 r9495  
    292292         &             nn_it000, nn_itend , nn_date0    , nn_time0     , nn_leapy  , nn_istate ,     & 
    293293         &             nn_stock, nn_write , ln_mskland  , ln_clobber   , nn_chunksz, nn_euler  ,     & 
    294          &             ln_cfmeta, ln_iscpl, ln_xios_read, nn_wxios, ln_rstdate 
     294         &             ln_cfmeta, ln_iscpl, ln_xios_read, nn_wxios 
    295295      NAMELIST/namdom/ ln_linssh, rn_isfhmin, rn_rdt, rn_atfp, ln_crs, ln_meshmask 
    296296#if defined key_netcdf4 
     
    338338         WRITE(numout,*) '      frequency of output file        nn_write        = ', nn_write 
    339339         WRITE(numout,*) '      mask land points                ln_mskland      = ', ln_mskland 
    340          WRITE(numout,*) '      date-stamp restart files        ln_rstdate = ', ln_rstdate 
    341340         WRITE(numout,*) '      additional CF standard metadata ln_cfmeta       = ', ln_cfmeta 
    342341         WRITE(numout,*) '      overwrite an existing file      ln_clobber      = ', ln_clobber 
  • branches/UKMO/dev_merge_2017_restart_datestamp_GO6_mixing/NEMOGCM/NEMO/OPA_SRC/ICB/icbrst.F90

    r9493 r9495  
    2525   USE netcdf         ! netcdf routines for IO 
    2626   USE iom 
    27    USE ioipsl, ONLY : ju2ymds    ! for calendar 
    2827   USE icb_oce        ! define iceberg arrays 
    2928   USE icbutl         ! iceberg utility routines 
     
    158157      INTEGER ::   jn   ! dummy loop index 
    159158      INTEGER ::   ix_dim, iy_dim, ik_dim, in_dim 
    160       INTEGER             ::   iyear, imonth, iday 
    161       REAL (wp)           ::   zsec 
    162       REAL (wp)           ::   zfjulday 
    163       CHARACTER(len=256)  :: cl_path 
    164       CHARACTER(len=256)  :: cl_filename 
    165       CHARACTER(LEN=20)   ::   clkt     ! ocean time-step deine as a character 
     159      CHARACTER(len=256)     :: cl_path 
     160      CHARACTER(len=256)     :: cl_filename 
    166161      TYPE(iceberg), POINTER :: this 
    167162      TYPE(point)  , POINTER :: pt 
     
    171166      cl_path = TRIM(cn_ocerst_outdir) 
    172167      IF( cl_path(LEN_TRIM(cl_path):) /= '/' ) cl_path = TRIM(cl_path) // '/' 
    173       IF ( ln_rstdate ) THEN 
    174          zfjulday = fjulday + rdt / rday 
    175          IF( ABS(zfjulday - REAL(NINT(zfjulday),wp)) < 0.1 / rday )   zfjulday = REAL(NINT(zfjulday),wp)   ! avoid truncation error 
    176          CALL ju2ymds( zfjulday, iyear, imonth, iday, zsec )            
    177          WRITE(clkt, '(i4.4,2i2.2)') iyear, imonth, iday 
     168      IF( lk_mpp ) THEN 
     169         WRITE(cl_filename,'(A,"_icebergs_",I8.8,"_restart_",I4.4,".nc")') TRIM(cexper), kt, narea-1 
    178170      ELSE 
    179          IF( kt > 999999999 ) THEN   ;   WRITE(clkt, *       ) kt 
    180          ELSE                        ;   WRITE(clkt, '(i8.8)') kt 
    181          ENDIF 
    182       ENDIF 
    183       IF( lk_mpp ) THEN 
    184          WRITE(cl_filename,'(A,"_icebergs_",A,"_restart_",I4.4,".nc")') TRIM(cexper), TRIM(ADJUSTL(clkt)), narea-1 
    185       ELSE 
    186          WRITE(cl_filename,'(A,"_icebergs_",A,"_restart.nc")') TRIM(cexper), TRIM(ADJUSTL(clkt)) 
     171         WRITE(cl_filename,'(A,"_icebergs_",I8.8,"_restart.nc")') TRIM(cexper), kt 
    187172      ENDIF 
    188173      IF (nn_verbose_level >= 0) WRITE(numout,'(2a)') 'icebergs, write_restart: creating ',TRIM(cl_path)//TRIM(cl_filename) 
  • branches/UKMO/dev_merge_2017_restart_datestamp_GO6_mixing/NEMOGCM/NEMO/OPA_SRC/IOM/in_out_manager.F90

    r9493 r9495  
    4040   INTEGER, DIMENSION(10) :: nn_stocklist  !: restart dump times 
    4141   LOGICAL       ::   ln_mskland       !: mask land points in NetCDF outputs (costly: + ~15%) 
    42    LOGICAL       ::   ln_rstdate       !: T=> stamp output restart files with date instead of timestep 
    4342   LOGICAL       ::   ln_cfmeta        !: output additional data to netCDF files required for compliance with the CF metadata standard 
    4443   LOGICAL       ::   ln_clobber       !: clobber (overwrite) an existing file 
  • branches/UKMO/dev_merge_2017_restart_datestamp_GO6_mixing/NEMOGCM/NEMO/OPA_SRC/IOM/restart.F90

    r9493 r9495  
    2727   USE in_out_manager  ! I/O manager 
    2828   USE iom             ! I/O module 
    29    USE ioipsl, ONLY : ju2ymds    ! for calendar 
    3029   USE diurnal_bulk 
    3130 
     
    5958      INTEGER, INTENT(in) ::   kt     ! ocean time-step 
    6059      !! 
    61       INTEGER             ::   iyear, imonth, iday 
    62       REAL (wp)           ::   zsec 
    63       REAL (wp)           ::   zfjulday 
    6460      CHARACTER(LEN=20)   ::   clkt     ! ocean time-step deine as a character 
    6561      CHARACTER(LEN=50)   ::   clname   ! ocean output restart file name 
     
    9086         IF( nitrst <= nitend .AND. nitrst > 0 ) THEN  
    9187            ! beware of the format used to write kt (default is i8.8, that should be large enough...) 
    92             IF ( ln_rstdate ) THEN 
    93                zfjulday = fjulday + rdt / rday 
    94                IF( ABS(zfjulday - REAL(NINT(zfjulday),wp)) < 0.1 / rday )   zfjulday = REAL(NINT(zfjulday),wp)   ! avoid truncation error 
    95                CALL ju2ymds( zfjulday, iyear, imonth, iday, zsec )            
    96                WRITE(clkt, '(i4.4,2i2.2)') iyear, imonth, iday 
    97             ELSE 
    98                IF( nitrst > 999999999 ) THEN   ;   WRITE(clkt, *       ) nitrst 
    99                ELSE                            ;   WRITE(clkt, '(i8.8)') nitrst 
    100                ENDIF 
     88            IF( nitrst > 999999999 ) THEN   ;   WRITE(clkt, *       ) nitrst 
     89            ELSE                            ;   WRITE(clkt, '(i8.8)') nitrst 
    10190            ENDIF 
    10291            ! create the file 
  • branches/UKMO/dev_merge_2017_restart_datestamp_GO6_mixing/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfmxl.F90

    r9493 r9495  
    1515   USE trc_oce  , ONLY: l_offline         ! ocean space and time domain variables 
    1616   USE zdf_oce        ! ocean vertical physics 
    17    USE eosbn2         ! for zdf_mxl_zint 
    1817   ! 
    1918   USE in_out_manager ! I/O manager 
     
    3231   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hmlp    !: mixed layer depth  (rho=rho0+zdcrit) [m]   (used by LDF) 
    3332   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hmlpt   !: depth of the last T-point inside the mixed layer [m] (used by LDF) 
    34    REAL(wp), PUBLIC, ALLOCATABLE,       DIMENSION(:,:) ::   hmld_zint  !: vertically-interpolated mixed layer depth   [m]  
    35    REAL(wp), PUBLIC, ALLOCATABLE,       DIMENSION(:,:) ::   htc_mld    ! Heat content of hmld_zint 
    36    LOGICAL, PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)    :: ll_found   ! Is T_b to be found by interpolation ?  
    37    LOGICAL, PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  :: ll_belowml ! Flag points below mixed layer when ll_found=F 
    3833 
    3934   REAL(wp), PUBLIC ::   rho_c = 0.01_wp    !: density criterion for mixed layer depth 
    4035   REAL(wp)         ::   avt_c = 5.e-4_wp   ! Kz criterion for the turbocline depth 
    41  
    42    TYPE, PUBLIC :: MXL_ZINT   !: Structure for MLD defs 
    43       INTEGER   :: mld_type   ! mixed layer type      
    44       REAL(wp)  :: zref       ! depth of initial T_ref 
    45       REAL(wp)  :: dT_crit    ! Critical temp diff 
    46       REAL(wp)  :: iso_frac   ! Fraction of rn_dT_crit used 
    47    END TYPE MXL_ZINT 
    4836 
    4937   !!---------------------------------------------------------------------- 
     
    6048      zdf_mxl_alloc = 0      ! set to zero if no array to be allocated 
    6149      IF( .NOT. ALLOCATED( nmln ) ) THEN 
    62          ALLOCATE( nmln(jpi,jpj), hmld(jpi,jpj), hmlp(jpi,jpj), hmlpt(jpi,jpj), hmld_zint(jpi,jpj),     & 
    63         &          htc_mld(jpi,jpj),                                                                    & 
    64         &          ll_found(jpi,jpj), ll_belowml(jpi,jpj,jpk), STAT= zdf_mxl_alloc ) 
     50         ALLOCATE( nmln(jpi,jpj), hmld(jpi,jpj), hmlp(jpi,jpj), hmlpt(jpi,jpj), STAT= zdf_mxl_alloc ) 
    6551         ! 
    6652         IF( lk_mpp             )   CALL mpp_sum ( zdf_mxl_alloc ) 
     
    151137      ENDIF 
    152138      ! 
    153       ! Vertically-interpolated mixed-layer depth diagnostic 
    154       CALL zdf_mxl_zint( kt ) 
    155       ! 
    156139      IF(ln_ctl)   CALL prt_ctl( tab2d_1=REAL(nmln,wp), clinfo1=' nmln : ', tab2d_2=hmlp, clinfo2=' hmlp : ' ) 
    157140      ! 
    158141   END SUBROUTINE zdf_mxl 
    159142 
    160    SUBROUTINE zdf_mxl_zint_mld( sf )  
    161       !!----------------------------------------------------------------------------------  
    162       !!                    ***  ROUTINE zdf_mxl_zint_mld  ***  
    163       !                                                                         
    164       !   Calculate vertically-interpolated mixed layer depth diagnostic.  
    165       !             
    166       !   This routine can calculate the mixed layer depth diagnostic suggested by 
    167       !   Kara et al, 2000, JGR, 105, 16803, but is more general and can calculate 
    168       !   vertically-interpolated mixed-layer depth diagnostics with other parameter 
    169       !   settings set in the namzdf_mldzint namelist.   
    170       !  
    171       !   If mld_type=1 the mixed layer depth is calculated as the depth at which the   
    172       !   density has increased by an amount equivalent to a temperature difference of   
    173       !   0.8C at the surface.  
    174       !  
    175       !   For other values of mld_type the mixed layer is calculated as the depth at   
    176       !   which the temperature differs by 0.8C from the surface temperature.   
    177       !                                                                         
    178       !   David Acreman, Daley Calvert                                       
    179       !  
    180       !!-----------------------------------------------------------------------------------  
    181  
    182       TYPE(MXL_ZINT), INTENT(in)  :: sf 
    183  
    184       ! Diagnostic criteria 
    185       INTEGER   :: nn_mld_type   ! mixed layer type      
    186       REAL(wp)  :: rn_zref       ! depth of initial T_ref 
    187       REAL(wp)  :: rn_dT_crit    ! Critical temp diff 
    188       REAL(wp)  :: rn_iso_frac   ! Fraction of rn_dT_crit used 
    189  
    190       ! Local variables 
    191       REAL(wp), PARAMETER :: zepsilon = 1.e-30          ! local small value 
    192       INTEGER, DIMENSION(jpi,jpj) :: ikmt          ! number of active tracer levels  
    193       INTEGER, DIMENSION(jpi,jpj) :: ik_ref        ! index of reference level  
    194       INTEGER, DIMENSION(jpi,jpj) :: ik_iso        ! index of last uniform temp level  
    195       REAL, DIMENSION(jpi,jpj,jpk)  :: zT            ! Temperature or density  
    196       REAL, DIMENSION(jpi,jpj)    :: ppzdep        ! depth for use in calculating d(rho)  
    197       REAL, DIMENSION(jpi,jpj)    :: zT_ref        ! reference temperature  
    198       REAL    :: zT_b                                   ! base temperature  
    199       REAL, DIMENSION(jpi,jpj,jpk)  :: zdTdz         ! gradient of zT  
    200       REAL, DIMENSION(jpi,jpj,jpk)  :: zmoddT        ! Absolute temperature difference  
    201       REAL    :: zdz                                    ! depth difference  
    202       REAL    :: zdT                                    ! temperature difference  
    203       REAL, DIMENSION(jpi,jpj)    :: zdelta_T      ! difference critereon  
    204       REAL, DIMENSION(jpi,jpj)    :: zRHO1, zRHO2  ! Densities  
    205       INTEGER :: ji, jj, jk                             ! loop counter  
    206  
    207       !!-------------------------------------------------------------------------------------  
    208       !   
    209       ! Unpack structure 
    210       nn_mld_type = sf%mld_type 
    211       rn_zref     = sf%zref 
    212       rn_dT_crit  = sf%dT_crit 
    213       rn_iso_frac = sf%iso_frac 
    214  
    215       ! Set the mixed layer depth criterion at each grid point  
    216       IF( nn_mld_type == 0 ) THEN 
    217          zdelta_T(:,:) = rn_dT_crit 
    218          zT(:,:,:) = rhop(:,:,:) 
    219       ELSE IF( nn_mld_type == 1 ) THEN 
    220          ppzdep(:,:)=0.0  
    221          call eos ( tsn(:,:,1,:), ppzdep(:,:), zRHO1(:,:) )  
    222 ! Use zT temporarily as a copy of tsn with rn_dT_crit added to SST  
    223 ! [assumes number of tracers less than number of vertical levels]  
    224          zT(:,:,1:jpts)=tsn(:,:,1,1:jpts)  
    225          zT(:,:,jp_tem)=zT(:,:,1)+rn_dT_crit  
    226          CALL eos( zT(:,:,1:jpts), ppzdep(:,:), zRHO2(:,:) )  
    227          zdelta_T(:,:) = abs( zRHO1(:,:) - zRHO2(:,:) ) * rau0  
    228          ! RHO from eos (2d version) doesn't calculate north or east halo:  
    229          CALL lbc_lnk( zdelta_T, 'T', 1. )  
    230          zT(:,:,:) = rhop(:,:,:)  
    231       ELSE  
    232          zdelta_T(:,:) = rn_dT_crit                       
    233          zT(:,:,:) = tsn(:,:,:,jp_tem)                            
    234       END IF  
    235  
    236       ! Calculate the gradient of zT and absolute difference for use later  
    237       DO jk = 1 ,jpk-2  
    238          zdTdz(:,:,jk)  =    ( zT(:,:,jk+1) - zT(:,:,jk) ) / e3w_n(:,:,jk+1)  
    239          zmoddT(:,:,jk) = abs( zT(:,:,jk+1) - zT(:,:,jk) )  
    240       END DO  
    241  
    242       ! Find density/temperature at the reference level (Kara et al use 10m).           
    243       ! ik_ref is the index of the box centre immediately above or at the reference level  
    244       ! Find rn_zref in the array of model level depths and find the ref     
    245       ! density/temperature by linear interpolation.                                    
    246       DO jk = jpkm1, 2, -1  
    247          WHERE ( gdept_n(:,:,jk) > rn_zref )  
    248            ik_ref(:,:) = jk - 1  
    249            zT_ref(:,:) = zT(:,:,jk-1) + zdTdz(:,:,jk-1) * ( rn_zref - gdept_n(:,:,jk-1) )  
    250          END WHERE  
    251       END DO  
    252  
    253       ! If the first grid box centre is below the reference level then use the  
    254       ! top model level to get zT_ref  
    255       WHERE ( gdept_n(:,:,1) > rn_zref )   
    256          zT_ref = zT(:,:,1)  
    257          ik_ref = 1  
    258       END WHERE  
    259  
    260       ! The number of active tracer levels is 1 less than the number of active w levels  
    261       ikmt(:,:) = mbkt(:,:) - 1  
    262  
    263       ! Initialize / reset 
    264       ll_found(:,:) = .false. 
    265  
    266       IF ( rn_iso_frac - zepsilon > 0. ) THEN 
    267          ! Search for a uniform density/temperature region where adjacent levels           
    268          ! differ by less than rn_iso_frac * deltaT.                                       
    269          ! ik_iso is the index of the last level in the uniform layer   
    270          ! ll_found indicates whether the mixed layer depth can be found by interpolation  
    271          ik_iso(:,:)   = ik_ref(:,:)  
    272          DO jj = 1, nlcj  
    273             DO ji = 1, nlci  
    274 !CDIR NOVECTOR  
    275                DO jk = ik_ref(ji,jj), ikmt(ji,jj)-1  
    276                   IF ( zmoddT(ji,jj,jk) > ( rn_iso_frac * zdelta_T(ji,jj) ) ) THEN  
    277                      ik_iso(ji,jj)   = jk  
    278                      ll_found(ji,jj) = ( zmoddT(ji,jj,jk) > zdelta_T(ji,jj) )  
    279                      EXIT  
    280                   END IF  
    281                END DO  
    282             END DO  
    283          END DO  
    284  
    285          ! Use linear interpolation to find depth of mixed layer base where possible  
    286          hmld_zint(:,:) = rn_zref  
    287          DO jj = 1, jpj  
    288             DO ji = 1, jpi  
    289                IF (ll_found(ji,jj) .and. tmask(ji,jj,1) == 1.0) THEN  
    290                   zdz =  abs( zdelta_T(ji,jj) / zdTdz(ji,jj,ik_iso(ji,jj)) )  
    291                   hmld_zint(ji,jj) = gdept_n(ji,jj,ik_iso(ji,jj)) + zdz  
    292                END IF  
    293             END DO  
    294          END DO  
    295       END IF 
    296  
    297       ! If ll_found = .false. then calculate MLD using difference of zdelta_T     
    298       ! from the reference density/temperature  
    299   
    300 ! Prevent this section from working on land points  
    301       WHERE ( tmask(:,:,1) /= 1.0 )  
    302          ll_found = .true.  
    303       END WHERE  
    304   
    305       DO jk=1, jpk  
    306          ll_belowml(:,:,jk) = abs( zT(:,:,jk) - zT_ref(:,:) ) >= zdelta_T(:,:)   
    307       END DO  
    308   
    309 ! Set default value where interpolation cannot be used (ll_found=false)   
    310       DO jj = 1, jpj  
    311          DO ji = 1, jpi  
    312             IF ( .not. ll_found(ji,jj) )  hmld_zint(ji,jj) = gdept_n(ji,jj,ikmt(ji,jj))  
    313          END DO  
    314       END DO  
    315  
    316       DO jj = 1, jpj  
    317          DO ji = 1, jpi  
    318 !CDIR NOVECTOR  
    319             DO jk = ik_ref(ji,jj)+1, ikmt(ji,jj)  
    320                IF ( ll_found(ji,jj) ) EXIT  
    321                IF ( ll_belowml(ji,jj,jk) ) THEN                 
    322                   zT_b = zT_ref(ji,jj) + zdelta_T(ji,jj) * SIGN(1.0, zdTdz(ji,jj,jk-1) )  
    323                   zdT  = zT_b - zT(ji,jj,jk-1)                                       
    324                   zdz  = zdT / zdTdz(ji,jj,jk-1)                                        
    325                   hmld_zint(ji,jj) = gdept_n(ji,jj,jk-1) + zdz  
    326                   EXIT                                                    
    327                END IF  
    328             END DO  
    329          END DO  
    330       END DO  
    331  
    332       hmld_zint(:,:) = hmld_zint(:,:)*tmask(:,:,1)  
    333       !   
    334    END SUBROUTINE zdf_mxl_zint_mld 
    335  
    336    SUBROUTINE zdf_mxl_zint_htc( kt ) 
    337       !!---------------------------------------------------------------------- 
    338       !!                  ***  ROUTINE zdf_mxl_zint_htc  *** 
    339       !!  
    340       !! ** Purpose :    
    341       !! 
    342       !! ** Method  :    
    343       !!---------------------------------------------------------------------- 
    344  
    345       INTEGER, INTENT(in) ::   kt   ! ocean time-step index 
    346  
    347       INTEGER :: ji, jj, jk 
    348       INTEGER :: ikmax 
    349       REAL(wp) :: zc, zcoef 
    350       ! 
    351       INTEGER,  ALLOCATABLE, DIMENSION(:,:) ::   ilevel 
    352       REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::   zthick_0, zthick 
    353  
    354       !!---------------------------------------------------------------------- 
    355  
    356       IF( .NOT. ALLOCATED(ilevel) ) THEN 
    357          ALLOCATE( ilevel(jpi,jpj), zthick_0(jpi,jpj), & 
    358          &         zthick(jpi,jpj), STAT=ji ) 
    359          IF( lk_mpp  )   CALL mpp_sum(ji) 
    360          IF( ji /= 0 )   CALL ctl_stop( 'STOP', 'zdf_mxl_zint_htc : unable to allocate arrays' ) 
    361       ENDIF 
    362  
    363       ! Find last whole model T level above the MLD 
    364       ilevel(:,:)   = 0 
    365       zthick_0(:,:) = 0._wp 
    366  
    367       DO jk = 1, jpkm1   
    368          DO jj = 1, jpj 
    369             DO ji = 1, jpi                     
    370                zthick_0(ji,jj) = zthick_0(ji,jj) + e3t_n(ji,jj,jk) 
    371                IF( zthick_0(ji,jj) < hmld_zint(ji,jj) )   ilevel(ji,jj) = jk 
    372             END DO 
    373          END DO 
    374          WRITE(numout,*) 'zthick_0(jk =',jk,') =',zthick_0(2,2) 
    375          WRITE(numout,*) 'gdepw_n(jk+1 =',jk+1,') =',gdepw_n(2,2,jk+1) 
    376       END DO 
    377  
    378       ! Surface boundary condition 
    379       IF( ln_linssh ) THEN  ;   zthick(:,:) = sshn(:,:)   ;   htc_mld(:,:) = tsn(:,:,1,jp_tem) * sshn(:,:) * tmask(:,:,1)    
    380       ELSE                  ;   zthick(:,:) = 0._wp       ;   htc_mld(:,:) = 0._wp                                    
    381       ENDIF 
    382  
    383       ! Deepest whole T level above the MLD 
    384       ikmax = MIN( MAXVAL( ilevel(:,:) ), jpkm1 ) 
    385  
    386       ! Integration down to last whole model T level 
    387       DO jk = 1, ikmax 
    388          DO jj = 1, jpj 
    389             DO ji = 1, jpi 
    390                zc = e3t_n(ji,jj,jk) * REAL( MIN( MAX( 0, ilevel(ji,jj) - jk + 1 ) , 1  )  )    ! 0 below ilevel 
    391                zthick(ji,jj) = zthick(ji,jj) + zc 
    392                htc_mld(ji,jj) = htc_mld(ji,jj) + zc * tsn(ji,jj,jk,jp_tem) * tmask(ji,jj,jk) 
    393             END DO 
    394          END DO 
    395       END DO 
    396  
    397       ! Subsequent partial T level 
    398       zthick(:,:) = hmld_zint(:,:) - zthick(:,:)   !   remaining thickness to reach MLD 
    399  
    400       DO jj = 1, jpj 
    401          DO ji = 1, jpi 
    402             htc_mld(ji,jj) = htc_mld(ji,jj) + tsn(ji,jj,ilevel(ji,jj)+1,jp_tem)  &  
    403       &                      * MIN( e3t_n(ji,jj,ilevel(ji,jj)+1), zthick(ji,jj) ) * tmask(ji,jj,ilevel(ji,jj)+1) 
    404          END DO 
    405       END DO 
    406  
    407       WRITE(numout,*) 'htc_mld(after) =',htc_mld(2,2) 
    408  
    409       ! Convert to heat content 
    410       zcoef = rau0 * rcp 
    411       htc_mld(:,:) = zcoef * htc_mld(:,:) 
    412  
    413    END SUBROUTINE zdf_mxl_zint_htc 
    414  
    415    SUBROUTINE zdf_mxl_zint( kt ) 
    416       !!---------------------------------------------------------------------- 
    417       !!                  ***  ROUTINE zdf_mxl_zint  *** 
    418       !!  
    419       !! ** Purpose :    
    420       !! 
    421       !! ** Method  :    
    422       !!---------------------------------------------------------------------- 
    423  
    424       INTEGER, INTENT(in) ::   kt   ! ocean time-step index 
    425  
    426       INTEGER :: ios 
    427       INTEGER :: jn 
    428  
    429       INTEGER :: nn_mld_diag = 0    ! number of diagnostics 
    430  
    431       CHARACTER(len=1) :: cmld 
    432  
    433       TYPE(MXL_ZINT) :: sn_mld1, sn_mld2, sn_mld3, sn_mld4, sn_mld5 
    434       TYPE(MXL_ZINT), SAVE, DIMENSION(5) ::   mld_diags 
    435  
    436       NAMELIST/namzdf_mldzint/ nn_mld_diag, sn_mld1, sn_mld2, sn_mld3, sn_mld4, sn_mld5 
    437  
    438       !!---------------------------------------------------------------------- 
    439        
    440       IF( kt == nit000 ) THEN 
    441          REWIND( numnam_ref )              ! Namelist namzdf_mldzint in reference namelist  
    442          READ  ( numnam_ref, namzdf_mldzint, IOSTAT = ios, ERR = 901) 
    443 901      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namzdf_mldzint in reference namelist', lwp ) 
    444  
    445          REWIND( numnam_cfg )              ! Namelist namzdf_mldzint in configuration namelist  
    446          READ  ( numnam_cfg, namzdf_mldzint, IOSTAT = ios, ERR = 902 ) 
    447 902      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namzdf_mldzint in configuration namelist', lwp ) 
    448          IF(lwm) WRITE ( numond, namzdf_mldzint ) 
    449  
    450          IF( nn_mld_diag > 5 )   CALL ctl_stop( 'STOP', 'zdf_mxl_ini: Specify no more than 5 MLD definitions' ) 
    451  
    452          mld_diags(1) = sn_mld1 
    453          mld_diags(2) = sn_mld2 
    454          mld_diags(3) = sn_mld3 
    455          mld_diags(4) = sn_mld4 
    456          mld_diags(5) = sn_mld5 
    457  
    458          IF( lwp .AND. (nn_mld_diag > 0) ) THEN 
    459             WRITE(numout,*) '=============== Vertically-interpolated mixed layer ================' 
    460             WRITE(numout,*) '(Diagnostic number, nn_mld_type, rn_zref, rn_dT_crit, rn_iso_frac)' 
    461             DO jn = 1, nn_mld_diag 
    462                WRITE(numout,*) 'MLD criterion',jn,':' 
    463                WRITE(numout,*) '    nn_mld_type =', mld_diags(jn)%mld_type 
    464                WRITE(numout,*) '    rn_zref ='    , mld_diags(jn)%zref 
    465                WRITE(numout,*) '    rn_dT_crit =' , mld_diags(jn)%dT_crit 
    466                WRITE(numout,*) '    rn_iso_frac =', mld_diags(jn)%iso_frac 
    467             END DO 
    468             WRITE(numout,*) '====================================================================' 
    469          ENDIF 
    470       ENDIF 
    471  
    472       IF( nn_mld_diag > 0 ) THEN 
    473          DO jn = 1, nn_mld_diag 
    474             WRITE(cmld,'(I1)') jn 
    475             IF( iom_use( "mldzint_"//cmld ) .OR. iom_use( "mldhtc_"//cmld ) ) THEN 
    476                CALL zdf_mxl_zint_mld( mld_diags(jn) ) 
    477  
    478                IF( iom_use( "mldzint_"//cmld ) ) THEN 
    479                   CALL iom_put( "mldzint_"//cmld, hmld_zint(:,:) ) 
    480                ENDIF 
    481  
    482                IF( iom_use( "mldhtc_"//cmld ) )  THEN 
    483                   CALL zdf_mxl_zint_htc( kt ) 
    484                   CALL iom_put( "mldhtc_"//cmld , htc_mld(:,:)   ) 
    485                ENDIF 
    486             ENDIF 
    487          END DO 
    488       ENDIF 
    489  
    490    END SUBROUTINE zdf_mxl_zint 
    491  
    492143   !!====================================================================== 
    493144END MODULE zdfmxl 
  • branches/UKMO/dev_merge_2017_restart_datestamp_GO6_mixing/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftke.F90

    r9493 r9495  
    709709      IF( nn_mxl  < 0   .OR.  nn_mxl  > 3 )   CALL ctl_stop( 'bad flag: nn_mxl is  0, 1 or 2 ' ) 
    710710      IF( nn_pdl  < 0   .OR.  nn_pdl  > 1 )   CALL ctl_stop( 'bad flag: nn_pdl is  0 or 1    ' ) 
    711       IF( nn_htau < 0 .OR. nn_htau == 2 .OR. nn_htau > 5 )   CALL ctl_stop( 'bad flag: nn_htau is 0, 1, 3, 4 or 5 ' ) 
     711      IF( nn_htau < 0   .OR.  nn_htau > 1 )   CALL ctl_stop( 'bad flag: nn_htau is 0, 1 or 2 ' ) 
    712712      IF( nn_etau == 3 .AND. .NOT. ln_cpl )   CALL ctl_stop( 'nn_etau == 3 : HF taum only known in coupled mode' ) 
    713713      ! 
     
    722722      !                               !* depth of penetration of surface tke 
    723723      IF( nn_etau /= 0 ) THEN       
    724          htau(:,:) = 0._wp 
    725724         SELECT CASE( nn_htau )             ! Choice of the depth of penetration 
    726725         CASE( 0 )                                 ! constant depth penetration (here 10 meters) 
     
    728727         CASE( 1 )                                 ! F(latitude) : 0.5m to 30m poleward of 40 degrees 
    729728            htau(:,:) = MAX(  0.5_wp, MIN( 30._wp, 45._wp* ABS( SIN( rpi/180._wp * gphit(:,:) ) ) )   )             
    730          CASE( 3 )                                 ! F(latitude) : 0.5m to 15m poleward of 20 degrees 
    731             htau(:,:) = MAX(  0.5_wp, MIN( 15._wp, 45._wp* ABS( SIN( rpi/180._wp * gphit(:,:) ) ) )   ) 
    732          CASE( 4 )                                 ! F(latitude) : 0.5m to 10m/30m poleward of 13/40 degrees north/south 
    733             DO jj = 2, jpjm1 
    734                DO ji = fs_2, fs_jpim1   ! vector opt. 
    735                   IF( gphit(ji,jj) <= 0._wp ) THEN 
    736                      htau(ji,jj) = MAX(  0.5_wp, MIN( 30._wp, 45._wp* ABS( SIN( rpi/180._wp * gphit(ji,jj) ) ) )   ) 
    737                   ELSE 
    738                      htau(ji,jj) = MAX(  0.5_wp, MIN( 10._wp, 45._wp* ABS( SIN( rpi/180._wp * gphit(ji,jj) ) ) )   ) 
    739                   ENDIF 
    740                END DO 
    741             END DO 
    742          CASE ( 5 )                                ! F(latitude) : 0.5m to 10m poleward of 13 degrees north/south, 
    743             DO jj = 2, jpjm1                       !               10m to 30m between 30/45 degrees south 
    744                DO ji = fs_2, fs_jpim1   ! vector opt. 
    745                   IF( gphit(ji,jj) <= -30._wp ) THEN 
    746                      htau(ji,jj) = MAX(  10._wp, MIN( 30._wp, 55._wp* ABS( SIN( rpi/120._wp * ( gphit(ji,jj) + 23._wp ) ) ) )   ) 
    747                   ELSE 
    748                      htau(ji,jj) = MAX(  0.5_wp, MIN( 10._wp, 45._wp* ABS( SIN( rpi/180._wp * gphit(ji,jj) ) ) )   ) 
    749                   ENDIF 
    750                END DO 
    751             END DO 
    752729         END SELECT 
    753          ! 
    754730      ENDIF 
    755731      !                                !* read or initialize all required files 
Note: See TracChangeset for help on using the changeset viewer.