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 3837 for branches/2011/DEV_r2739_STFC_dCSE/NEMOGCM/NEMO/OPA_SRC/DOM/domzgr.F90 – NEMO

Ignore:
Timestamp:
2013-03-12T15:55:32+01:00 (11 years ago)
Author:
trackstand2
Message:

Merge of finiss

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2011/DEV_r2739_STFC_dCSE/NEMOGCM/NEMO/OPA_SRC/DOM/domzgr.F90

    r3432 r3837  
    4343 
    4444   PUBLIC   dom_zgr        ! called by dom_init.F90 
     45   PUBLIC   zgr_z, zgr_bat, zgr_zco, zgr_zps ! called by nemogcm::recursive_partition 
     46   PUBLIC   fssig1         ! called by partition_mod::smooth_bathy 
    4547 
    4648   !                                       !!* Namelist namzgr_sco * 
     
    5456   !                                        ! ( rn_bb=0; top only, rn_bb =1; top and bottom) 
    5557   REAL(wp) ::   rn_hc       =  150._wp     ! Critical depth for s-sigma coordinates 
    56    PUBLIC rn_sbot_min, rn_sbot_max, rn_theta, rn_thetb, rn_rmax, ln_s_sigma, rn_bb,rn_hc 
     58   PUBLIC rn_sbot_min, rn_sbot_max, rn_theta, rn_thetb, rn_rmax, & 
     59          ln_s_sigma, rn_bb, rn_hc 
     60   PUBLIC ln_zco, ln_zps, ln_sco 
     61 
    5762   !! * Control permutation of array indices 
    5863#  include "oce_ftrans.h90" 
     
    6267#  include "domzgr_substitute.h90" 
    6368#  include "vectopt_loop_substitute.h90" 
     69 
     70   NAMELIST/namzgr/ ln_zco, ln_zps, ln_sco 
     71   NAMELIST/namzgr_sco/ rn_sbot_max, rn_sbot_min, rn_theta, rn_thetb, & 
     72                        rn_rmax, ln_s_sigma, rn_bb, rn_hc 
     73   PUBLIC namzgr, namzgr_sco 
    6474   !!---------------------------------------------------------------------- 
    6575   !! NEMO/OPA 3.3.1 , NEMO Consortium (2011) 
     
    8898      INTEGER ::   ioptio = 0   ! temporary integer 
    8999      ! 
    90       NAMELIST/namzgr/ ln_zco, ln_zps, ln_sco 
    91       !!---------------------------------------------------------------------- 
    92  
    93       REWIND( numnam )                 ! Read Namelist namzgr : vertical coordinate' 
     100      !NAMELIST/namzgr/ ln_zco, ln_zps, ln_sco 
     101      !!---------------------------------------------------------------------- 
     102 
     103      REWIND( numnam )                 ! Read Namelist namzgr : vertical coordinate 
    94104      READ  ( numnam, namzgr ) 
    95105 
     
    287297 
    288298 
    289    SUBROUTINE zgr_bat 
     299   SUBROUTINE zgr_bat(global_domain) 
    290300      !!---------------------------------------------------------------------- 
    291301      !!                    ***  ROUTINE zgr_bat  *** 
     
    317327      !!              - bathy : meter bathymetry (in meters) 
    318328      !!---------------------------------------------------------------------- 
     329      LOGICAL, OPTIONAL, INTENT(in) :: global_domain ! Whether dealing with 
     330                                                     ! whole domain (T) or a  
     331                                                     ! sub-domain after domain 
     332                                                     ! decomposition 
     333      ! Locals 
    319334      INTEGER  ::   ji, jj, jl, jk            ! dummy loop indices 
    320335      INTEGER  ::   inum                      ! temporary logical unit 
     
    325340      INTEGER , DIMENSION(jpidta,jpjdta) ::   idta   ! global domain integer data 
    326341      REAL(wp), DIMENSION(jpidta,jpjdta) ::   zdta   ! global domain scalar data 
     342      LOGICAL  ::   is_global  
    327343      !!---------------------------------------------------------------------- 
    328344 
     
    330346      IF(lwp) WRITE(numout,*) '    zgr_bat : defines level and meter bathymetry' 
    331347      IF(lwp) WRITE(numout,*) '    ~~~~~~~' 
     348 
     349      ! Set local flag to signal whether we're dealing with the global domain 
     350      ! (pre decomposition) or a local part of it. Required by the  
     351      ! recursive k-section partitioning. 
     352      is_global = .FALSE. 
     353      IF( PRESENT(global_domain) )THEN 
     354         IF( global_domain ) is_global = .TRUE. 
     355      END IF 
    332356 
    333357      !                                               ! ================== !  
     
    347371            ii_bump = jpidta / 2                           ! i-index of the bump center 
    348372            ij_bump = jpjdta / 2                           ! j-index of the bump center 
    349             r_bump  = 50000._wp                            ! bump radius (meters)        
    350             h_bump  =  2700._wp                            ! bump height (meters) 
     373            r_bump  = 0.165*MIN(jpidta,jpjdta)             ! bump radius (grid cells)        
     374            h_bump  =  3000._wp                            ! bump height (meters) 
    351375            h_oce   = gdepw_0(jpk)                         ! background ocean depth (meters) 
    352376            IF(lwp) WRITE(numout,*) '            bump characteristics: ' 
    353             IF(lwp) WRITE(numout,*) '               bump center (i,j)   = ', ii_bump, ii_bump 
     377            IF(lwp) WRITE(numout,*) '               bump center (i,j)   = ', ii_bump, ij_bump 
    354378            IF(lwp) WRITE(numout,*) '               bump height         = ', h_bump , ' meters' 
    355             IF(lwp) WRITE(numout,*) '               bump radius         = ', r_bump , ' index' 
     379            IF(lwp) WRITE(numout,*) '               bump radius         = ', r_bump , ' cells' 
    356380            IF(lwp) WRITE(numout,*) '            background ocean depth = ', h_oce  , ' meters' 
    357381            !                                         
    358382            DO jj = 1, jpjdta                              ! zdta : 
    359383               DO ji = 1, jpidta 
    360                   zi = FLOAT( ji - ii_bump ) * ppe1_m / r_bump 
    361                   zj = FLOAT( jj - ij_bump ) * ppe2_m / r_bump 
     384                  !zi = FLOAT( ji - ii_bump ) * ppe1_m / r_bump 
     385                  !zj = FLOAT( jj - ij_bump ) * ppe2_m / r_bump 
     386                  zi = FLOAT( ji - ii_bump ) / r_bump 
     387                  zj = FLOAT( jj - ij_bump ) / r_bump 
    362388                  zdta(ji,jj) = h_oce - h_bump * EXP( -( zi*zi + zj*zj ) ) 
    363389               END DO 
     
    412438         ! 
    413439         IF( ln_zco )   THEN                          ! zco : read level bathymetry  
    414             CALL iom_open ( 'bathy_level.nc', inum )   
    415             CALL iom_get  ( inum, jpdom_data, 'Bathy_level', bathy ) 
     440            CALL iom_open ( 'bathy_level.nc', inum )  
     441            IF(is_global)THEN  
     442               CALL iom_get ( inum, jpdom_unknown, 'Bathymetry' , bathy, & 
     443                              kstart=(/jpizoom,jpjzoom/),                & 
     444                              kcount=(/jpiglo,jpjglo/) ) 
     445            ELSE 
     446               CALL iom_get  ( inum, jpdom_data, 'Bathy_level', bathy ) 
     447            END IF 
     448 
    416449            CALL iom_close( inum ) 
    417450            mbathy(:,:) = INT( bathy(:,:) ) 
     
    446479         IF( ln_zps .OR. ln_sco )   THEN              ! zps or sco : read meter bathymetry 
    447480            CALL iom_open ( 'bathy_meter.nc', inum )  
    448             CALL iom_get  ( inum, jpdom_data, 'Bathymetry', bathy ) 
     481            IF(is_global)THEN  
     482               CALL iom_get ( inum, jpdom_unknown, 'Bathymetry' , bathy, & 
     483                              kstart=(/jpizoom,jpjzoom/),                & 
     484                              kcount=(/jpiglo,jpjglo/) ) 
     485            ELSE 
     486               CALL iom_get  ( inum, jpdom_data, 'Bathymetry', bathy ) 
     487            END IF 
    449488            CALL iom_close( inum ) 
    450489            !                                                ! ===================== 
     
    516555         zhmin = gdepw_0(ik+1)                                                         ! minimum depth = ik+1 w-levels  
    517556         WHERE( bathy(:,:) <= 0._wp )   ;   bathy(:,:) = 0._wp                         ! min=0     over the lands 
    518          ELSE WHERE                     ;   bathy(:,:) = MAX(  zhmin , bathy(:,:)  )   ! min=zhmin over the oceans 
     557         ELSEWHERE                      ;   bathy(:,:) = MAX(  zhmin , bathy(:,:)  )   ! min=zhmin over the oceans 
    519558         END WHERE 
    520559         IF(lwp) write(numout,*) 'Minimum ocean depth: ', zhmin, ' minimum number of ocean levels : ', ik 
     
    730769      !!                                     (min value = 1 over land) 
    731770      !!---------------------------------------------------------------------- 
     771      !USE arpdebugging, ONLY: dump_array 
    732772      USE wrk_nemo, ONLY:   wrk_in_use, wrk_not_released 
    733773      USE wrk_nemo, ONLY:   zmbk => wrk_2d_1 
     
    744784      IF(lwp) WRITE(numout,*) '    ~~~~~~~~~~~~~' 
    745785      ! 
     786      !CALL dump_array(0, 'mbathy', mbathy, withHalos=.TRUE.) 
     787 
    746788      mbkt(:,:) = MAX( mbathy(:,:) , 1 )    ! bottom k-index of T-level (=1 over land) 
    747789      !                                     ! bottom k-index of W-level = mbkt+1 
     
    755797      zmbk(:,:) = REAL( mbku(:,:), wp )   ;   CALL lbc_lnk(zmbk,'U',1.)   ;   mbku  (:,:) = MAX( INT( zmbk(:,:) ), 1 ) 
    756798      zmbk(:,:) = REAL( mbkv(:,:), wp )   ;   CALL lbc_lnk(zmbk,'V',1.)   ;   mbkv  (:,:) = MAX( INT( zmbk(:,:) ), 1 ) 
     799      ! 
     800      ! Compute and store the deepest bottom level of any grid-type at each grid point 
     801      ! For use in removing data below ocean floor from halo exchanges. 
     802      mbkmax(:,:) = MAX(mbkt(:,:)+1, mbku(:,:), mbkv(:,:)) 
    757803      ! 
    758804      IF( wrk_not_released(2, 1) )   CALL ctl_stop('zgr_bot_level: failed to release workspace array') 
     
    809855 
    810856 
    811    SUBROUTINE zgr_zps 
     857   SUBROUTINE zgr_zps(pre_domain_decomp) 
    812858      !!---------------------------------------------------------------------- 
    813859      !!                  ***  ROUTINE zgr_zps  *** 
     
    858904      !! DCSE_NEMO: wrk_3d_1 renamed, need additional directive 
    859905!FTRANS zprt :I :I :z 
     906      LOGICAL, INTENT(in), OPTIONAL :: pre_domain_decomp 
    860907      !! 
    861908      INTEGER  ::   ji, jj, jk       ! dummy loop indices 
     
    892939      bathy(:,:) = MIN( zmax ,  bathy(:,:) )    ! bounded value of bathy (min already set at the end of zgr_bat) 
    893940      WHERE( bathy(:,:) == 0._wp )   ;   mbathy(:,:) = 0       ! land  : set mbathy to 0 
    894       ELSE WHERE                     ;   mbathy(:,:) = jpkm1   ! ocean : initialize mbathy to the max ocean level 
     941      ELSEWHERE                     ;   mbathy(:,:) = jpkm1   ! ocean : initialize mbathy to the max ocean level 
    895942      END WHERE 
    896943 
     
    903950         WHERE( 0._wp < bathy(:,:) .AND. bathy(:,:) <= zdepth )   mbathy(:,:) = jk-1 
    904951      END DO 
     952 
     953      ! If we've been called before domain decomposition then we only want to compute 
     954      ! mbathy and the return. 
     955      IF( PRESENT(pre_domain_decomp) )THEN 
     956         IF( pre_domain_decomp )RETURN 
     957      ENDIF 
    905958 
    906959      ! Scale factors and depth at T- and W-points 
     
    12081261      USE mapcomm_mod, ONLY: trimmed, cyclic_bc 
    12091262      USE mapcomm_mod, ONLY: nidx, eidx, sidx, widx 
    1210 !      USE arpdebugging, ONLY: dump_array 
     1263 
    12111264      !! DCSE_NEMO: wrk_nemo module variables renamed, need additional directives 
    12121265!FTRANS gsigw3 :I :I :z 
     
    12271280      ! 
    12281281 
    1229       NAMELIST/namzgr_sco/ rn_sbot_max, rn_sbot_min, rn_theta, rn_thetb, rn_rmax, ln_s_sigma, rn_bb, rn_hc 
     1282!      NAMELIST/namzgr_sco/ rn_sbot_max, rn_sbot_min, rn_theta, rn_thetb, rn_rmax, ln_s_sigma, rn_bb, rn_hc 
    12301283      !!---------------------------------------------------------------------- 
    12311284 
     
    13191372                                                         jl, zrmax, INT( SUM(zmsk(:,:) ) ) 
    13201373         ! 
    1321 !!$         IF(jl < 6)THEN ! .OR. (MOD(jl,1000) == 0) )THEN 
    1322 !!$            CALL dump_array(jl, 'zenv_before', zenv, withHalos=.TRUE.) 
    1323 !!$            CALL dump_array(jl, 'ztmp_before', ztmp, withHalos=.TRUE.) 
    1324 !!$            CALL dump_array(jl, 'zmsk_before', zmsk, withHalos=.TRUE.) 
    1325 !!$         END IF 
    13261374 
    13271375         ! Copy current surface before next smoothing iteration  
     
    13611409         ! Apply lateral boundary condition but do not zero on closed boundaries 
    13621410         CALL lbc_lnk( zenv, 'T', 1._wp, lzero=.FALSE. ) 
    1363  
    1364 !!$         IF(jl < 6)THEN ! .OR. (MOD(jl,1000) == 0) )THEN 
    1365 !!$            CALL dump_array(jl, 'zenv', zenv, withHalos=.TRUE.) 
    1366 !!$            CALL dump_array(jl, 'ztmp', ztmp, withHalos=.TRUE.) 
    1367 !!$            CALL dump_array(jl, 'zmsk', zmsk, withHalos=.TRUE.) 
    1368 !!$         END IF 
    13691411 
    13701412         !                                                  ! ================ ! 
     
    17091751                  CALL ctl_stop( ctmp1 ) 
    17101752               ENDIF 
     1753#if defined key_vvl 
    17111754               IF( gdepw_1(ji,jj,jk) < 0._wp .OR. gdept_1(ji,jj,jk) < 0._wp ) THEN 
    17121755                  WRITE(ctmp1,*) 'zgr_sco :   gdepw or gdept =< 0  at point (i,j,k)= ', ji, jj, jk 
    17131756                  CALL ctl_stop( ctmp1 ) 
    17141757               ENDIF 
     1758#endif 
    17151759            END DO 
    17161760         END DO 
Note: See TracChangeset for help on using the changeset viewer.