Changeset 7233


Ignore:
Timestamp:
2016-11-15T17:44:18+01:00 (4 years ago)
Author:
jpaul
Message:

see ticket #1781

Location:
branches/2016/dev_r6999_CONFIGMAN_1/NEMOGCM/TOOLS/SIREN/src
Files:
7 edited

Legend:

Unmodified
Added
Removed
  • branches/2016/dev_r6999_CONFIGMAN_1/NEMOGCM/TOOLS/SIREN/src/create_meshmask.f90

    r7153 r7233  
    6666!> 
    6767!>    * _input files namelist (namin)_:<br/> 
    68 !>       - cn_bathy  : Bathymetry file 
    69 !>       - cn_coord  : coordinate file (in_mshhgr=0) 
    70 !>       - cn_isfdep : Iceshelf draft  (ln_isfcav=true) 
    71 !>       - in_perio  : NEMO periodicity 
    72 !>       - ln_closea : 
     68!>       - cn_bathy     : Bathymetry file 
     69!>       - cn_varbathy  : Bathymetry variable name 
     70!>       - cn_coord     : coordinate file (in_mshhgr=0) 
     71!>       - cn_isfdep    : Iceshelf draft  (ln_isfcav=true) 
     72!>       - cn_varisfdep : Iceshelf draft variable name (ln_isfcav=true) 
     73!>       - in_perio     : NEMO periodicity 
     74!>       - ln_closea    : 
    7375!> 
    7476!>    * _horizontal grid namelist (namhgr)_:<br/> 
     
    164166!       - ln_zoom   : use zoom (namzoom) 
    165167!>       - ln_c1d    : use configuration 1D 
     168!>       - ln_e3_dep : vertical scale factors =T: e3.=dk[depth] =F: old definition 
    166169! 
    167170!    * _zoom namelist (namzoom)_:<br/> 
     
    175178! 
    176179!>    * _output namelist (namout)_:<br/> 
     180!>       - cn_domcfg    : output file name 
    177181!>       - in_msh       : number of output file and contain (0-9) 
    178182!>       - in_nproc     : number of processor to be used 
     
    204208!> - do not use anymore special case for ORCA grid 
    205209!> - allow to write domain_cfg file 
     210!> @date November, 2016 
     211!> - choose vertical scale factors (e3.=dk[depth] or old definition) 
    206212!> 
    207213!> @note Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
     
    298304   ! namin 
    299305   CHARACTER(LEN=lc) :: cn_bathy    = '' 
     306   CHARACTER(LEN=lc) :: cn_varbathy = '' 
    300307   CHARACTER(LEN=lc) :: cn_coord    = '' 
    301308   CHARACTER(LEN=lc) :: cn_isfdep   = '' 
     309   CHARACTER(LEN=lc) :: cn_varisfdep= 'isfdraft' 
    302310   INTEGER(i4)       :: in_perio    = -1 
    303311   LOGICAL           :: ln_closea   = .TRUE. 
     
    313321 
    314322   ! namlbc 
    315    REAL(sp)          :: rn_shlat    = 2. 
     323   REAL(sp)          :: rn_shlat    = -100. !2. 
    316324 
    317325   ! namout 
     326   CHARACTER(LEN=lc) :: cn_domcfg   = 'domain_cfg.nc' 
    318327   INTEGER(i4)       :: in_msh      = 0  
    319328   CHARACTER(LEN=lc) :: cn_type     = 'cdf' 
     
    334343   NAMELIST /namin/  &  !< input namelist 
    335344   &  cn_bathy,      &  !< Bathymetry file 
     345   &  cn_varbathy,   &  !< Bathymetry variable name 
    336346   &  cn_coord,      &  !< Coordinate file (in_mshhgr=0) 
    337347   &  cn_isfdep,     &  !< Iceshelf draft  (ln_isfcav=true) 
     348   &  cn_varisfdep,  &  !< Iceshelf draft variable name (ln_isfcav=true) 
    338349   &  in_perio,      &  !< NEMO periodicity 
    339350   &  ln_closea 
     
    352363 
    353364   NAMELIST /namout/ &  !< output namlist 
     365   &  cn_domcfg,     &  !< output file name 
    354366   &  in_msh,        &  !< number of output file (1,2,3) 
    355367   &  in_nproc,      &  !< number of processor to be used 
     
    449461   WRITE(*,*) 'DIMENSION TO BE USED :',jpi,jpj,jpk 
    450462 
    451    IF( ln_zco )   THEN 
    452       WRITE(*,*) 'VARIABLE READ : Bathy_level' 
    453       tl_bathy=iom_mpp_read_var(tl_mpp, 'Bathy_level') 
    454    ELSEIF( ln_zps .OR. ln_sco )   THEN 
    455       IF ( ln_isfcav ) THEN 
    456          WRITE(*,*) 'VARIABLE READ : Bathymetry_isf' 
    457          tl_bathy=iom_mpp_read_var(tl_mpp, 'Bathymetry_isf' ) 
    458       ELSE 
    459          WRITE(*,*) 'VARIABLE READ : Bathymetry' 
    460          tl_bathy=iom_mpp_read_var(tl_mpp, 'Bathymetry' ) 
    461       END IF 
     463   ! read variable 
     464   IF( TRIM(cn_varbathy) == '' )THEN 
     465      IF( ln_zco )THEN 
     466         cn_varbathy='Bathy_level' 
     467      ELSEIF( ln_zps .OR. ln_sco )THEN 
     468         IF( ln_isfcav )THEN 
     469            cn_varbathy='Bathymetry_isf' 
     470         ELSE 
     471            cn_varbathy='Bathymetry' 
     472         ENDIF 
     473      ENDIF 
    462474   ENDIF 
     475   WRITE(*,*) 'VARIABLE READ : '//TRIM(cn_varbathy) 
     476   tl_bathy=iom_mpp_read_var(tl_mpp, TRIM(cn_varbathy)) 
    463477   CALL iom_mpp_close(tl_mpp) 
    464478 
     
    469483  
    470484   IF ( ln_isfcav ) THEN 
     485      WRITE(*,*) 'ICESHELF DRAFT FILE TO BE USED:',TRIM(cn_isfdep) 
     486      WRITE(*,*) 'ICESHELF VARIABLE READ : '//TRIM(cn_varisfdep) 
    471487      ! open Iceshelf draft 
    472488      IF( cn_isfdep /= '' )THEN 
     
    481497      CALL iom_mpp_open(tl_mpp) 
    482498      IF( ln_zps .OR. ln_sco )   THEN 
    483          tl_risfdep=iom_mpp_read_var(tl_mpp, 'isf_draft') 
     499         tl_risfdep=iom_mpp_read_var(tl_mpp, cn_varisfdep) 
    484500      ENDIF 
    485501      CALL iom_mpp_close(tl_mpp) 
     
    488504      dl_tmp2D(:,:)=0._dp 
    489505 
    490       tl_risfdep=var_init('isf_draft',dl_tmp2D(:,:), id_type=NF90_FLOAT) 
     506      tl_risfdep=var_init(cn_varisfdep, dl_tmp2D(:,:), id_type=NF90_DOUBLE) 
    491507 
    492508      DEALLOCATE(dl_tmp2D) 
     
    561577         !                                  ! create 'domain_cfg.nc' file 
    562578         !                                  ! ============================ 
    563          tl_mppout0=mpp_init( 'domain_cfg', tg_tmask, & 
     579         tl_mppout0=mpp_init( cn_domcfg, tg_tmask, & 
    564580         &                    in_niproc, in_njproc, in_nproc, & 
    565581         &                    cd_type=cn_type ) 
     
    789805   &                            tg_mbkt%d_value(:,:,:,:) 
    790806   ! 
    791    IF( ll_domcfg ) tg_mbathy%c_name='bottom_level' 
     807   IF( ll_domcfg )THEN 
     808      tg_mbathy%c_name='bottom_level' 
     809   ENDIF 
    792810   CALL mpp_add_var(tl_mppzgr, tg_mbathy) 
    793811   CALL var_clean(tg_mbathy) 
     
    862880      CALL mpp_add_var(tl_mppzgr, tg_e3v_0) 
    863881      CALL var_clean(tg_e3v_0) 
     882      ! e3f_0 
     883      CALL mpp_add_var(tl_mppzgr, tg_e3f_0) 
     884      CALL var_clean(tg_e3f_0) 
    864885      ! e3w_0 
    865886      CALL mpp_add_var(tl_mppzgr, tg_e3w_0) 
    866887      CALL var_clean(tg_e3w_0) 
     888      ! e3uw_0 
     889      CALL mpp_add_var(tl_mppzgr, tg_e3uw_0) 
     890      CALL var_clean(tg_e3uw_0) 
     891      ! e3vw_0 
     892      CALL mpp_add_var(tl_mppzgr, tg_e3vw_0) 
     893      CALL var_clean(tg_e3vw_0) 
    867894 
    868895      ! Max. grid stiffness ratio 
     
    873900 
    874901      ! stretched system 
    875       ! gdept_1d 
    876       CALL mpp_add_var(tl_mppzgr, tg_gdept_1d) 
    877       CALL var_clean(tg_gdept_1d) 
    878       ! gdepw_1d 
    879       CALL mpp_add_var(tl_mppzgr, tg_gdepw_1d) 
    880       CALL var_clean(tg_gdepw_1d) 
    881        
    882       ! gdept_0 
    883       CALL mpp_add_var(tl_mppzgr, tg_gdept_0)       
    884       CALL var_clean(tg_gdept_0) 
    885       ! gdepw_0 
    886       CALL mpp_add_var(tl_mppzgr, tg_gdepw_0)       
    887       CALL var_clean(tg_gdepw_0) 
     902      IF( .NOT. tl_namz%l_e3_dep )THEN 
     903         ! gdept_1d 
     904         CALL mpp_add_var(tl_mppzgr, tg_gdept_1d) 
     905         CALL var_clean(tg_gdept_1d) 
     906         ! gdepw_1d 
     907         CALL mpp_add_var(tl_mppzgr, tg_gdepw_1d) 
     908         CALL var_clean(tg_gdepw_1d) 
     909          
     910         ! gdept_0 
     911         CALL mpp_add_var(tl_mppzgr, tg_gdept_0)       
     912         CALL var_clean(tg_gdept_0) 
     913         ! gdepw_0 
     914         CALL mpp_add_var(tl_mppzgr, tg_gdepw_0)       
     915         CALL var_clean(tg_gdepw_0) 
     916      ENDIF 
    888917 
    889918   ENDIF 
     
    928957      IF( ll_domcfg .OR. in_msh <= 3 ) THEN ! 3D depth 
    929958          
    930          ! gdept_0 
    931          CALL mpp_add_var(tl_mppzgr, tg_gdept_0) 
     959         IF( .NOT. tl_namz%l_e3_dep )THEN 
    932960       
    933          ! gdepu, gdepv 
    934          IF( .NOT. ll_domcfg )THEN 
    935             ALLOCATE(dl_tmp3D(jpi,jpj,jpk)) 
    936             dl_tmp3D(:,:,:)=dp_fill 
    937  
    938             tl_gdepu=var_init('gdepu',dl_tmp3D(:,:,:), id_type=NF90_FLOAT) 
    939             tl_gdepv=var_init('gdepv',dl_tmp3D(:,:,:), id_type=NF90_FLOAT) 
    940  
    941             DEALLOCATE(dl_tmp3D) 
    942             DO jk = 1,jpk    
    943                DO jj = 1, jpj-1    
    944                   DO ji = 1, jpi-1   ! vector opt. 
    945                      tl_gdepu%d_value(ji,jj,jk,1) = MIN( tg_gdept_0%d_value(ji  ,jj  ,jk,1) , & 
    946                         &                                tg_gdept_0%d_value(ji+1,jj  ,jk,1) ) 
    947  
    948                      tl_gdepv%d_value(ji,jj,jk,1) = MIN( tg_gdept_0%d_value(ji  ,jj  ,jk,1) , & 
    949                         &                                tg_gdept_0%d_value(ji  ,jj+1,jk,1) ) 
     961            ! gdepu, gdepv 
     962            IF( .NOT. ll_domcfg )THEN 
     963               ALLOCATE(dl_tmp3D(jpi,jpj,jpk)) 
     964               dl_tmp3D(:,:,:)=dp_fill 
     965 
     966               tl_gdepu=var_init('gdepu',dl_tmp3D(:,:,:), id_type=NF90_FLOAT) 
     967               tl_gdepv=var_init('gdepv',dl_tmp3D(:,:,:), id_type=NF90_FLOAT) 
     968 
     969               DEALLOCATE(dl_tmp3D) 
     970               DO jk = 1,jpk    
     971                  DO jj = 1, jpj-1    
     972                     DO ji = 1, jpi-1   ! vector opt. 
     973                        tl_gdepu%d_value(ji,jj,jk,1) = MIN( tg_gdept_0%d_value(ji  ,jj  ,jk,1) , & 
     974                           &                                tg_gdept_0%d_value(ji+1,jj  ,jk,1) ) 
     975 
     976                        tl_gdepv%d_value(ji,jj,jk,1) = MIN( tg_gdept_0%d_value(ji  ,jj  ,jk,1) , & 
     977                           &                                tg_gdept_0%d_value(ji  ,jj+1,jk,1) ) 
     978                     END DO    
    950979                  END DO    
    951                END DO    
    952             END DO          
    953             CALL lbc_lnk( tl_gdepu%d_value(:,:,:,1), 'U', in_perio, 1._dp ) 
    954             CALL lbc_lnk( tl_gdepv%d_value(:,:,:,1), 'V', in_perio, 1._dp ) 
    955  
    956             ! gdepu 
    957             CALL mpp_add_var(tl_mppzgr, tl_gdepu) 
    958             CALL var_clean(tl_gdepu) 
    959             ! gdepv 
    960             CALL mpp_add_var(tl_mppzgr, tl_gdepv) 
    961             CALL var_clean(tl_gdepv) 
     980               END DO          
     981               CALL lbc_lnk( tl_gdepu%d_value(:,:,:,1), 'U', in_perio, 1._dp ) 
     982               CALL lbc_lnk( tl_gdepv%d_value(:,:,:,1), 'V', in_perio, 1._dp ) 
     983 
     984               ! gdepu 
     985               CALL mpp_add_var(tl_mppzgr, tl_gdepu) 
     986               CALL var_clean(tl_gdepu) 
     987               ! gdepv 
     988               CALL mpp_add_var(tl_mppzgr, tl_gdepv) 
     989               CALL var_clean(tl_gdepv) 
     990            ENDIF 
     991 
     992            ! gdept_0 
     993            CALL mpp_add_var(tl_mppzgr, tg_gdept_0) 
     994            CALL var_clean(tg_gdept_0) 
     995 
     996            ! gdepw_0 
     997            CALL mpp_add_var(tl_mppzgr, tg_gdepw_0) 
     998            CALL var_clean(tg_gdepw_0) 
    962999         ENDIF 
    963  
    964          ! clean 
    965          CALL var_clean(tg_gdept_0) 
    966  
    967          ! gdepw_0 
    968          CALL mpp_add_var(tl_mppzgr, tg_gdepw_0) 
    969          CALL var_clean(tg_gdepw_0) 
    9701000 
    9711001      ELSE ! 2D bottom depth 
     
    9981028   ENDIF 
    9991029 
    1000    IF( ln_zps .OR. ln_zco )THEN ! z-coordinate  
    1001       ! depth 
    1002       ! gdept_1d 
    1003       CALL mpp_add_var(tl_mppzgr, tg_gdept_1d) 
    1004       CALL var_clean(tg_gdept_1d) 
    1005       ! gdepw_1d 
    1006       CALL mpp_add_var(tl_mppzgr, tg_gdepw_1d) 
    1007       CALL var_clean(tg_gdepw_1d) 
    1008       ! scale factors 
     1030   ! scale factors 
     1031   IF( ll_domcfg )THEN 
    10091032      ! e3t_1d 
    10101033      CALL mpp_add_var(tl_mppzgr, tg_e3t_1d) 
     
    10151038   ENDIF 
    10161039 
     1040   IF( ln_zps .OR. ln_zco )THEN ! z-coordinate  
     1041      IF( .NOT. tl_namz%l_e3_dep )THEN 
     1042         ! depth 
     1043         ! gdept_1d 
     1044         CALL mpp_add_var(tl_mppzgr, tg_gdept_1d) 
     1045         CALL var_clean(tg_gdept_1d) 
     1046         ! gdepw_1d 
     1047         CALL mpp_add_var(tl_mppzgr, tg_gdepw_1d) 
     1048         CALL var_clean(tg_gdepw_1d) 
     1049      ENDIF 
     1050   ENDIF 
     1051 
    10171052   ! define global attributes 
    10181053   ALLOCATE(tl_gatt(ip_maxatt)) 
     
    10211056 
    10221057 
     1058   IF( in_msh == 0 ) in_msh=1 
    10231059   SELECT CASE ( MOD(in_msh, 3) ) 
    10241060      CASE ( 1 ) 
     
    14041440!      ENDIF 
    14051441 
    1406       IF( .NOT. ld_domcfg )THEN 
    14071442      ! 2. Ocean/land mask at u-,  v-, and z-points (computed from tmask) 
    14081443      ! ------------------------------------------- 
     
    14171452            !END DO 
    14181453            !DO ji = 1, jpi-1  ! NO vector opt. 
    1419                tg_fmask%d_value(ji,jj,jk,1) = tg_tmask%d_value(ji  ,jj  ,jk,1) * & 
    1420                                             & tg_tmask%d_value(ji+1,jj  ,jk,1) * & 
    1421                                             & tg_tmask%d_value(ji  ,jj+1,jk,1) * & 
    1422                                             & tg_tmask%d_value(ji+1,jj+1,jk,1) 
     1454               IF( .NOT. ld_domcfg )THEN 
     1455                  tg_fmask%d_value(ji,jj,jk,1) = tg_tmask%d_value(ji  ,jj  ,jk,1) * & 
     1456                                               & tg_tmask%d_value(ji+1,jj  ,jk,1) * & 
     1457                                               & tg_tmask%d_value(ji  ,jj+1,jk,1) * & 
     1458                                               & tg_tmask%d_value(ji+1,jj+1,jk,1) 
     1459               ENDIF 
    14231460            ENDDO 
    14241461         ENDDO 
     
    14451482      CALL lbc_lnk( tg_umask%d_value  (:,:,:,1), 'U', td_nam%i_perio, 1._dp )      ! Lateral boundary conditions 
    14461483      CALL lbc_lnk( tg_vmask%d_value  (:,:,:,1), 'V', td_nam%i_perio, 1._dp ) 
    1447       CALL lbc_lnk( tg_fmask%d_value  (:,:,:,1), 'F', td_nam%i_perio, 1._dp ) 
     1484      IF( .NOT. ld_domcfg )THEN 
     1485         CALL lbc_lnk( tg_fmask%d_value  (:,:,:,1), 'F', td_nam%i_perio, 1._dp ) 
     1486      ENDIF 
    14481487!      CALL lbc_lnk( tg_ssumask%d_value(:,:,:,1), 'U', td_nam%i_perio, 1._dp )      ! Lateral boundary conditions 
    14491488!      CALL lbc_lnk( tg_ssvmask%d_value(:,:,:,1), 'V', td_nam%i_perio, 1._dp ) 
     
    14671506      ! Lateral boundary conditions on velocity (modify fmask) 
    14681507      ! ---------------------------------------      
    1469       ALLOCATE( zwf(jpi,jpj) ) 
    1470       DO jk = 1, jpk 
    1471          zwf(:,:) = tg_fmask%d_value(:,:,jk,1)          
    1472          DO jj = 2, jpj-1 
    1473             DO ji = 2, jpi-1   ! vector opt. 
    1474                IF( tg_fmask%d_value(ji,jj,jk,1) == 0._dp )THEN 
    1475                   tg_fmask%d_value(ji,jj,jk,1) = rn_shlat * & 
    1476                                                & MIN(1._dp , MAX(zwf(ji+1,jj), zwf(ji,jj+1), & 
    1477                                                &                 zwf(ji-1,jj), zwf(ji,jj-1)) ) 
     1508      IF( .NOT. ld_domcfg )THEN 
     1509         ALLOCATE( zwf(jpi,jpj) ) 
     1510         DO jk = 1, jpk 
     1511            zwf(:,:) = tg_fmask%d_value(:,:,jk,1)          
     1512            DO jj = 2, jpj-1 
     1513               DO ji = 2, jpi-1   ! vector opt. 
     1514                  IF( tg_fmask%d_value(ji,jj,jk,1) == 0._dp )THEN 
     1515                     tg_fmask%d_value(ji,jj,jk,1) = rn_shlat * & 
     1516                                                  & MIN(1._dp , MAX(zwf(ji+1,jj), zwf(ji,jj+1), & 
     1517                                                  &                 zwf(ji-1,jj), zwf(ji,jj-1)) ) 
     1518                  ENDIF 
     1519               END DO 
     1520            END DO 
     1521            DO jj = 2, jpj-1 
     1522               IF( tg_fmask%d_value(1,jj,jk,1) == 0._dp )THEN 
     1523                  tg_fmask%d_value(1  ,jj,jk,1) = rn_shlat * & 
     1524                                                & MIN(1._dp, MAX(zwf(2,jj), zwf(1,jj+1), zwf(1,jj-1))) 
     1525               ENDIF 
     1526               IF( tg_fmask%d_value(jpi,jj,jk,1) == 0._dp )THEN 
     1527                  tg_fmask%d_value(jpi,jj,jk,1) = rn_shlat * & 
     1528                                                & MIN(1._wp, MAX(zwf(jpi,jj+1), zwf(jpi-1,jj), zwf(jpi,jj-1))) 
     1529               ENDIF 
     1530            END DO          
     1531            DO ji = 2, jpi-1 
     1532               IF( tg_fmask%d_value(ji,1,jk,1) == 0._dp )THEN 
     1533                  tg_fmask%d_value(ji, 1 ,jk,1) = rn_shlat * & 
     1534                                                & MIN(1._dp, MAX(zwf(ji+1,1), zwf(ji,2), zwf(ji-1,1))) 
     1535               ENDIF 
     1536               IF( tg_fmask%d_value(ji,jpj,jk,1) == 0._dp )THEN 
     1537                  tg_fmask%d_value(ji,jpj,jk,1) = rn_shlat * & 
     1538                                                & MIN(1._dp, MAX(zwf(ji+1,jpj), zwf(ji-1,jpj), zwf(ji,jpj-1))) 
    14781539               ENDIF 
    14791540            END DO 
    14801541         END DO 
    1481          DO jj = 2, jpj-1 
    1482             IF( tg_fmask%d_value(1,jj,jk,1) == 0._dp )THEN 
    1483                tg_fmask%d_value(1  ,jj,jk,1) = rn_shlat * & 
    1484                                              & MIN(1._dp, MAX(zwf(2,jj), zwf(1,jj+1), zwf(1,jj-1))) 
    1485             ENDIF 
    1486             IF( tg_fmask%d_value(jpi,jj,jk,1) == 0._dp )THEN 
    1487                tg_fmask%d_value(jpi,jj,jk,1) = rn_shlat * & 
    1488                                              & MIN(1._wp, MAX(zwf(jpi,jj+1), zwf(jpi-1,jj), zwf(jpi,jj-1))) 
    1489             ENDIF 
    1490          END DO          
    1491          DO ji = 2, jpi-1 
    1492             IF( tg_fmask%d_value(ji,1,jk,1) == 0._dp )THEN 
    1493                tg_fmask%d_value(ji, 1 ,jk,1) = rn_shlat * & 
    1494                                              & MIN(1._dp, MAX(zwf(ji+1,1), zwf(ji,2), zwf(ji-1,1))) 
    1495             ENDIF 
    1496             IF( tg_fmask%d_value(ji,jpj,jk,1) == 0._dp )THEN 
    1497                tg_fmask%d_value(ji,jpj,jk,1) = rn_shlat * & 
    1498                                              & MIN(1._dp, MAX(zwf(ji+1,jpj), zwf(ji-1,jpj), zwf(ji,jpj-1))) 
    1499             ENDIF 
    1500          END DO 
    1501       END DO 
    1502       DEALLOCATE( zwf ) 
     1542         DEALLOCATE( zwf ) 
    15031543 
    15041544!      IF( td_nam%c_cfg == "orca" .AND. td_nam%i_cfg == 2 )THEN   ! ORCA_R2 configuration 
     
    15761616!      ENDIF 
    15771617      ! 
    1578       CALL lbc_lnk( tg_fmask%d_value(:,:,:,1), 'F', td_nam%i_perio, 1._dp )      ! Lateral boundary conditions on fmask 
     1618         CALL lbc_lnk( tg_fmask%d_value(:,:,:,1), 'F', td_nam%i_perio, 1._dp )      ! Lateral boundary conditions on fmask 
     1619      ENDIF 
    15791620 
    15801621!      DEALLOCATE( imsk ) 
    1581       ENDIF ! ld_domcfg 
    15821622 
    15831623!      DEALLOCATE( dl_tpol ) 
     
    17001740            ji=ji+1 ; td_att(ji)=att_init("bb",td_namz%d_bb) 
    17011741         ELSEIF( td_namz%l_s_sf12 )THEN 
    1702             IF( td_namz%l_sigcrit ) ji=ji+1 ; td_att(ji)=att_init("sigma below critical depth","activated") 
     1742            IF( td_namz%l_sigcrit ) ji=ji+1 ; td_att(ji)=att_init("sigma_below_critical_depth","activated") 
    17031743            ji=ji+1 ; td_att(ji)=att_init("alpha",td_namz%d_alpha) 
    17041744            ji=ji+1 ; td_att(ji)=att_init("efold",td_namz%d_efold) 
  • branches/2016/dev_r6999_CONFIGMAN_1/NEMOGCM/TOOLS/SIREN/src/grid_hgr.f90

    r7153 r7233  
    273273!      tg_ssfmask  = var_init('ssfmask',dl_tmp2D(:,:)  , dd_fill=dp_fill_i1, id_type=NF90_BYTE) 
    274274 
    275       dl_tmp2D(:,:)=dp_fill_sp 
    276  
    277       tg_glamt    = var_init('glamt',dl_tmp2D(:,:)   , dd_fill=dp_fill_sp, id_type=NF90_FLOAT) 
    278       tg_glamu    = var_init('glamu',dl_tmp2D(:,:)   , dd_fill=dp_fill_sp, id_type=NF90_FLOAT) 
    279       tg_glamv    = var_init('glamv',dl_tmp2D(:,:)   , dd_fill=dp_fill_sp, id_type=NF90_FLOAT) 
    280       tg_glamf    = var_init('glamf',dl_tmp2D(:,:)   , dd_fill=dp_fill_sp, id_type=NF90_FLOAT) 
    281  
    282       tg_gphit    = var_init('gphit',dl_tmp2D(:,:)   , dd_fill=dp_fill_sp, id_type=NF90_FLOAT) 
    283       tg_gphiu    = var_init('gphiu',dl_tmp2D(:,:)   , dd_fill=dp_fill_sp, id_type=NF90_FLOAT) 
    284       tg_gphiv    = var_init('gphiv',dl_tmp2D(:,:)   , dd_fill=dp_fill_sp, id_type=NF90_FLOAT) 
    285       tg_gphif    = var_init('gphif',dl_tmp2D(:,:)   , dd_fill=dp_fill_sp, id_type=NF90_FLOAT) 
    286  
    287275      dl_tmp2D(:,:)=dp_fill 
     276 
     277      tg_glamt    = var_init('glamt',dl_tmp2D(:,:)   , dd_fill=dp_fill, id_type=NF90_DOUBLE) 
     278      tg_glamu    = var_init('glamu',dl_tmp2D(:,:)   , dd_fill=dp_fill, id_type=NF90_DOUBLE) 
     279      tg_glamv    = var_init('glamv',dl_tmp2D(:,:)   , dd_fill=dp_fill, id_type=NF90_DOUBLE) 
     280      tg_glamf    = var_init('glamf',dl_tmp2D(:,:)   , dd_fill=dp_fill, id_type=NF90_DOUBLE) 
     281 
     282      tg_gphit    = var_init('gphit',dl_tmp2D(:,:)   , dd_fill=dp_fill, id_type=NF90_DOUBLE) 
     283      tg_gphiu    = var_init('gphiu',dl_tmp2D(:,:)   , dd_fill=dp_fill, id_type=NF90_DOUBLE) 
     284      tg_gphiv    = var_init('gphiv',dl_tmp2D(:,:)   , dd_fill=dp_fill, id_type=NF90_DOUBLE) 
     285      tg_gphif    = var_init('gphif',dl_tmp2D(:,:)   , dd_fill=dp_fill, id_type=NF90_DOUBLE) 
    288286 
    289287      tg_e1t      = var_init('e1t'  ,dl_tmp2D(:,:)   , dd_fill=dp_fill, id_type=NF90_DOUBLE) 
     
    316314 
    317315      tg_tmask   = var_init('tmask'  ,dl_tmp3D(:,:,:), dd_fill=dp_fill_i1, id_type=NF90_BYTE) 
     316      tg_umask   = var_init('umask'  ,dl_tmp3D(:,:,:), dd_fill=dp_fill_i1, id_type=NF90_BYTE) 
     317      tg_vmask   = var_init('vmask'  ,dl_tmp3D(:,:,:), dd_fill=dp_fill_i1, id_type=NF90_BYTE) 
    318318      IF( .NOT. ld_domcfg )THEN 
    319          tg_umask   = var_init('umask'  ,dl_tmp3D(:,:,:), dd_fill=dp_fill_i1, id_type=NF90_BYTE) 
    320          tg_vmask   = var_init('vmask'  ,dl_tmp3D(:,:,:), dd_fill=dp_fill_i1, id_type=NF90_BYTE) 
    321319         tg_fmask   = var_init('fmask'  ,dl_tmp3D(:,:,:), dd_fill=dp_fill_i1, id_type=NF90_BYTE) 
    322320      ENDIF 
     
    380378 
    381379      CALL var_clean(tg_tmask   ) 
     380      CALL var_clean(tg_umask   ) 
     381      CALL var_clean(tg_vmask   ) 
    382382      IF( .NOT. ld_domcfg )THEN 
    383          CALL var_clean(tg_umask   ) 
    384          CALL var_clean(tg_vmask   ) 
    385383         CALL var_clean(tg_fmask   ) 
    386384      ENDIF 
     
    533531      ! loop indices 
    534532      !---------------------------------------------------------------- 
    535       CALL logger_info('GRID HGR FILL : define the horizontal mesh from ithe'//& 
     533      CALL logger_info('GRID HGR FILL : define the horizontal mesh from the'//& 
    536534      &  ' type of horizontal mesh mshhgr = '//TRIM(fct_str(td_nam%i_mshhgr))) 
    537535      IF( td_nam%i_mshhgr == 1 )THEN 
     
    657655 
    658656      ! force output type 
    659       tg_glamt%i_type=NF90_FLOAT 
    660       tg_glamu%i_type=NF90_FLOAT 
    661       tg_glamv%i_type=NF90_FLOAT 
    662       tg_glamf%i_type=NF90_FLOAT 
    663  
    664       tg_gphit%i_type=NF90_FLOAT 
    665       tg_gphiu%i_type=NF90_FLOAT 
    666       tg_gphiv%i_type=NF90_FLOAT 
    667       tg_gphif%i_type=NF90_FLOAT 
     657      tg_glamt%i_type=NF90_DOUBLE 
     658      tg_glamu%i_type=NF90_DOUBLE 
     659      tg_glamv%i_type=NF90_DOUBLE 
     660      tg_glamf%i_type=NF90_DOUBLE 
     661 
     662      tg_gphit%i_type=NF90_DOUBLE 
     663      tg_gphiu%i_type=NF90_DOUBLE 
     664      tg_gphiv%i_type=NF90_DOUBLE 
     665      tg_gphif%i_type=NF90_DOUBLE 
    668666 
    669667      tg_e1t  =iom_mpp_read_var(tl_coord, 'e1t') 
  • branches/2016/dev_r6999_CONFIGMAN_1/NEMOGCM/TOOLS/SIREN/src/grid_zgr.f90

    r7153 r7233  
    5555!> - J, Paul : update from trunk (revision 6961): add wetting and drying, ice sheet coupling.. 
    5656!> - J, Paul : do not use anymore special case for ORCA grid. 
     57!> @date November, 2016 
     58!> - J, Paul : vertical scale factors e3. = dk[gdep] or old definition 
    5759!> 
    5860!> @note Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
     
    99101   PUBLIC :: tg_e3v_0 
    100102   PUBLIC :: tg_e3w_0 
    101 !   PUBLIC :: tg_e3f_0     !useless to create meshmask 
    102 !   PUBLIC :: tg_e3uw_0    !useless to create meshmask 
    103 !   PUBLIC :: tg_e3vw_0    !useless to create meshmask 
     103   PUBLIC :: tg_e3f_0     !useless to create meshmask 
     104   PUBLIC :: tg_e3uw_0    !useless to create meshmask 
     105   PUBLIC :: tg_e3vw_0    !useless to create meshmask 
    104106    
    105107   PUBLIC :: tg_mbkt 
     
    144146   PRIVATE :: grid_zgr__isf_fill 
    145147!   PRIVATE :: grid_zgr__isf_fill_e3x 
    146 !   PRIVATE :: grid_zgr__isf_fill_e3uw 
     148   PRIVATE :: grid_zgr__isf_fill_e3uw 
    147149!   PRIVATE :: grid_zgr__isf_fill_gdep3w_0 
    148150   PRIVATE :: grid_zgr__sco_fill 
     
    155157 
    156158   TYPE TNAMZ 
     159 
    157160      CHARACTER(LEN=lc) :: c_coord    
    158161      INTEGER(i4)       :: i_perio    
     
    214217!      LOGICAL           :: l_zoom 
    215218      LOGICAL           :: l_c1d 
     219      LOGICAL           :: l_e3_dep 
    216220                 
    217221!      CHARACTER(LEN=lc) :: c_cfz       
     
    244248   TYPE(TVAR), SAVE :: tg_e3v_0     !      zps & sco 
    245249   TYPE(TVAR), SAVE :: tg_e3w_0     !      zps & sco 
    246    !TYPE(TVAR), SAVE :: tg_e3f_0 
    247    !TYPE(TVAR), SAVE :: tg_e3uw_0 
    248    !TYPE(TVAR), SAVE :: tg_e3vw_0 
     250   TYPE(TVAR), SAVE :: tg_e3f_0 
     251   TYPE(TVAR), SAVE :: tg_e3uw_0 
     252   TYPE(TVAR), SAVE :: tg_e3vw_0 
    249253    
    250254   TYPE(TVAR), SAVE :: tg_mbkt      !zco & zps & sco 
     
    315319      dl_tmp2D(:,:)  =dp_fill_i2 
    316320 
    317       tg_mbathy  =var_init('mbathy  ',dl_tmp2D(:,:)  , dd_fill=dp_fill_i2, id_type=NF90_SHORT) 
    318       tg_misfdep =var_init('misfdep ',dl_tmp2D(:,:)  , dd_fill=dp_fill_i2, id_type=NF90_SHORT) 
    319  
    320321      tg_mbkt    =var_init('mbkt    ',dl_tmp2D(:,:)  , dd_fill=dp_fill_i2, id_type=NF90_SHORT) 
    321322      !tg_mbku    =var_init('mbku    ',dl_tmp2D(:,:)  , dd_fill=dp_fill_i2, id_type=NF90_SHORT) 
     
    326327      !tg_mikf    =var_init('mikf    ',dl_tmp2D(:,:)  , dd_fill=dp_fill_i2, id_type=NF90_SHORT) 
    327328 
     329      dl_tmp2D(:,:)  =dp_fill_i4 
     330 
     331      tg_mbathy  =var_init('mbathy  ',dl_tmp2D(:,:)  , dd_fill=dp_fill_i4, id_type=NF90_INT) 
     332      tg_misfdep =var_init('misfdep ',dl_tmp2D(:,:)  , dd_fill=dp_fill_i4, id_type=NF90_INT) 
     333 
    328334      dl_tmp2D(:,:)  =dp_fill 
    329335 
     
    337343 
    338344      ! variable 3D 
    339       dl_tmp3D(:,:,:)=dp_fill_sp 
    340  
    341       tg_gdept_0 =var_init('gdept_0 ',dl_tmp3D(:,:,:), dd_fill=dp_fill_sp, id_type=NF90_FLOAT) 
    342       tg_gdepw_0 =var_init('gdepw_0 ',dl_tmp3D(:,:,:), dd_fill=dp_fill_sp, id_type=NF90_FLOAT) 
     345      dl_tmp3D(:,:,:)=dp_fill 
     346 
     347      tg_gdept_0 =var_init('gdept_0 ',dl_tmp3D(:,:,:), dd_fill=dp_fill, id_type=NF90_DOUBLE) 
     348      tg_gdepw_0 =var_init('gdepw_0 ',dl_tmp3D(:,:,:), dd_fill=dp_fill, id_type=NF90_DOUBLE) 
    343349      !tg_gdep3w_0=var_init('gdep3w_0',dl_tmp3D(:,:,:), dd_fill=dp_fill_sp, id_type=NF90_FLOAT) 
    344350 
    345       dl_tmp3D(:,:,:)=dp_fill 
    346        
    347351      tg_e3t_0   =var_init('e3t_0   ',dl_tmp3D(:,:,:), dd_fill=dp_fill, id_type=NF90_DOUBLE) 
    348352      tg_e3u_0   =var_init('e3u_0   ',dl_tmp3D(:,:,:), dd_fill=dp_fill, id_type=NF90_DOUBLE) 
    349353      tg_e3v_0   =var_init('e3v_0   ',dl_tmp3D(:,:,:), dd_fill=dp_fill, id_type=NF90_DOUBLE) 
    350354      tg_e3w_0   =var_init('e3w_0   ',dl_tmp3D(:,:,:), dd_fill=dp_fill, id_type=NF90_DOUBLE) 
    351       !tg_e3f_0   =var_init('e3f_0   ',dl_tmp3D(:,:,:), dd_fill=dp_fill, id_type=NF90_DOUBLE) 
    352       !tg_e3uw_0  =var_init('e3uw_0  ',dl_tmp3D(:,:,:), dd_fill=dp_fill, id_type=NF90_DOUBLE) 
    353       !tg_e3vw_0  =var_init('e3vw_0  ',dl_tmp3D(:,:,:), dd_fill=dp_fill, id_type=NF90_DOUBLE) 
     355      tg_e3f_0   =var_init('e3f_0   ',dl_tmp3D(:,:,:), dd_fill=dp_fill, id_type=NF90_DOUBLE) 
     356      tg_e3uw_0  =var_init('e3uw_0  ',dl_tmp3D(:,:,:), dd_fill=dp_fill, id_type=NF90_DOUBLE) 
     357      tg_e3vw_0  =var_init('e3vw_0  ',dl_tmp3D(:,:,:), dd_fill=dp_fill, id_type=NF90_DOUBLE) 
    354358 
    355359   END SUBROUTINE grid_zgr_init 
     
    409413      CALL var_clean(tg_e3u_0   ) 
    410414      CALL var_clean(tg_e3v_0   ) 
    411       !CALL var_clean(tg_e3f_0   ) 
    412       !CALL var_clean(tg_e3uw_0  ) 
    413       !CALL var_clean(tg_e3vw_0  ) 
     415      CALL var_clean(tg_e3f_0   ) 
     416      CALL var_clean(tg_e3uw_0  ) 
     417      CALL var_clean(tg_e3vw_0  ) 
    414418       
    415419   END SUBROUTINE grid_zgr_clean 
     
    493497      REAL(dp)          :: dn_zb_b     = NF90_FILL_DOUBLE 
    494498 
    495       ! namcla 
    496       INTEGER(i4)       :: in_cla      = 0 
     499!      ! namcla 
     500!      INTEGER(i4)       :: in_cla      = 0 
    497501 
    498502      ! namwd 
     
    507511!      LOGICAL           :: ln_zoom     = .FALSE. 
    508512      LOGICAL           :: ln_c1d      = .FALSE. 
     513      LOGICAL           :: ln_e3_dep   = .FALSE. 
    509514 
    510515!      ! namzoom 
     
    568573      &   dn_zb_b         !<  offset for calculating Zb 
    569574 
    570       NAMELIST /namcla/ & 
    571       &  in_cla            !< =1 cross land advection for exchanges through some straits (ORCA2) 
     575!      NAMELIST /namcla/ & 
     576!      &  in_cla            !< =1 cross land advection for exchanges through some straits (ORCA2) 
    572577 
    573578      NAMELIST /namwd/ &  !< wetting and drying 
     
    581586!      &  in_bench,      &  !< benchmark parameter (in_mshhgr = 5 ) 
    582587!      &  ln_zoom,       &  !< use zoom 
    583       &  ln_c1d            !< use configuration 1D 
     588      &  ln_c1d,         &  !< use configuration 1D 
     589      &  ln_e3_dep          !< new vertical scale factors [T, F:old definition] 
    584590 
    585591!      NAMELIST /namzoom/& 
     
    616622      IF( ln_zps    ) READ( il_fileid, NML = namzps  ) 
    617623      IF( ln_sco    ) READ( il_fileid, NML = namsco  ) 
    618       READ( il_fileid, NML = namcla  ) 
     624!      READ( il_fileid, NML = namcla  ) 
    619625      READ( il_fileid, NML = namwd   ) 
    620626      READ( il_fileid, NML = namgrd  ) 
     
    676682      grid_zgr_nam%d_zb_b     = dn_zb_b 
    677683 
    678       grid_zgr_nam%i_cla      = in_cla 
     684!      grid_zgr_nam%i_cla      = in_cla 
    679685 
    680686      grid_zgr_nam%d_wdmin1   = dn_wdmin1   
     
    687693!      grid_zgr_nam%l_zoom     = ln_zoom 
    688694      grid_zgr_nam%l_c1d      = ln_c1d 
     695      grid_zgr_nam%l_e3_dep   = ln_e3_dep 
    689696 
    690697!      grid_zgr_nam%c_cfz      = cn_cfz    
     
    743750      CALL logger_info(' s- or hybrid z-s-coordinate    ln_sco    = '//TRIM(fct_str(td_nam%l_sco))) 
    744751      CALL logger_info(' ice shelf cavities             ln_isfcav = '//TRIM(fct_str(td_nam%l_isfcav))) 
     752      CALL logger_info(' vertical scale factors         ln_e3_dep = '//TRIM(fct_str(td_nam%l_e3_dep))) 
    745753 
    746754      il_count=0 
     
    759767      ENDIF 
    760768 
     769      IF(.NOT. td_nam%l_e3_dep )THEN 
     770         CALL logger_info("Obsolescent definition of e3 scale factors is used") 
     771      ENDIF 
    761772      ! Build the vertical coordinate system 
    762773      ! ------------------------------------ 
     
    773784      ! z-coordinate 
    774785      IF( td_nam%l_zco ) CALL grid_zgr__zco(jpk) 
    775        
     786 
    776787      ! Partial step z-coordinate 
    777788      IF( td_nam%l_zps ) CALL grid_zgr__zps_fill( td_nam,jpi,jpj,jpk,td_bathy,td_risfdep ) 
    778        
     789 
    779790      ! s-coordinate or hybrid z-s coordinate 
    780791      IF( td_nam%l_sco ) CALL grid_zgr__sco_fill( td_nam,jpi,jpj,jpk,td_bathy ) 
     
    807818         !&   ' 3w '//TRIM(fct_str(MINVAL( tg_gdep3w_0%d_value(:,:,:,1) )))//& 
    808819         &   '  t '//TRIM(fct_str(MINVAL( tg_e3t_0%d_value(:,:,:,1) )))//& 
    809          !&   '  f '//TRIM(fct_str(MINVAL( tg_e3f_0%d_value(:,:,:,1) )))//& 
     820         &   '  f '//TRIM(fct_str(MINVAL( tg_e3f_0%d_value(:,:,:,1) )))//& 
    810821         &   '  u '//TRIM(fct_str(MINVAL( tg_e3u_0%d_value(:,:,:,1) )))//& 
    811822         &   '  v '//TRIM(fct_str(MINVAL( tg_e3v_0%d_value(:,:,:,1) )))//& 
    812          !&   ' uw '//TRIM(fct_str(MINVAL( tg_e3uw_0%d_value(:,:,:,1) )))//& 
    813          !&   ' vw '//TRIM(fct_str(MINVAL( tg_e3vw_0%d_value(:,:,:,1) )))//& 
     823         &   ' uw '//TRIM(fct_str(MINVAL( tg_e3uw_0%d_value(:,:,:,1) )))//& 
     824         &   ' vw '//TRIM(fct_str(MINVAL( tg_e3vw_0%d_value(:,:,:,1) )))//& 
    814825         &   '  w '//TRIM(fct_str(MINVAL( tg_e3w_0%d_value(:,:,:,1) ))) ) 
    815826      CALL logger_info(' MAX val depth t '//TRIM(fct_str(MAXVAL( tg_gdept_0%d_value(:,:,:,1) )))//& 
     
    817828         !&   ' 3w '//TRIM(fct_str(MAXVAL( tg_gdep3w_0%d_value(:,:,:,1) ))) ) 
    818829      CALL logger_info(' MAX val e3    t '//TRIM(fct_str(MAXVAL( tg_e3t_0%d_value(:,:,:,1) )))//& 
    819          !&   ' f '//TRIM(fct_str(MAXVAL( tg_e3f_0%d_value(:,:,:,1) )))//& 
     830         &   ' f '//TRIM(fct_str(MAXVAL( tg_e3f_0%d_value(:,:,:,1) )))//& 
    820831         &   ' u '//TRIM(fct_str(MAXVAL( tg_e3u_0%d_value(:,:,:,1) )))//& 
    821832         &   ' v '//TRIM(fct_str(MAXVAL( tg_e3v_0%d_value(:,:,:,1) )))//& 
    822          !&   ' uw '//TRIM(fct_str(MAXVAL( tg_e3uw_0%d_value(:,:,:,1) )))//& 
    823          !&   ' vw '//TRIM(fct_str(MAXVAL( tg_e3vw_0%d_value(:,:,:,1) )))//& 
     833         &   ' uw '//TRIM(fct_str(MAXVAL( tg_e3uw_0%d_value(:,:,:,1) )))//& 
     834         &   ' vw '//TRIM(fct_str(MAXVAL( tg_e3vw_0%d_value(:,:,:,1) )))//& 
    824835         &   ' w '//TRIM(fct_str(MAXVAL( tg_e3w_0%d_value(:,:,:,1) ))) ) 
    825836 
     
    889900         &  td_nam%d_ppsur == dp_pp_to_be_computed           ) THEN 
    890901         ! 
    891          za1  = (  zdzmin - zhmax / FLOAT(jpk-1)  )                                                 & 
    892             & / ( TANH((1-zkth)/zacr) - zacr/FLOAT(jpk-1) * (  LOG( COSH( (jpk - zkth) / zacr) )    & 
    893             &                                                - LOG( COSH( ( 1  - zkth) / zacr) )  ) ) 
     902         za1  = (  zdzmin - zhmax / REAL(jpk-1,dp)  )                                               & 
     903            & / ( TANH((1-zkth)/zacr) - zacr/REAL(jpk-1,dp) * ( LOG( COSH( (jpk - zkth) / zacr) )   & 
     904            &                                               - LOG( COSH( ( 1  - zkth) / zacr) ) ) ) 
    894905         za0  = zdzmin  - za1 *             TANH( (1-zkth) / zacr ) 
    895906         zsur =   - za0 - za1 * zacr * LOG( COSH( (1-zkth) / zacr )  ) 
     
    904915      CALL logger_info(' GRID ZGR Z : Reference vertical z-coordinates') 
    905916      CALL logger_info('~~~~~~~~~~~') 
    906       IF( zkth == 0._wp ) THEN  
     917      IF( zkth == 0._dp ) THEN  
    907918         CALL logger_info('Uniform grid with '//TRIM(fct_str(jpk-1))//' layers') 
    908919         CALL logger_info('Total depth    :'//TRIM(fct_str(zhmax))) 
    909920         CALL logger_info('Layer thickness:'//TRIM(fct_str(zhmax/(jpk-1))))          
    910921      ELSE 
    911          IF( za1 == 0._wp .AND. za0 == 0._wp .AND. zsur == 0._wp ) THEN 
     922         IF( za1 == 0._dp .AND. za0 == 0._dp .AND. zsur == 0._dp ) THEN 
    912923            CALL logger_info('zsur, za0, za1 computed from ') 
    913924            CALL logger_info('        zdzmin = '//TRIM(fct_str(zdzmin))) 
     
    931942      ! init 
    932943      IF( zkth == 0._dp ) THEN            !  uniform vertical grid        
    933          za1 = zhmax / FLOAT(jpk-1)  
     944         za1 = zhmax / REAL(jpk-1,dp)  
    934945         DO jk = 1, jpk 
    935             zw = FLOAT( jk ) 
    936             zt = FLOAT( jk ) + 0.5_dp 
     946            zw = REAL( jk, dp ) 
     947            zt = REAL( jk, dp ) + 0.5_dp 
    937948            tg_gdepw_1d%d_value(1,1,jk,1) = ( zw - 1 ) * za1 
    938949            tg_gdept_1d%d_value(1,1,jk,1) = ( zt - 1 ) * za1 
     
    943954         IF( .NOT. td_nam%l_dbletanh ) THEN 
    944955            DO jk = 1, jpk 
    945                zw = REAL( jk , wp ) 
    946                zt = REAL( jk , wp ) + 0.5_dp 
     956               zw = REAL( jk , dp ) 
     957               zt = REAL( jk , dp ) + 0.5_dp 
    947958               tg_gdepw_1d%d_value(1,1,jk,1) = ( zsur + za0 * zw + za1 * zacr * LOG ( COSH( (zw-zkth) / zacr ) )  ) 
    948959               tg_gdept_1d%d_value(1,1,jk,1) = ( zsur + za0 * zt + za1 * zacr * LOG ( COSH( (zt-zkth) / zacr ) )  ) 
     
    952963         ELSE 
    953964            DO jk = 1, jpk 
    954                zw = FLOAT( jk ) 
    955                zt = FLOAT( jk ) + 0.5_dp 
     965               zw = REAL( jk, dp ) 
     966               zt = REAL( jk, dp ) + 0.5_dp 
    956967               ! Double tanh function 
    957968               tg_gdepw_1d%d_value(1,1,jk,1) = ( zsur + za0 * zw + za1 * zacr * LOG ( COSH( (zw-zkth ) / zacr  ) )    & 
     
    968979      ENDIF       
    969980 
    970       IF ( td_nam%l_isfcav ) THEN 
     981      IF ( td_nam%l_isfcav .OR. td_nam%l_e3_dep ) THEN 
    971982         ! need to be like this to compute the pressure gradient with ISF.  
    972983         ! If not, level beneath the ISF are not aligned (sum(e3t) /= depth) 
     
    9891000!       
    9901001!      ! ref. depth with tolerance (10% of minimum layer thickness) 
    991 !      zrefdep = 10._wp - 0.1_wp * MINVAL( tg_e3w_1d%d_value(1,1,:,1) ) 
     1002!      zrefdep = 10._dp - 0.1_dp * MINVAL( tg_e3w_1d%d_value(1,1,:,1) ) 
    9921003!       
    9931004!      ! shallowest W level Below ~10m 
     
    11261137            WHERE ( td_bathy%d_value(:,:,1,1) <= td_risfdep%d_value(:,:,1,1) + td_nam%d_isfhmin ) 
    11271138               tg_misfdep%d_value(:,:,1,1) = 0 
    1128                td_risfdep%d_value(:,:,1,1) = 0._wp 
     1139               td_risfdep%d_value(:,:,1,1) = 0._dp 
    11291140               tg_mbathy%d_value (:,:,1,1) = 0 
    1130                td_bathy%d_value  (:,:,1,1) = 0._wp 
     1141               td_bathy%d_value  (:,:,1,1) = 0._dp 
    11311142            END WHERE 
    11321143         END IF 
     
    12111222         tg_e3u_0%d_value   (:,:,jk,1) = tg_e3t_1d%d_value  (1,1,jk,1) 
    12121223         tg_e3v_0%d_value   (:,:,jk,1) = tg_e3t_1d%d_value  (1,1,jk,1) 
    1213          !tg_e3f_0%d_value   (:,:,jk,1) = tg_e3t_1d%d_value  (1,1,jk,1) 
     1224         tg_e3f_0%d_value   (:,:,jk,1) = tg_e3t_1d%d_value  (1,1,jk,1) 
    12141225         tg_e3w_0%d_value   (:,:,jk,1) = tg_e3w_1d%d_value  (1,1,jk,1) 
    1215          !tg_e3uw_0%d_value  (:,:,jk,1) = tg_e3w_1d%d_value  (1,1,jk,1) 
    1216          !tg_e3vw_0%d_value  (:,:,jk,1) = tg_e3w_1d%d_value  (1,1,jk,1) 
     1226         tg_e3uw_0%d_value  (:,:,jk,1) = tg_e3w_1d%d_value  (1,1,jk,1) 
     1227         tg_e3vw_0%d_value  (:,:,jk,1) = tg_e3w_1d%d_value  (1,1,jk,1) 
    12171228      END DO 
    12181229 
     
    17601771                     tg_gdept_0%d_value(ji,jj,ik+1,1) = tg_gdept_0%d_value(ji,jj,ik,1) + tg_e3t_0%d_value(ji,jj,ik,1) 
    17611772                  ENDIF 
    1762                   !gm Bug?  check the gdepw_1d 
    1763                   !       ... on ik 
    1764                   tg_gdept_0%d_value(ji,jj,ik,1) = tg_gdepw_1d%d_value ( 1, 1,ik  ,1)   & 
    1765                                                 &    + ( tg_gdepw_0%d_value  (ji,jj,ik+1,1) - tg_gdepw_1d%d_value(1,1,ik,1) )& 
    1766                                                 &    * ( (tg_gdept_1d%d_value( 1, 1,ik  ,1) - tg_gdepw_1d%d_value(1,1,ik,1)) & 
    1767                                                 &      / (tg_gdepw_1d%d_value( 1, 1,ik+1,1) - tg_gdepw_1d%d_value(1,1,ik,1)) ) 
    1768  
    1769                   tg_e3t_0%d_value  (ji,jj,ik,1) = tg_e3t_1d%d_value  ( 1, 1,ik  ,1)   & 
    1770                                                 &    * ( tg_gdepw_0%d_value (ji,jj,ik+1,1) - tg_gdepw_1d%d_value(1,1,ik,1) ) & 
    1771                                                 &    / ( tg_gdepw_1d%d_value( 1, 1,ik+1,1) - tg_gdepw_1d%d_value(1,1,ik,1) )  
    1772  
    1773                   tg_e3w_0%d_value  (ji,jj,ik,1) = 0.5_dp   & 
    1774                                                 &     *  ( tg_gdepw_0%d_value (ji,jj,ik+1,1)   & 
    1775                                                 &        + tg_gdepw_1d%d_value( 1, 1,ik+1,1)   & 
    1776                                                 &        - tg_gdepw_1d%d_value( 1, 1,ik  ,1)*2._dp) & 
    1777                                                 &     *  ( tg_e3w_1d%d_value  ( 1, 1,ik  ,1)   & 
    1778                                                 &        / ( tg_gdepw_1d%d_value( 1, 1,ik+1,1)   & 
    1779                                                 &          - tg_gdepw_1d%d_value( 1, 1,ik  ,1) ) ) 
    1780  
    1781                   !       ... on ik+1 
    1782                   tg_e3w_0%d_value  (ji,jj,ik+1,1) = tg_e3t_0%d_value  (ji,jj,ik,1) 
    1783                   tg_e3t_0%d_value  (ji,jj,ik+1,1) = tg_e3t_0%d_value  (ji,jj,ik,1) 
    1784                   tg_gdept_0%d_value(ji,jj,ik+1,1) = tg_gdept_0%d_value(ji,jj,ik,1) + tg_e3t_0%d_value(ji,jj,ik,1) 
    17851773               ENDIF 
    17861774            END DO 
     
    18241812         tg_e3u_0%d_value (:,:,jk,1) = tg_e3t_1d%d_value(1,1,jk,1) 
    18251813         tg_e3v_0%d_value (:,:,jk,1) = tg_e3t_1d%d_value(1,1,jk,1) 
    1826          !tg_e3uw_0%d_value(:,:,jk,1) = tg_e3w_1d%d_value(1,1,jk,1) 
    1827          !tg_e3vw_0%d_value(:,:,jk,1) = tg_e3w_1d%d_value(1,1,jk,1) 
     1814         tg_e3uw_0%d_value(:,:,jk,1) = tg_e3w_1d%d_value(1,1,jk,1) 
     1815         tg_e3vw_0%d_value(:,:,jk,1) = tg_e3w_1d%d_value(1,1,jk,1) 
    18281816      END DO 
    18291817 
     
    18341822               tg_e3u_0%d_value (ji,jj,jk,1) = MIN( tg_e3t_0%d_value(ji,jj,jk,1), tg_e3t_0%d_value(ji+1,jj  ,jk,1) ) 
    18351823               tg_e3v_0%d_value (ji,jj,jk,1) = MIN( tg_e3t_0%d_value(ji,jj,jk,1), tg_e3t_0%d_value(ji  ,jj+1,jk,1) ) 
    1836                !tg_e3uw_0%d_value(ji,jj,jk,1) = MIN( tg_e3w_0%d_value(ji,jj,jk,1), tg_e3w_0%d_value(ji+1,jj  ,jk,1) ) 
    1837                !tg_e3vw_0%d_value(ji,jj,jk,1) = MIN( tg_e3w_0%d_value(ji,jj,jk,1), tg_e3w_0%d_value(ji  ,jj+1,jk,1) ) 
     1824               tg_e3uw_0%d_value(ji,jj,jk,1) = MIN( tg_e3w_0%d_value(ji,jj,jk,1), tg_e3w_0%d_value(ji+1,jj  ,jk,1) ) 
     1825               tg_e3vw_0%d_value(ji,jj,jk,1) = MIN( tg_e3w_0%d_value(ji,jj,jk,1), tg_e3w_0%d_value(ji  ,jj+1,jk,1) ) 
    18381826            END DO 
    18391827         END DO 
    18401828      END DO 
    18411829 
    1842       !IF ( td_nam%l_isfcav ) THEN 
    1843       !   ! (ISF) define e3uw (adapted for 2 cells in the water column) 
    1844       !   CALL grid_zgr__isf_fill_e3uw(jpi,jpj) 
    1845       !END IF 
     1830      IF ( td_nam%l_isfcav ) THEN 
     1831         ! (ISF) define e3uw (adapted for 2 cells in the water column) 
     1832         CALL grid_zgr__isf_fill_e3uw(jpi,jpj) 
     1833      END IF 
    18461834 
    18471835      ! lateral boundary conditions 
    18481836      CALL lbc_lnk( tg_e3u_0%d_value (:,:,:,1), 'U', td_nam%i_perio, 1._dp ) 
    18491837      CALL lbc_lnk( tg_e3v_0%d_value (:,:,:,1), 'V', td_nam%i_perio, 1._dp ) 
    1850       !CALL lbc_lnk( tg_e3uw_0%d_value(:,:,:,1), 'U', td_nam%i_perio, 1._dp ) 
    1851       !CALL lbc_lnk( tg_e3vw_0%d_value(:,:,:,1), 'V', td_nam%i_perio, 1._dp ) 
     1838      CALL lbc_lnk( tg_e3uw_0%d_value(:,:,:,1), 'U', td_nam%i_perio, 1._dp ) 
     1839      CALL lbc_lnk( tg_e3vw_0%d_value(:,:,:,1), 'V', td_nam%i_perio, 1._dp ) 
    18521840 
    18531841      ! set to z-scale factor if zero (i.e. along closed boundaries) 
     
    18551843         WHERE( tg_e3u_0%d_value (:,:,jk,1) == 0._dp )   tg_e3u_0%d_value (:,:,jk,1) = tg_e3t_1d%d_value(1,1,jk,1) 
    18561844         WHERE( tg_e3v_0%d_value (:,:,jk,1) == 0._dp )   tg_e3v_0%d_value (:,:,jk,1) = tg_e3t_1d%d_value(1,1,jk,1) 
    1857          !WHERE( tg_e3uw_0%d_value(:,:,jk,1) == 0._dp )   tg_e3uw_0%d_value(:,:,jk,1) = tg_e3w_1d%d_value(1,1,jk,1) 
    1858          !WHERE( tg_e3vw_0%d_value(:,:,jk,1) == 0._dp )   tg_e3vw_0%d_value(:,:,jk,1) = tg_e3w_1d%d_value(1,1,jk,1) 
     1845         WHERE( tg_e3uw_0%d_value(:,:,jk,1) == 0._dp )   tg_e3uw_0%d_value(:,:,jk,1) = tg_e3w_1d%d_value(1,1,jk,1) 
     1846         WHERE( tg_e3vw_0%d_value(:,:,jk,1) == 0._dp )   tg_e3vw_0%d_value(:,:,jk,1) = tg_e3w_1d%d_value(1,1,jk,1) 
    18591847      END DO 
    18601848 
     
    20001988 
    20011989      ! remove very shallow ice shelf (less than ~ 10m if 75L) 
    2002       WHERE( td_risfdep%d_value(:,:,1,1) <= 10._wp .AND. & 
     1990      WHERE( td_risfdep%d_value(:,:,1,1) <= 10._dp .AND. & 
    20031991           & tg_misfdep%d_value(:,:,1,1) > 1 ) 
    20041992         tg_misfdep%d_value(:,:,1,1) = 0 
    2005          td_risfdep%d_value(:,:,1,1) = 0.0_wp 
     1993         td_risfdep%d_value(:,:,1,1) = 0.0_dp 
    20061994         tg_mbathy%d_value (:,:,1,1) = 0 
    2007          td_bathy%d_value  (:,:,1,1) = 0.0_wp 
     1995         td_bathy%d_value  (:,:,1,1) = 0.0_dp 
    20081996      END WHERE 
    2009       WHERE( td_bathy%d_value(:,:,1,1) <= 30.0_wp .AND. & 
    2010            & tg_gphit%d_value(:,:,1,1) < -60._wp ) 
     1997      WHERE( td_bathy%d_value(:,:,1,1) <= 30.0_dp .AND. & 
     1998           & tg_gphit%d_value(:,:,1,1) < -60._dp ) 
    20111999         tg_misfdep%d_value(:,:,1,1) = 0 
    2012          td_risfdep%d_value(:,:,1,1) = 0.0_wp 
     2000         td_risfdep%d_value(:,:,1,1) = 0.0_dp 
    20132001         tg_mbathy%d_value (:,:,1,1) = 0 
    2014          td_bathy%d_value  (:,:,1,1) = 0.0_wp 
     2002         td_bathy%d_value  (:,:,1,1) = 0.0_dp 
    20152003      END WHERE 
    20162004 
     
    26692657                  zdepwp = td_bathy%d_value(ji,jj,1,1) 
    26702658                  ze3tp  = td_bathy%d_value(ji,jj,1,1) - tg_gdepw_1d%d_value(1,1,ik,1) 
    2671                   ze3wp  = 0.5_wp * tg_e3w_1d%d_value(1,1,ik,1) & 
    2672                      &            * ( 1._wp + ( ze3tp/tg_e3t_1d%d_value(1,1,ik,1) ) ) 
     2659                  ze3wp  = 0.5_dp * tg_e3w_1d%d_value(1,1,ik,1) & 
     2660                     &            * ( 1._dp + ( ze3tp/tg_e3t_1d%d_value(1,1,ik,1) ) ) 
    26732661                  tg_e3t_0%d_value  (ji,jj,ik  ,1) = ze3tp 
    26742662                  tg_e3t_0%d_value  (ji,jj,ik+1,1) = ze3tp 
     
    28952883! 
    28962884!   END SUBROUTINE grid_zgr__isf_fill_e3x 
    2897 !   !------------------------------------------------------------------- 
    2898 !   !> @brief This subroutine define e3uw  
    2899 !   !>    (adapted for 2 cells in the water column) for ISF case  
    2900 !   !> 
    2901 !   !> @details 
    2902 !   !> 
    2903 !   !> @author J.Paul 
    2904 !   !> @date September, 2015 - rewrite from zgr_zps 
    2905 !   !> 
    2906 !   !> @param[in] jpi 
    2907 !   !> @param[in] jpj 
    2908 !   !------------------------------------------------------------------- 
    2909 !   SUBROUTINE grid_zgr__isf_fill_e3uw(jpi,jpj)  
    2910 !      IMPLICIT NONE 
    2911 !      ! Argument       
    2912 !      INTEGER(i4)   , INTENT(IN   ) :: jpi 
    2913 !      INTEGER(i4)   , INTENT(IN   ) :: jpj 
    2914 ! 
    2915 !      ! local variable 
    2916 !      INTEGER(i4) :: ikb, ikt 
    2917 ! 
    2918 !      ! loop indices 
    2919 !      INTEGER(i4) :: ji 
    2920 !      INTEGER(i4) :: jj 
    2921 !      !---------------------------------------------------------------- 
    2922 ! 
    2923 !      DO jj = 2, jpj - 1  
    2924 !         DO ji = 2, jpi - 1  
    2925 ! 
    2926 !            ikb = MAX(tg_mbathy%d_value (ji,jj,1,1), tg_mbathy%d_value (ji+1,jj,1,1)) 
    2927 !            ikt = MAX(tg_misfdep%d_value(ji,jj,1,1), tg_misfdep%d_value(ji+1,jj,1,1)) 
    2928 !            IF( ikb == ikt+1 )THEN 
    2929 !               tg_e3uw_0%d_value(ji,jj,ikb,1) = MIN( tg_gdept_0%d_value(ji,jj,ikb  ,1), tg_gdept_0%d_value(ji+1,jj  ,ikb  ,1) ) - & 
    2930 !               &                                MAX( tg_gdept_0%d_value(ji,jj,ikb-1,1), tg_gdept_0%d_value(ji+1,jj  ,ikb-1,1) ) 
    2931 !            ENDIF 
    2932 !             
    2933 !            ikb = MAX( tg_mbathy%d_value (ji,jj,1,1), tg_mbathy%d_value (ji,jj+1,1,1)) 
    2934 !            ikt = MAX( tg_misfdep%d_value(ji,jj,1,1), tg_misfdep%d_value(ji,jj+1,1,1)) 
    2935 !            IF( ikb == ikt+1 )THEN 
    2936 !               tg_e3vw_0%d_value(ji,jj,ikb,1) = MIN( tg_gdept_0%d_value(ji,jj,ikb  ,1), tg_gdept_0%d_value(ji  ,jj+1,ikb  ,1) ) - & 
    2937 !               &                                MAX( tg_gdept_0%d_value(ji,jj,ikb-1,1), tg_gdept_0%d_value(ji  ,jj+1,ikb-1,1) ) 
    2938 !            ENDIF 
    2939 !         END DO 
    2940 !      END DO 
    2941 ! 
    2942 !   END SUBROUTINE grid_zgr__isf_fill_e3uw 
     2885   !------------------------------------------------------------------- 
     2886   !> @brief This subroutine define e3uw  
     2887   !>    (adapted for 2 cells in the water column) for ISF case  
     2888   !> 
     2889   !> @details 
     2890   !> 
     2891   !> @author J.Paul 
     2892   !> @date September, 2015 - rewrite from zgr_zps 
     2893   !> 
     2894   !> @param[in] jpi 
     2895   !> @param[in] jpj 
     2896   !------------------------------------------------------------------- 
     2897   SUBROUTINE grid_zgr__isf_fill_e3uw(jpi,jpj)  
     2898      IMPLICIT NONE 
     2899      ! Argument       
     2900      INTEGER(i4)   , INTENT(IN   ) :: jpi 
     2901      INTEGER(i4)   , INTENT(IN   ) :: jpj 
     2902 
     2903      ! local variable 
     2904      INTEGER(i4) :: ikb, ikt 
     2905 
     2906      ! loop indices 
     2907      INTEGER(i4) :: ji 
     2908      INTEGER(i4) :: jj 
     2909      !---------------------------------------------------------------- 
     2910 
     2911      DO jj = 2, jpj - 1  
     2912         DO ji = 2, jpi - 1  
     2913 
     2914            ikb = MAX(tg_mbathy%d_value (ji,jj,1,1), tg_mbathy%d_value (ji+1,jj,1,1)) 
     2915            ikt = MAX(tg_misfdep%d_value(ji,jj,1,1), tg_misfdep%d_value(ji+1,jj,1,1)) 
     2916            IF( ikb == ikt+1 )THEN 
     2917               tg_e3uw_0%d_value(ji,jj,ikb,1) = MIN( tg_gdept_0%d_value(ji,jj,ikb  ,1), tg_gdept_0%d_value(ji+1,jj  ,ikb  ,1) ) - & 
     2918               &                                MAX( tg_gdept_0%d_value(ji,jj,ikb-1,1), tg_gdept_0%d_value(ji+1,jj  ,ikb-1,1) ) 
     2919            ENDIF 
     2920             
     2921            ikb = MAX( tg_mbathy%d_value (ji,jj,1,1), tg_mbathy%d_value (ji,jj+1,1,1)) 
     2922            ikt = MAX( tg_misfdep%d_value(ji,jj,1,1), tg_misfdep%d_value(ji,jj+1,1,1)) 
     2923            IF( ikb == ikt+1 )THEN 
     2924               tg_e3vw_0%d_value(ji,jj,ikb,1) = MIN( tg_gdept_0%d_value(ji,jj,ikb  ,1), tg_gdept_0%d_value(ji  ,jj+1,ikb  ,1) ) - & 
     2925               &                                MAX( tg_gdept_0%d_value(ji,jj,ikb-1,1), tg_gdept_0%d_value(ji  ,jj+1,ikb-1,1) ) 
     2926            ENDIF 
     2927         END DO 
     2928      END DO 
     2929 
     2930   END SUBROUTINE grid_zgr__isf_fill_e3uw 
    29432931!   !------------------------------------------------------------------- 
    29442932!   !> @brief This subroutine compute gdep3w_0 (vertical sum of e3w)  
     
    30123000      INTEGER(i4), INTENT(IN) :: jpi 
    30133001      INTEGER(i4), INTENT(IN) :: jpj 
    3014 !      LOGICAL    , INTENT(IN) :: ld_domcfg 
    30153002 
    30163003      ! local variable 
     
    31223109      REAL(dp), DIMENSION(:,:), ALLOCATABLE :: zrj 
    31233110      REAL(dp), DIMENSION(:,:), ALLOCATABLE :: zhbat 
    3124       REAL(wp), DIMENSION(:,:), ALLOCATABLE :: ztmpi1 
    3125       REAL(wp), DIMENSION(:,:), ALLOCATABLE :: ztmpi2 
    3126       REAL(wp), DIMENSION(:,:), ALLOCATABLE :: ztmpj1 
    3127       REAL(wp), DIMENSION(:,:), ALLOCATABLE :: ztmpj2 
     3111      REAL(dp), DIMENSION(:,:), ALLOCATABLE :: ztmpi1 
     3112      REAL(dp), DIMENSION(:,:), ALLOCATABLE :: ztmpi2 
     3113      REAL(dp), DIMENSION(:,:), ALLOCATABLE :: ztmpj1 
     3114      REAL(dp), DIMENSION(:,:), ALLOCATABLE :: ztmpj2 
    31283115 
    31293116      ! loop indices 
     
    31703157      ! set maximum ocean depth 
    31713158      td_bathy%d_value(:,:,1,1) = MIN( td_nam%d_sbot_max, td_bathy%d_value(:,:,1,1) ) 
    3172       IF( .NOT.td_nam%l_wd ) THEN 
     3159      IF( .NOT. td_nam%l_wd )THEN 
    31733160         DO jj = 1, jpj 
    31743161            DO ji = 1, jpi 
     
    31873174      zenv(:,:) = td_bathy%d_value(:,:,1,1) 
    31883175 
    3189       IF( .NOT. td_nam%l_wd ) THEN 
     3176      IF( .NOT. td_nam%l_wd )THEN 
    31903177      ! set first land point adjacent to a wet cell to sbot_min as this needs to be included in smoothing 
    31913178         DO jj = 1, jpj 
    31923179            DO ji = 1, jpi 
    3193               IF( td_bathy%d_value(ji,jj,1,1) == 0._dp ) THEN 
     3180              IF( td_bathy%d_value(ji,jj,1,1) == 0._dp )THEN 
    31943181                iip1 = MIN( ji+1, jpi ) 
    31953182                ijp1 = MIN( jj+1, jpj ) 
     
    32033190                   &  td_bathy%d_value(iim1,ijm1,1,1) + & 
    32043191                   &  td_bathy%d_value(ji  ,ijm1,1,1) + & 
    3205                    &  td_bathy%d_value(iip1,ijp1,1,1) ) > 0._dp ) THEN 
     3192                   &  td_bathy%d_value(iip1,ijp1,1,1) ) > 0._dp )THEN 
    32063193                  zenv(ji,jj) = td_nam%d_sbot_min 
    32073194                ENDIF 
     
    32123199 
    32133200      ! apply lateral boundary condition   CAUTION: keep the value when the lbc field is zero 
    3214       CALL lbc_lnk( zenv(:,:), 'T', td_nam%i_perio, 1._dp) ! ?? , 'no0' ) 
     3201      ! this is only in mpp case, so here just do nothing 
     3202      !! CALL lbc_lnk( zenv(:,:), 'T', td_nam%i_perio, 1._dp, 'no0' )  
    32153203 
    32163204      ! smooth the bathymetry (if required) 
     
    32263214      zrfact = ( 1._dp - td_nam%d_rmax ) / ( 1._dp + td_nam%d_rmax )       
    32273215 
    3228       ! initialise temporary evelope depth arrays 
     3216      ! initialise temporary envelope depth arrays 
    32293217      ALLOCATE(ztmpi1(jpi,jpj)) 
    32303218      ALLOCATE(ztmpi2(jpi,jpj)) 
     
    32363224      ztmpj1(:,:) = zenv(:,:) 
    32373225      ztmpj2(:,:) = zenv(:,:) 
    3238        
     3226 
    32393227      ! initialise temporary r-value arrays 
    32403228      ALLOCATE(zri(jpi,jpj)) 
     
    32663254               IF( (zenv(ji  ,jj) > 0._dp) .AND. & 
    32673255                 & (zenv(iip1,jj) > 0._dp) )THEN 
    3268                   zri(ji,jj) = ( zenv(iip1,jj) - zenv(ji,jj) ) / ( zenv(iip1,jj) + zenv(ji,jj) ) 
     3256                  zri(ji,jj) = ( zenv(iip1,  jj) - zenv(ji,jj) ) / ( zenv(iip1,jj  ) + zenv(ji,jj) ) 
    32693257               END IF 
    32703258               IF( (zenv(ji,jj  ) > 0._dp) .AND. & 
    32713259                 & (zenv(ji,ijp1) > 0._dp) )THEN 
    3272                   zrj(ji,jj) = ( zenv(ji,ijp1) - zenv(ji,jj) ) / ( zenv(ji,ijp1) + zenv(ji,jj) ) 
     3260                  zrj(ji,jj) = ( zenv(ji  ,ijp1) - zenv(ji,jj) ) / ( zenv(ji  ,ijp1) + zenv(ji,jj) ) 
    32733261               END IF 
    32743262               IF( zri(ji,jj) >  td_nam%d_rmax ) ztmpi1(ji  ,jj  ) = zenv(iip1,jj  ) * zrfact 
     
    32933281         END DO 
    32943282         ! apply lateral boundary condition   CAUTION: keep the value when the lbc field is zero 
    3295          CALL lbc_lnk( zenv, 'T', td_nam%i_perio, 1._dp) ! ?? , 'no0' ) 
     3283         ! this is only in mpp case, so here just do nothing 
     3284         !!CALL lbc_lnk( zenv, 'T', td_nam%i_perio, 1._dp, 'no0' ) 
    32963285      ENDDO 
    32973286      !     End loop     ! 
     
    33453334      tg_hbatv%d_value(:,:,1,1) = td_nam%d_sbot_min 
    33463335      tg_hbatf%d_value(:,:,1,1) = td_nam%d_sbot_min 
    3347        
     3336 
    33483337      DO jj = 1, jpj-1 
    33493338        DO ji = 1, jpi-1   ! NO vector opt. 
     
    33633352          DO ji = 1, jpi 
    33643353            IF( ABS(tg_hbatt%d_value(ji,jj,1,1)) < td_nam%d_wdmin1 )THEN 
    3365                tg_hbatt%d_value(ji,jj,1,1) = SIGN(1._wp, tg_hbatt%d_value(ji,jj,1,1)) * td_nam%d_wdmin1 
     3354               tg_hbatt%d_value(ji,jj,1,1) = SIGN(1._dp, tg_hbatt%d_value(ji,jj,1,1)) * td_nam%d_wdmin1 
    33663355            ENDIF 
    33673356            IF( ABS(tg_hbatu%d_value(ji,jj,1,1)) < td_nam%d_wdmin1 )THEN 
    3368                tg_hbatu%d_value(ji,jj,1,1) = SIGN(1._wp, tg_hbatu%d_value(ji,jj,1,1)) * td_nam%d_wdmin1 
     3357               tg_hbatu%d_value(ji,jj,1,1) = SIGN(1._dp, tg_hbatu%d_value(ji,jj,1,1)) * td_nam%d_wdmin1 
    33693358            ENDIF 
    33703359            IF( ABS(tg_hbatv%d_value(ji,jj,1,1)) < td_nam%d_wdmin1 )THEN 
    3371                tg_hbatv%d_value(ji,jj,1,1) = SIGN(1._wp, tg_hbatv%d_value(ji,jj,1,1)) * td_nam%d_wdmin1 
     3360               tg_hbatv%d_value(ji,jj,1,1) = SIGN(1._dp, tg_hbatv%d_value(ji,jj,1,1)) * td_nam%d_wdmin1 
    33723361            ENDIF 
    33733362            IF( ABS(tg_hbatf%d_value(ji,jj,1,1)) < td_nam%d_wdmin1 )THEN 
    3374                tg_hbatf%d_value(ji,jj,1,1) = SIGN(1._wp, tg_hbatf%d_value(ji,jj,1,1)) * td_nam%d_wdmin1 
     3363               tg_hbatf%d_value(ji,jj,1,1) = SIGN(1._dp, tg_hbatf%d_value(ji,jj,1,1)) * td_nam%d_wdmin1 
    33753364            ENDIF 
    33763365          END DO 
     
    33863375      DO jj = 1, jpj 
    33873376         DO ji = 1, jpi 
    3388             IF( tg_hbatu%d_value(ji,jj,1,1) == 0._dp ) THEN 
     3377            IF( tg_hbatu%d_value(ji,jj,1,1) == 0._dp )THEN 
    33893378               !No worries about the following line when l_wd == .true. 
    33903379               IF( zhbat(ji,jj) == 0._dp ) tg_hbatu%d_value(ji,jj,1,1) = td_nam%d_sbot_min 
     
    34193408 
    34203409!!bug:  key_helsinki a verifer 
    3421       IF( .NOT.td_nam%l_wd ) THEN 
     3410      IF( .NOT.td_nam%l_wd )THEN 
    34223411         dl_hift(:,:) = MIN( dl_hift(:,:), tg_hbatt%d_value(:,:,1,1) ) 
    34233412         dl_hifu(:,:) = MIN( dl_hifu(:,:), tg_hbatu%d_value(:,:,1,1) ) 
    34243413         dl_hifv(:,:) = MIN( dl_hifv(:,:), tg_hbatv%d_value(:,:,1,1) ) 
    34253414         dl_hiff(:,:) = MIN( dl_hiff(:,:), tg_hbatf%d_value(:,:,1,1) ) 
    3426       ENDIf 
     3415      ENDIF 
    34273416 
    34283417      CALL logger_info(' MAX val hif   t '//TRIM(fct_str(MAXVAL( dl_hift(:,:) )))//& 
     
    34643453         CALL grid_zgr__sco_s_tanh( td_nam,jpi,jpj,jpk, & 
    34653454         &                          dl_scosrf, & 
    3466          &                          dl_hift, dl_hifu, dl_hifv)!, dl_hiff ) 
     3455         &                          dl_hift, dl_hifu, dl_hifv, dl_hiff ) 
    34673456      ENDIF 
    34683457 
     
    34773466      CALL lbc_lnk( tg_e3u_0%d_value(:,:,:,1) , 'U', td_nam%i_perio, 1._dp ) 
    34783467      CALL lbc_lnk( tg_e3v_0%d_value(:,:,:,1) , 'V', td_nam%i_perio, 1._dp ) 
    3479       !CALL lbc_lnk( tg_e3f_0%d_value(:,:,:,1) , 'F', td_nam%i_perio, 1._dp ) 
     3468      CALL lbc_lnk( tg_e3f_0%d_value(:,:,:,1) , 'F', td_nam%i_perio, 1._dp ) 
    34803469      CALL lbc_lnk( tg_e3w_0%d_value(:,:,:,1) , 'W', td_nam%i_perio, 1._dp ) 
    3481       !CALL lbc_lnk( tg_e3uw_0%d_value(:,:,:,1), 'U', td_nam%i_perio, 1._dp ) 
    3482       !CALL lbc_lnk( tg_e3vw_0%d_value(:,:,:,1), 'V', td_nam%i_perio, 1._dp ) 
     3470      CALL lbc_lnk( tg_e3uw_0%d_value(:,:,:,1), 'U', td_nam%i_perio, 1._dp ) 
     3471      CALL lbc_lnk( tg_e3vw_0%d_value(:,:,:,1), 'V', td_nam%i_perio, 1._dp ) 
    34833472 
    34843473      IF( .NOT. td_nam%l_wd ) THEN 
    3485          WHERE( tg_e3t_0%d_value(:,:,:,1) ==0_dp )  tg_e3t_0%d_value(:,:,:,1) = 1._dp 
    3486          WHERE( tg_e3u_0%d_value(:,:,:,1) ==0_dp )  tg_e3u_0%d_value(:,:,:,1) = 1._dp 
    3487          WHERE( tg_e3v_0%d_value(:,:,:,1) ==0_dp )  tg_e3v_0%d_value(:,:,:,1) = 1._dp 
    3488          !WHERE( tg_e3f_0%d_value(:,:,:,1) ==0_dp )  tg_e3f_0%d_value(:,:,:,1) = 1._dp 
    3489          WHERE( tg_e3w_0%d_value(:,:,:,1) ==0_dp )  tg_e3w_0%d_value(:,:,:,1) = 1._dp 
    3490          !WHERE( tg_e3uw_0%d_value(:,:,:,1)==0_dp )  tg_e3uw_0%d_value(:,:,:,1)= 1._dp 
    3491          !WHERE( tg_e3vw_0%d_value(:,:,:,1)==0_dp )  tg_e3vw_0%d_value(:,:,:,1)= 1._dp 
     3474         WHERE( tg_e3t_0%d_value(:,:,:,1) == 0_dp )  tg_e3t_0%d_value(:,:,:,1) = 1._dp 
     3475         WHERE( tg_e3u_0%d_value(:,:,:,1) == 0_dp )  tg_e3u_0%d_value(:,:,:,1) = 1._dp 
     3476         WHERE( tg_e3v_0%d_value(:,:,:,1) == 0_dp )  tg_e3v_0%d_value(:,:,:,1) = 1._dp 
     3477         WHERE( tg_e3f_0%d_value(:,:,:,1) == 0_dp )  tg_e3f_0%d_value(:,:,:,1) = 1._dp 
     3478         WHERE( tg_e3w_0%d_value(:,:,:,1) == 0_dp )  tg_e3w_0%d_value(:,:,:,1) = 1._dp 
     3479         WHERE( tg_e3uw_0%d_value(:,:,:,1)== 0_dp )  tg_e3uw_0%d_value(:,:,:,1)= 1._dp 
     3480         WHERE( tg_e3vw_0%d_value(:,:,:,1)== 0_dp )  tg_e3vw_0%d_value(:,:,:,1)= 1._dp 
    34923481      ENDIF 
    34933482 
     
    35293518 
    35303519      CALL logger_info(' MIN val e3    t '//TRIM(fct_str(MINVAL( tg_e3t_0%d_value(:,:,:,1) )))//& 
    3531          !&                           ' f '//TRIM(fct_str(MINVAL( tg_e3f_0%d_value(:,:,:,1) )))//& 
     3520         &                           ' f '//TRIM(fct_str(MINVAL( tg_e3f_0%d_value(:,:,:,1) )))//& 
    35323521         &                           ' u '//TRIM(fct_str(MINVAL( tg_e3u_0%d_value(:,:,:,1) )))//& 
    35333522         &                           ' v '//TRIM(fct_str(MINVAL( tg_e3v_0%d_value(:,:,:,1) )))//& 
    3534          !&                          ' uw '//TRIM(fct_str(MINVAL( tg_e3uw_0%d_value(:,:,:,1) )))//& 
    3535          !&                          ' vw '//TRIM(fct_str(MINVAL( tg_e3vw_0%d_value(:,:,:,1) )))//& 
     3523         &                          ' uw '//TRIM(fct_str(MINVAL( tg_e3uw_0%d_value(:,:,:,1) )))//& 
     3524         &                          ' vw '//TRIM(fct_str(MINVAL( tg_e3vw_0%d_value(:,:,:,1) )))//& 
    35363525         &                           ' w '//TRIM(fct_str(MINVAL( tg_e3w_0%d_value(:,:,:,1) ))) ) 
    35373526 
     
    35413530 
    35423531      CALL logger_info(' MAX val e3    t '//TRIM(fct_str(MAXVAL( tg_e3t_0%d_value(:,:,:,1) )))//& 
    3543          !&                           ' f '//TRIM(fct_str(MAXVAL( tg_e3f_0%d_value(:,:,:,1) )))//& 
     3532         &                           ' f '//TRIM(fct_str(MAXVAL( tg_e3f_0%d_value(:,:,:,1) )))//& 
    35443533         &                           ' u '//TRIM(fct_str(MAXVAL( tg_e3u_0%d_value(:,:,:,1) )))//& 
    35453534         &                           ' v '//TRIM(fct_str(MAXVAL( tg_e3v_0%d_value(:,:,:,1) )))//& 
    3546          !&                          ' uw '//TRIM(fct_str(MAXVAL( tg_e3uw_0%d_value(:,:,:,1) )))//& 
    3547          !&                          ' vw '//TRIM(fct_str(MAXVAL( tg_e3vw_0%d_value(:,:,:,1) )))//& 
     3535         &                          ' uw '//TRIM(fct_str(MAXVAL( tg_e3uw_0%d_value(:,:,:,1) )))//& 
     3536         &                          ' vw '//TRIM(fct_str(MAXVAL( tg_e3vw_0%d_value(:,:,:,1) )))//& 
    35483537         &                           ' w '//TRIM(fct_str(MAXVAL( tg_e3w_0%d_value(:,:,:,1) ))) ) 
    35493538 
     
    35543543         DO jj = 1, jpj 
    35553544 
    3556             IF( tg_hbatt%d_value(ji,jj,1,1) > 0._dp) THEN 
     3545            IF( tg_hbatt%d_value(ji,jj,1,1) > 0._dp )THEN 
    35573546               DO jk = 1, INT(tg_mbathy%d_value(ji,jj,1,1),i4) 
    35583547                  ! check coordinate is monotonically increasing 
     
    35963585                  ENDIF 
    35973586               ENDDO 
    3598  
    35993587            ENDIF 
    36003588 
     
    36373625      REAL(dp) :: zcoefw 
    36383626 
    3639       REAL(wp) ::   ztmpu 
    3640       REAL(wp) ::   ztmpv 
    3641       REAL(wp) ::   ztmpf 
    3642       REAL(wp) ::   ztmpu1 
    3643       REAL(wp) ::   ztmpv1 
    3644       REAL(wp) ::   ztmpf1 
     3627      REAL(dp) ::   ztmpu 
     3628      REAL(dp) ::   ztmpv 
     3629      REAL(dp) ::   ztmpf 
     3630      REAL(dp) ::   ztmpu1 
     3631      REAL(dp) ::   ztmpv1 
     3632      REAL(dp) ::   ztmpf1 
    36453633 
    36463634      REAL(dp), DIMENSION(:,:,:), ALLOCATABLE :: z_gsigw3 
     
    37473735 
    37483736               IF( td_nam%l_wd .AND. & 
    3749                  & ( ztmpu1 < 0._wp .OR. ABS(ztmpu) < td_nam%d_wdmin1 ) )THEN 
    3750                   z_esigtu3(ji,jj,jk) = 0.5_wp * ( z_esigt3(ji,jj,jk) + z_esigt3(ji+1,jj,jk) ) 
    3751                   z_esigwu3(ji,jj,jk) = 0.5_wp * ( z_esigw3(ji,jj,jk) + z_esigw3(ji+1,jj,jk) ) 
     3737                 & ( ztmpu1 < 0._dp .OR. ABS(ztmpu) < td_nam%d_wdmin1 ) )THEN 
     3738                  z_esigtu3(ji,jj,jk) = 0.5_dp * ( z_esigt3(ji,jj,jk) + z_esigt3(ji+1,jj,jk) ) 
     3739                  z_esigwu3(ji,jj,jk) = 0.5_dp * ( z_esigw3(ji,jj,jk) + z_esigw3(ji+1,jj,jk) ) 
    37523740               ELSE 
    37533741                  z_esigtu3(ji,jj,jk) = ( tg_hbatt%d_value(ji  ,jj,1,1)*z_esigt3(ji  ,jj,jk) & 
     
    37583746 
    37593747               IF( td_nam%l_wd .AND. & 
    3760                  & ( ztmpv1 < 0._wp .OR. ABS(ztmpv) < td_nam%d_wdmin1 ) )THEN 
    3761                   z_esigtv3(ji,jj,jk) = 0.5_wp * ( z_esigt3(ji,jj,jk) + z_esigt3(ji,jj+1,jk) ) 
    3762                   z_esigwv3(ji,jj,jk) = 0.5_wp * ( z_esigw3(ji,jj,jk) + z_esigw3(ji,jj+1,jk) ) 
     3748                 & ( ztmpv1 < 0._dp .OR. ABS(ztmpv) < td_nam%d_wdmin1 ) )THEN 
     3749                  z_esigtv3(ji,jj,jk) = 0.5_dp * ( z_esigt3(ji,jj,jk) + z_esigt3(ji,jj+1,jk) ) 
     3750                  z_esigwv3(ji,jj,jk) = 0.5_dp * ( z_esigw3(ji,jj,jk) + z_esigw3(ji,jj+1,jk) ) 
    37633751               ELSE 
    37643752                  z_esigtv3(ji,jj,jk) = ( tg_hbatt%d_value(ji,jj  ,1,1)*z_esigt3(ji,jj  ,jk) & 
     
    37693757 
    37703758               IF( td_nam%l_wd .AND. & 
    3771                  & ( ztmpf1 < 0._wp .OR. ABS(ztmpf) < td_nam%d_wdmin1 ) )THEN 
    3772                   z_esigtf3(ji,jj,jk) = 0.25_wp * ( z_esigt3(ji  ,jj  ,jk) & 
     3759                 & ( ztmpf1 < 0._dp .OR. ABS(ztmpf) < td_nam%d_wdmin1 ) )THEN 
     3760                  z_esigtf3(ji,jj,jk) = 0.25_dp * ( z_esigt3(ji  ,jj  ,jk) & 
    37733761                                      &           + z_esigt3(ji+1,jj  ,jk) & 
    37743762                                      &           + z_esigt3(ji  ,jj+1,jk) & 
     
    37873775               tg_e3v_0%d_value(ji,jj,jk,1) = ( ( tg_hbatv%d_value(ji,jj,1,1) - td_nam%d_hc )*z_esigtv3(ji,jj,jk) & 
    37883776                  &                            + td_nam%d_hc / REAL(jpk-1,dp) ) 
    3789                !tg_e3f_0%d_value(ji,jj,jk,1) = ( ( tg_hbatf%d_value(ji,jj,1,1) - td_nam%d_hc ) *z_esigtf3(ji,jj,jk) & 
    3790                !   &                            + td_nam%d_hc/REAL(jpk-1,dp) ) 
    3791                ! 
     3777               tg_e3f_0%d_value(ji,jj,jk,1) = ( ( tg_hbatf%d_value(ji,jj,1,1) - td_nam%d_hc ) *z_esigtf3(ji,jj,jk) & 
     3778                  &                            + td_nam%d_hc/REAL(jpk-1,dp) ) 
     3779                
    37923780               tg_e3w_0%d_value (ji,jj,jk,1)= ( ( tg_hbatt%d_value(ji,jj,1,1) - td_nam%d_hc )*z_esigw3 (ji,jj,jk) & 
    37933781                  &                            + td_nam%d_hc / REAL(jpk-1,dp) ) 
    3794                !tg_e3uw_0%d_value(ji,jj,jk,1)= ( ( tg_hbatu%d_value(ji,jj,1,1) - td_nam%d_hc)*z_esigwu3(ji,jj,jk) & 
    3795                !  &                             + td_nam%d_hc/REAL(jpk-1,dp) ) 
    3796                !tg_e3vw_0%d_value(ji,jj,jk,1)= ( ( tg_hbatv%d_value(ji,jj,1,1) - td_nam%d_hc)*z_esigwv3(ji,jj,jk) & 
    3797                !  &                             + td_nam%d_hc/REAL(jpk-1,dp) ) 
     3782               tg_e3uw_0%d_value(ji,jj,jk,1)= ( ( tg_hbatu%d_value(ji,jj,1,1) - td_nam%d_hc)*z_esigwu3(ji,jj,jk) & 
     3783                 &                             + td_nam%d_hc/REAL(jpk-1,dp) ) 
     3784               tg_e3vw_0%d_value(ji,jj,jk,1)= ( ( tg_hbatv%d_value(ji,jj,1,1) - td_nam%d_hc)*z_esigwv3(ji,jj,jk) & 
     3785                 &                             + td_nam%d_hc/REAL(jpk-1,dp) ) 
    37983786            END DO 
    37993787        END DO 
     
    38503838      REAL(dp) ::   zzs, zzb  ! Surface and bottom cell thickness in sigma space 
    38513839 
    3852       REAL(wp) ::   ztmpu 
    3853       REAL(wp) ::   ztmpv 
    3854       REAL(wp) ::   ztmpf 
    3855       REAL(wp) ::   ztmpu1 
    3856       REAL(wp) ::   ztmpv1 
    3857       REAL(wp) ::   ztmpf1 
     3840      REAL(dp) ::   ztmpu 
     3841      REAL(dp) ::   ztmpv 
     3842      REAL(dp) ::   ztmpf 
     3843      REAL(dp) ::   ztmpu1 
     3844      REAL(dp) ::   ztmpv1 
     3845      REAL(dp) ::   ztmpf1 
    38583846 
    38593847      REAL(dp), DIMENSION(:,:,:), ALLOCATABLE :: z_gsigw3 
     
    39093897               
    39103898              IF( td_nam%d_efold /= 0.0_dp )THEN 
    3911                 zsmth   = tanh( (tg_hbatt%d_value(ji,jj,1,1)- td_nam%d_hc ) / td_nam%d_efold ) 
     3899                zsmth = TANH( (tg_hbatt%d_value(ji,jj,1,1)-td_nam%d_hc ) / td_nam%d_efold ) 
    39123900              ELSE 
    39133901                zsmth = 1.0_dp  
     
    39323920            DO jk = 1, jpk 
    39333921              z_gsigw3(ji,jj,jk) =  REAL(jk-1,dp)        /REAL(jpk-1,dp)*(td_nam%d_hc/tg_hbatt%d_value(ji,jj,1,1)) 
    3934               z_gsigt3(ji,jj,jk) = (REAL(jk-1,dp)+0.5_wp)/REAL(jpk-1,dp)*(td_nam%d_hc/tg_hbatt%d_value(ji,jj,1,1)) 
     3922              z_gsigt3(ji,jj,jk) = (REAL(jk-1,dp)+0.5_dp)/REAL(jpk-1,dp)*(td_nam%d_hc/tg_hbatt%d_value(ji,jj,1,1)) 
    39353923            END DO 
    39363924 
     
    39513939 
    39523940          DO jk = 1, jpk 
    3953              tg_gdept_0%d_value (ji,jj,jk,1) = (dd_scosrf(ji,jj)+tg_hbatt%d_value(ji,jj,1,1))*z_gsigt3(ji,jj,jk) 
    3954              tg_gdepw_0%d_value (ji,jj,jk,1) = (dd_scosrf(ji,jj)+tg_hbatt%d_value(ji,jj,1,1))*z_gsigw3(ji,jj,jk) 
     3941             tg_gdept_0%d_value(ji,jj,jk,1) = (dd_scosrf(ji,jj)+tg_hbatt%d_value(ji,jj,1,1))*z_gsigt3(ji,jj,jk) 
     3942             tg_gdepw_0%d_value(ji,jj,jk,1) = (dd_scosrf(ji,jj)+tg_hbatt%d_value(ji,jj,1,1))*z_gsigw3(ji,jj,jk) 
    39553943             !tg_gdep3w_0%d_value(ji,jj,jk,1) = (dd_scosrf(ji,jj)+tg_hbatt%d_value(ji,jj,1,1))*z_gsi3w3(ji,jj,jk) 
    39563944          END DO 
     
    39783966 
    39793967            IF( td_nam%l_wd .AND. & 
    3980               & ( ztmpu1 < 0._wp .OR. ABS(ztmpu) < td_nam%d_wdmin1 ) )THEN 
    3981                z_esigtu3(ji,jj,jk) = 0.5_wp * ( z_esigt3(ji,jj,jk) + z_esigt3(ji+1,jj,jk) ) 
    3982                z_esigwu3(ji,jj,jk) = 0.5_wp * ( z_esigw3(ji,jj,jk) + z_esigw3(ji+1,jj,jk) ) 
     3968              & ( ztmpu1 < 0._dp .OR. ABS(ztmpu) < td_nam%d_wdmin1 ) )THEN 
     3969               z_esigtu3(ji,jj,jk) = 0.5_dp * ( z_esigt3(ji,jj,jk) + z_esigt3(ji+1,jj,jk) ) 
     3970               z_esigwu3(ji,jj,jk) = 0.5_dp * ( z_esigw3(ji,jj,jk) + z_esigw3(ji+1,jj,jk) ) 
    39833971            ELSE 
    39843972               z_esigtu3(ji,jj,jk) = ( tg_hbatt%d_value(ji  ,jj,1,1)*z_esigt3(ji  ,jj,jk) & 
    3985                                    & + tg_hbatt%d_value(ji+1,jj,1,1)*z_esigt3(ji+1,jj,jk)) / ztmpu 
     3973                                   & + tg_hbatt%d_value(ji+1,jj,1,1)*z_esigt3(ji+1,jj,jk) ) / ztmpu 
    39863974               z_esigwu3(ji,jj,jk) = ( tg_hbatt%d_value(ji  ,jj,1,1)*z_esigw3(ji  ,jj,jk) & 
    3987                                    & + tg_hbatt%d_value(ji+1,jj,1,1)*z_esigw3(ji+1,jj,jk)) / ztmpu 
     3975                                   & + tg_hbatt%d_value(ji+1,jj,1,1)*z_esigw3(ji+1,jj,jk) ) / ztmpu 
    39883976            ENDIF 
    39893977 
    39903978            IF( td_nam%l_wd .AND. & 
    3991               & ( ztmpv1 < 0._wp .OR. ABS(ztmpv) < td_nam%d_wdmin1 ) )THEN 
    3992                z_esigtv3(ji,jj,jk) = 0.5_wp * ( z_esigt3(ji,jj,jk) + z_esigt3(ji,jj+1,jk) ) 
    3993                z_esigwv3(ji,jj,jk) = 0.5_wp * ( z_esigw3(ji,jj,jk) + z_esigw3(ji,jj+1,jk) ) 
     3979              & ( ztmpv1 < 0._dp .OR. ABS(ztmpv) < td_nam%d_wdmin1 ) )THEN 
     3980               z_esigtv3(ji,jj,jk) = 0.5_dp * ( z_esigt3(ji,jj,jk) + z_esigt3(ji,jj+1,jk) ) 
     3981               z_esigwv3(ji,jj,jk) = 0.5_dp * ( z_esigw3(ji,jj,jk) + z_esigw3(ji,jj+1,jk) ) 
    39943982            ELSE 
    39953983               z_esigtv3(ji,jj,jk) = ( tg_hbatt%d_value(ji,jj  ,1,1)*z_esigt3(ji,jj  ,jk) & 
     
    40003988 
    40013989            IF( td_nam%l_wd .AND. & 
    4002               & ( ztmpf1 < 0._wp .OR. ABS(ztmpf) < td_nam%d_wdmin1 ) )THEN 
    4003                z_esigtf3(ji,jj,jk) = 0.25_wp * ( z_esigt3(ji  ,jj  ,jk) & 
     3990              & ( ztmpf1 < 0._dp .OR. ABS(ztmpf) < td_nam%d_wdmin1 ) )THEN 
     3991               z_esigtf3(ji,jj,jk) = 0.25_dp * ( z_esigt3(ji  ,jj  ,jk) & 
    40043992                                   &           + z_esigt3(ji+1,jj  ,jk) & 
    40053993                                   &           + z_esigt3(ji  ,jj+1,jk) & 
     
    40154003             tg_e3u_0%d_value(ji,jj,jk,1) = (dd_scosrf(ji,jj)+tg_hbatu%d_value(ji,jj,1,1))*z_esigtu3(ji,jj,jk) 
    40164004             tg_e3v_0%d_value(ji,jj,jk,1) = (dd_scosrf(ji,jj)+tg_hbatv%d_value(ji,jj,1,1))*z_esigtv3(ji,jj,jk) 
    4017              !tg_e3f_0%d_value(ji,jj,jk,1) = (dd_scosrf(ji,jj)+tg_hbatf%d_value(ji,jj,1,1))*z_esigtf3(ji,jj,jk) 
     4005             tg_e3f_0%d_value(ji,jj,jk,1) = (dd_scosrf(ji,jj)+tg_hbatf%d_value(ji,jj,1,1))*z_esigtf3(ji,jj,jk) 
    40184006             ! 
    4019              tg_e3w_0%d_value (ji,jj,jk,1) = tg_hbatt%d_value(ji,jj,1,1)*z_esigw3 (ji,jj,jk) 
    4020              !tg_e3uw_0%d_value(ji,jj,jk,1) = tg_hbatu%d_value(ji,jj,1,1)*z_esigwu3(ji,jj,jk) 
    4021              !tg_e3vw_0%d_value(ji,jj,jk,1) = tg_hbatv%d_value(ji,jj,1,1)*z_esigwv3(ji,jj,jk) 
     4007             tg_e3w_0%d_value(ji,jj,jk,1) = tg_hbatt%d_value(ji,jj,1,1)*z_esigw3 (ji,jj,jk) 
     4008             tg_e3uw_0%d_value(ji,jj,jk,1) = tg_hbatu%d_value(ji,jj,1,1)*z_esigwu3(ji,jj,jk) 
     4009             tg_e3vw_0%d_value(ji,jj,jk,1) = tg_hbatv%d_value(ji,jj,1,1)*z_esigwv3(ji,jj,jk) 
    40224010          END DO 
    40234011 
     
    40284016      CALL lbc_lnk(tg_e3u_0 %d_value(:,:,:,1),'T', td_nam%i_perio, 1._dp) 
    40294017      CALL lbc_lnk(tg_e3v_0 %d_value(:,:,:,1),'T', td_nam%i_perio, 1._dp)  
    4030       !CALL lbc_lnk(tg_e3f_0 %d_value(:,:,:,1),'T', td_nam%i_perio, 1._dp) 
     4018      CALL lbc_lnk(tg_e3f_0 %d_value(:,:,:,1),'T', td_nam%i_perio, 1._dp) 
    40314019      CALL lbc_lnk(tg_e3w_0 %d_value(:,:,:,1),'T', td_nam%i_perio, 1._dp) 
    4032       !CALL lbc_lnk(tg_e3uw_0%d_value(:,:,:,1),'T', td_nam%i_perio, 1._dp)  
    4033       !CALL lbc_lnk(tg_e3vw_0%d_value(:,:,:,1),'T', td_nam%i_perio, 1._dp) 
     4020      CALL lbc_lnk(tg_e3uw_0%d_value(:,:,:,1),'T', td_nam%i_perio, 1._dp)  
     4021      CALL lbc_lnk(tg_e3vw_0%d_value(:,:,:,1),'T', td_nam%i_perio, 1._dp) 
    40344022 
    40354023      DEALLOCATE( z_gsigw3  ) 
     
    40684056   SUBROUTINE grid_zgr__sco_s_tanh( td_nam,jpi,jpj,jpk, & 
    40694057         &                          dd_scosrf, & 
    4070          &                          dd_hift, dd_hifu, dd_hifv)!, dd_hiff )  
     4058         &                          dd_hift, dd_hifu, dd_hifv, dd_hiff )  
    40714059      IMPLICIT NONE 
    40724060      ! Argument       
     
    40794067      REAL(dp)   , DIMENSION(:,:), INTENT(IN   ) :: dd_hifu 
    40804068      REAL(dp)   , DIMENSION(:,:), INTENT(IN   ) :: dd_hifv 
    4081 !      REAL(dp)   , DIMENSION(:,:), INTENT(IN   ) :: dd_hiff 
     4069      REAL(dp)   , DIMENSION(:,:), INTENT(IN   ) :: dd_hiff 
    40824070 
    40834071      ! local variable 
     
    41404128              tg_e3v_0%d_value(ji,jj,jk,1) = ( (tg_hbatv%d_value(ji,jj,1 ,1) - dd_hifv(ji,jj)) & 
    41414129                 &                            * tg_esigt%d_value(1 ,1 ,jk,1) + dd_hifv(ji,jj)/REAL(jpk-1,dp) ) 
    4142               !tg_e3f_0%d_value(ji,jj,jk,1) = ( (tg_hbatf%d_value(ji,jj,1 ,1) - dd_hiff(ji,jj)) & 
    4143               !   &                            * tg_esigt%d_value(1 ,1, jk,1) + dd_hiff(ji,jj)/REAL(jpk-1,dp) ) 
     4130              tg_e3f_0%d_value(ji,jj,jk,1) = ( (tg_hbatf%d_value(ji,jj,1 ,1) - dd_hiff(ji,jj)) & 
     4131                 &                            * tg_esigt%d_value(1 ,1, jk,1) + dd_hiff(ji,jj)/REAL(jpk-1,dp) ) 
    41444132                
    41454133              tg_e3w_0%d_value (ji,jj,jk,1)= ( (tg_hbatt%d_value(ji,jj,1 ,1) - dd_hift(ji,jj)) & 
    41464134                 &                            * tg_esigw%d_value(1 ,1 ,jk,1) + dd_hift(ji,jj)/REAL(jpk-1,dp) ) 
    4147               !tg_e3uw_0%d_value(ji,jj,jk,1)= ( (tg_hbatu%d_value(ji,jj,1 ,1) - dd_hifu(ji,jj)) & 
    4148               !   &                            * tg_esigw%d_value(1 ,1 ,jk,1) + dd_hifu(ji,jj)/REAL(jpk-1,dp) ) 
    4149               !tg_e3vw_0%d_value(ji,jj,jk,1)= ( (tg_hbatv%d_value(ji,jj,1 ,1) - dd_hifv(ji,jj)) & 
    4150               !   &                            * tg_esigw%d_value(1 ,1 ,jk,1) + dd_hifv(ji,jj)/REAL(jpk-1,dp) ) 
     4135              tg_e3uw_0%d_value(ji,jj,jk,1)= ( (tg_hbatu%d_value(ji,jj,1 ,1) - dd_hifu(ji,jj)) & 
     4136                 &                            * tg_esigw%d_value(1 ,1 ,jk,1) + dd_hifu(ji,jj)/REAL(jpk-1,dp) ) 
     4137              tg_e3vw_0%d_value(ji,jj,jk,1)= ( (tg_hbatv%d_value(ji,jj,1 ,1) - dd_hifv(ji,jj)) & 
     4138                 &                            * tg_esigw%d_value(1 ,1 ,jk,1) + dd_hifv(ji,jj)/REAL(jpk-1,dp) ) 
    41514139            END DO 
    41524140         END DO 
     
    41794167      !!---------------------------------------------------------------------- 
    41804168      ! 
    4181       pf =   (   TANH( td_nam%d_theta * ( -(pk-0.5_dp) / REAL(jpk-1) + td_nam%d_thetb )  )  & 
    4182          &     - TANH( td_nam%d_thetb * td_nam%d_theta                                )  )  & 
    4183          & * (   COSH( td_nam%d_theta                           )                      & 
    4184          &     + COSH( td_nam%d_theta * ( 2._dp * td_nam%d_thetb - 1._dp ) )  )              & 
     4169      pf =   (   TANH( td_nam%d_theta * ( -(pk-0.5_dp) / REAL(jpk-1,dp) + td_nam%d_thetb )  ) & 
     4170         &     - TANH( td_nam%d_thetb * td_nam%d_theta                                   )  ) & 
     4171         & * (   COSH( td_nam%d_theta                                      )                  & 
     4172         &     + COSH( td_nam%d_theta * ( 2._dp * td_nam%d_thetb - 1._dp ) )  )               & 
    41854173         & / ( 2._dp * SINH( td_nam%d_theta ) ) 
    41864174      ! 
     
    42154203      ! 
    42164204      IF ( td_nam%d_theta == 0 ) then      ! uniform sigma 
    4217          pf1 = - ( pk1 - 0.5_dp ) / REAL( jpk-1 ) 
     4205         pf1 = - ( pk1 - 0.5_dp ) / REAL(jpk-1,dp) 
    42184206      ELSE                        ! stretched sigma 
    4219          pf1 =   ( 1._dp - pbb ) * ( SINH( td_nam%d_theta*(-(pk1-0.5_dp)/REAL(jpk-1)) ) ) / SINH( td_nam%d_theta )              & 
    4220             &  + pbb * (  (TANH( td_nam%d_theta*( (-(pk1-0.5_dp)/REAL(jpk-1)) + 0.5_dp) ) - TANH( 0.5_dp * td_nam%d_theta )  )  & 
     4207         pf1 =   ( 1._dp - pbb ) * ( SINH( td_nam%d_theta*(-(pk1-0.5_dp)/REAL(jpk-1,dp)) ) ) / SINH( td_nam%d_theta )              & 
     4208            &  + pbb * (  (TANH( td_nam%d_theta*( (-(pk1-0.5_dp)/REAL(jpk-1,dp)) + 0.5_dp) ) - TANH( 0.5_dp * td_nam%d_theta )  )  & 
    42214209            &        / ( 2._dp * TANH( 0.5_dp * td_nam%d_theta ) )  ) 
    42224210      ENDIF 
     
    42534241      TYPE(TNAMZ), INTENT(IN   ) :: td_nam 
    42544242      INTEGER(i4), INTENT(IN   ) :: jpk 
    4255       REAL(dp)   , INTENT(IN   ) :: pk1(jpk)       ! continuous "k" coordinate 
    4256       REAL(dp)                   :: p_gamma(jpk)   ! stretched coordinate 
     4243      REAL(dp)   , DIMENSION(:) , INTENT(IN   ) :: pk1       ! continuous "k" coordinate 
     4244      ! function 
     4245      REAL(dp)   , DIMENSION(jpk) :: p_gamma   ! stretched coordinate 
     4246      ! local variable 
    42574247      REAL(dp)   , INTENT(IN   ) :: pzb           ! Bottom box depth 
    42584248      REAL(dp)   , INTENT(IN   ) :: pzs           ! surface box depth 
     
    42664256      !!---------------------------------------------------------------------- 
    42674257      ! 
    4268       zn1  =  1./(jpk-1.) 
    4269       zn2  =  1. -  zn1 
     4258      zn1  =  1._dp / REAL(jpk-1,dp) 
     4259      zn2  =  1._dp -  zn1 
    42704260 
    42714261      za1 = (td_nam%d_alpha+2.0_dp)*zn1**(td_nam%d_alpha+1.0_dp)-(td_nam%d_alpha+1.0_dp)*zn1**(td_nam%d_alpha+2.0_dp)  
     
    43634353                  &                     + tg_gdepw_0%d_value(ji  ,jj  ,jk+1,1)& 
    43644354                  &                     - tg_gdepw_0%d_value(ji  ,jj-1,jk+1,1) ) & 
    4365                   &                 / ( tg_gdepw_0%d_value  (ji  ,jj  ,jk  ,1)& 
     4355                  &                   / ( tg_gdepw_0%d_value(ji  ,jj  ,jk  ,1)& 
    43664356                  &                     + tg_gdepw_0%d_value(ji  ,jj-1,jk  ,1)& 
    43674357                  &                     - tg_gdepw_0%d_value(ji  ,jj  ,jk+1,1)& 
  • branches/2016/dev_r6999_CONFIGMAN_1/NEMOGCM/TOOLS/SIREN/src/iom_cdf.f90

    r7153 r7233  
    21032103                  & 'gcost','gcosu','gcosv','gcosf', & 
    21042104                  & 'gsint','gsinu','gsinv','gsinf', & 
    2105                   & 'mbathy','misf','isf_draft',     & 
     2105                  & 'mbathy','misf','isf_draft','isfdraft',     & 
    21062106                  & 'hbatt','hbatu','hbatv','hbatf', & 
    21072107                  & 'gsigt','gsigu','gsigv','gsigf', & 
     
    22122212         CALL iom_cdf__check(il_status,"IOM CDF WRITE VAR DEF: ") 
    22132213      ENDIF 
     2214      CALL logger_debug("IOM CDF WRITE VAR DEF: type = "//TRIM(fct_str(tl_var%i_type))) 
    22142215 
    22152216      ! remove unuseful attribute 
  • branches/2016/dev_r6999_CONFIGMAN_1/NEMOGCM/TOOLS/SIREN/src/mpp.f90

    r7025 r7233  
    25502550   !> and the number of processors following I and J. 
    25512551   !> Then the number of sea/land processors is compute with mask 
    2552    ! 
     2552   !> 
    25532553   !> @author J.Paul 
    25542554   !> @date October, 2015 - Initial version 
    2555    ! 
     2555   !> @date October, 2016 
     2556   !> - compare index to td_lay number of proc instead of td_mpp (bug fix) 
     2557   !> 
    25562558   !> @param[in] td_mpp mpp strcuture 
    25572559   !> @param[in] id_mask   sub domain mask (sea=1, land=0) 
     
    26702672 
    26712673            ! east boundary 
    2672             IF( ji == td_mpp%i_niproc )THEN 
     2674            IF( ji == td_lay%i_niproc )THEN 
    26732675               il_lei = td_lay%i_lci(ji,jj) 
    26742676            ELSE 
     
    26772679 
    26782680            ! north boundary 
    2679             IF( jj == td_mpp%i_njproc )THEN 
     2681            IF( jj == td_lay%i_njproc )THEN 
    26802682               il_lej = td_lay%i_lcj(ji,jj) 
    26812683            ELSE 
  • branches/2016/dev_r6999_CONFIGMAN_1/NEMOGCM/TOOLS/SIREN/src/phycst.f90

    r7025 r7233  
    3434   PUBLIC :: dp_siday   !< sideral day                        [s] 
    3535 
    36    REAL(wp), PUBLIC ::   rday = 24.*60.*60.     !: day                                [s] 
    37    REAL(wp), PUBLIC ::   rsiyea                 !: sideral year                       [s] 
    38    REAL(wp), PUBLIC ::   rsiday                 !: sideral day                        [s] 
     36   REAL(dp), PUBLIC ::   rday = 24.*60.*60.     !: day                                [s] 
     37   REAL(dp), PUBLIC ::   rsiyea                 !: sideral year                       [s] 
     38   REAL(dp), PUBLIC ::   rsiday                 !: sideral day                        [s] 
    3939 
    4040   REAL(dp), PARAMETER :: dp_pi = 3.14159274101257_dp 
    41    REAL(dp), PARAMETER :: dp_eps = EPSILON(1._dp) 
     41   REAL(dp), PARAMETER :: dp_eps = 0.5 * EPSILON(1._dp) 
    4242   REAL(dp), PARAMETER :: dp_rearth = 6371229._dp 
    4343   REAL(dp), PARAMETER :: dp_deg2rad = dp_pi/180.0 
     
    4545 
    4646   REAL(dp), PARAMETER :: dp_day = 24.*60.*60.      
    47    REAL(dp), PARAMETER :: dp_siyea = 365.25_wp * dp_day * & 
    48       &  2._wp * dp_pi / 6.283076_dp 
    49    REAL(dp), PARAMETER :: dp_siday = dp_day / ( 1._wp + dp_day / dp_siyea ) 
     47   REAL(dp), PARAMETER :: dp_siyea = 365.25_dp * dp_day * & 
     48      &  2._dp * dp_pi / 6.283076_dp 
     49   REAL(dp), PARAMETER :: dp_siday = dp_day / ( 1._dp + dp_day / dp_siyea ) 
    5050 
    5151   REAL(dp), PARAMETER :: dp_delta=1.e-5 
  • branches/2016/dev_r6999_CONFIGMAN_1/NEMOGCM/TOOLS/SIREN/src/variable.f90

    r7153 r7233  
    11391139   !> @date July, 2015 
    11401140   !> - add unit factor (to change unit) 
     1141   !> @date November, 2016 
     1142   !> - allow to add scalar value 
    11411143   !> 
    11421144   !> @param[in] cd_name         variable name 
     
    12921294         dl_value(1,1,1,:) = dd_value(:) 
    12931295      ELSE 
    1294          CALL logger_fatal("VAR INIT: can not add value from variable "//& 
    1295          &  TRIM(cd_name)//". invalid dimension to be used") 
     1296         IF( SIZE(dd_value(:)) > 1 )THEN 
     1297            CALL logger_fatal("VAR INIT: can not add value from variable "//& 
     1298            &  TRIM(cd_name)//". invalid dimension to be used") 
     1299         ELSE 
     1300            dl_value(1,1,1,1) = dd_value(1) 
     1301            CALL logger_warn("VAR INIT: add scalar value for variable "//& 
     1302            &  TRIM(cd_name)) 
     1303 
     1304         ENDIF 
    12961305      ENDIF 
    12971306 
Note: See TracChangeset for help on using the changeset viewer.