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 9078 for branches/UKMO – NEMO

Changeset 9078 for branches/UKMO


Ignore:
Timestamp:
2017-12-15T15:31:35+01:00 (6 years ago)
Author:
davestorkey
Message:

UKMO/dev_r8600_closea_rewrite branch : slight redesign after discussion with Gurvan:

  1. Mask fields included in domain_cfg file.
  2. Control variables changed.
  3. Suppression of closed seas now done in dom_zgr.
Location:
branches/UKMO/dev_r8600_closea_rewrite/NEMOGCM
Files:
10 edited

Legend:

Unmodified
Added
Removed

  • branches/UKMO/dev_r8600_closea_rewrite/NEMOGCM/CONFIG/SHARED/namelist_ref

    r8995 r9078  
    8080   ln_use_jattr = .false.  !  use (T) the file attribute: open_ocean_jstart, if present 
    8181   !                       !  in netcdf input files, as the start j-row for reading 
     82   ln_closea = .true.      !  T => keep closed seas (defined by closea_mask field) in the domain and apply 
     83                           !       special treatment of freshwater fluxes. 
     84                           !  F => suppress closed seas (defined by closea_mask field) from the bathymetry 
     85                           !       at runtime. 
     86                           !  If there is no closea_mask field in the domain_cfg file or we do not use 
     87                           !  a domain_cfg file then this logical does nothing. 
    8288/ 
    8389!----------------------------------------------------------------------- 
    … …  
    8591!----------------------------------------------------------------------- 
    8692   ln_linssh   = .false.   !  =T  linear free surface  ==>>  model level are fixed in time 
    87    nn_closea   =    0      !  remove (=0) or keep (=1) closed seas and lakes (ORCA) 
    88    ! 
    8993   nn_msh      =    0      !  create (>0) a mesh file or not (=0) 
    9094   rn_isfhmin  =    1.00   !  treshold (m) to discriminate grounding ice to floating ice 

  • branches/UKMO/dev_r8600_closea_rewrite/NEMOGCM/NEMO/OPA_SRC/DOM/closea.F90

    r9021 r9078  
    1414 
    1515   !!---------------------------------------------------------------------- 
    16    !!   dom_clo    : modification of the ocean domain for closed seas cases 
    17    !!   sbc_clo    : Special handling of closed seas 
     16   !!   dom_clo    : read in masks which define closed seas and runoff areas 
     17   !!   sbc_clo    : Special handling of freshwater fluxes over closed seas 
    1818   !!   clo_rnf    : set close sea outflows as river mouths (see sbcrnf) 
    19    !!   clo_bat    : set to zero a field over closed sea (see domzrg) 
     19   !!   clo_bat    : set to zero a field over closed sea (see domzgr) 
    2020   !!---------------------------------------------------------------------- 
    2121   USE oce             ! dynamics and tracers 
    … …  
    3838   PUBLIC sbc_clo      ! called by sbcmod module 
    3939   PUBLIC clo_rnf      ! called by sbcrnf module 
    40    PUBLIC clo_bat      ! called in domain module 
    41  
    42    INTEGER, PUBLIC :: jncs !: number of closed seas (inferred from closea_mask field) 
    43    INTEGER, PUBLIC :: jncsr !: number of closed seas rnf mappings (inferred from closea_mask_rnf field) 
    44    INTEGER, PUBLIC :: jncse !: number of closed seas empmr mappings (inferred from closea_mask_empmr field) 
     40   PUBLIC clo_bat      ! called in domzgr module 
     41 
     42   LOGICAL, PUBLIC :: ln_closea  !:  T => keep closed seas (defined by closea_mask field) in the domain and apply 
     43                                 !:       special treatment of freshwater fluxes. 
     44                                 !:  F => suppress closed seas (defined by closea_mask field) from the bathymetry 
     45                                 !:       at runtime. 
     46                                 !:  If there is no closea_mask field in the domain_cfg file or we do not use 
     47                                 !:  a domain_cfg file then this logical does nothing. 
     48                                 !: 
     49   LOGICAL, PUBLIC :: l_sbc_clo  !: T => Closed seas defined, apply special treatment of freshwater fluxes. 
     50                                 !: F => No closed seas defined (closea_mask field not found). 
     51   LOGICAL, PUBLIC :: l_clo_rnf  !: T => Some closed seas output freshwater (RNF or EMPMR) to specified runoff points. 
     52   INTEGER, PUBLIC :: jncs       !: number of closed seas (inferred from closea_mask field) 
     53   INTEGER, PUBLIC :: jncsr      !: number of closed seas rnf mappings (inferred from closea_mask_rnf field) 
     54   INTEGER, PUBLIC :: jncse      !: number of closed seas empmr mappings (inferred from closea_mask_empmr field) 
    4555    
    4656   INTEGER, PUBLIC, ALLOCATABLE, DIMENSION(:,:) ::  closea_mask       !: mask of integers defining closed seas 
    … …  
    6171CONTAINS 
    6272 
    63    SUBROUTINE dom_clo( k_bot ) 
     73   SUBROUTINE dom_clo() 
    6474      !!--------------------------------------------------------------------- 
    6575      !!                  ***  ROUTINE dom_clo  *** 
    … …  
    7080      !!                just the thermodynamic processes are applied. 
    7181      !! 
    72       !! ** Action  :   Read closea mask fields from closea_masks.nc and infer 
     82      !! ** Action  :   Read closea_mask* fields (if they exist) from domain_cfg file and infer 
    7383      !!                number of closed seas from closea_mask field. 
    7484      !!                closea_mask       : integer values defining closed seas (or groups of closed seas) 
    … …  
    7888      !!                                    closed seas to a runoff area for net fluxes. 
    7989      !! 
    80       !!                Python code to generate the closea_masks.nc file from the old-style indices 
     90      !!                Python code to generate the closea_masks* fields from the old-style indices 
    8191      !!                definitions is available at TOOLS/DOMAINcfg/make_closea_masks.py 
    8292      !!---------------------------------------------------------------------- 
    83       INTEGER, DIMENSION(:,:), INTENT(in)  ::   k_bot   ! ocean last level index (for masking input fields)  
    8493      INTEGER ::   inum    ! input file identifier 
    8594      INTEGER ::   ierr    ! error code 
    8695      INTEGER ::   id      ! netcdf variable ID 
     96 
    8797      REAL(wp), POINTER, DIMENSION(:,:) :: zdata_in ! temporary real array for input 
    8898      !!---------------------------------------------------------------------- 
    8999      ! 
    90100      IF(lwp) WRITE(numout,*) 
    91       IF(lwp) WRITE(numout,*)'dom_clo : closed seas ' 
     101      IF(lwp) WRITE(numout,*)'dom_clo : read in masks to define closed seas ' 
    92102      IF(lwp) WRITE(numout,*)'~~~~~~~' 
    93103      ! 
    94104      CALL wrk_alloc( jpi,jpj, zdata_in) 
    95105      !  
    96       ! read the closed seas masks 
    97       ! -------------------------- 
    98       ! 
    99       ALLOCATE( closea_mask(jpi,jpj) , STAT=ierr ) 
    100       IF( ierr /= 0 )   CALL ctl_stop( 'STOP', 'dom_clo: failed to allocate closea_mask array') 
    101       ALLOCATE( closea_mask_rnf(jpi,jpj) , STAT=ierr ) 
    102       IF( ierr /= 0 )   CALL ctl_stop( 'STOP', 'dom_clo: failed to allocate closea_mask_rnf array') 
    103       ALLOCATE( closea_mask_empmr(jpi,jpj) , STAT=ierr ) 
    104       IF( ierr /= 0 )   CALL ctl_stop( 'STOP', 'dom_clo: failed to allocate closea_mask_empmr array') 
    105  
    106       closea_mask(:,:) = 0 
    107       closea_mask_rnf(:,:) = 0 
    108       closea_mask_empmr(:,:) = 0 
    109  
    110       CALL iom_open( 'closea_masks.nc', inum ) 
    111  
    112       zdata_in(:,:) = 0.0 
    113       CALL iom_get ( inum, jpdom_data, 'closea_mask', zdata_in ) 
    114       closea_mask(:,:) = NINT(zdata_in(:,:)) 
    115       ! necessary because fill values in input fields can confuse things 
    116       ! we can't multiply by tmask here because it isn't defined yet 
    117       WHERE( k_bot(:,:) == 0 ) closea_mask(:,:) = 0 
    118  
    119       id = iom_varid(inum, 'closea_mask_rnf', ldstop = .false.) 
    120       IF(lwp) WRITE(numout,*) 'closea_mask_rnf, id : ',id 
     106      ! read the closed seas masks (if they exist) from domain_cfg file 
     107      ! --------------------------------------------------------------- 
     108      ! 
     109      CALL iom_open( cn_domcfg, inum ) 
     110 
     111      id = iom_varid(inum, 'closea_mask', ldstop = .false.) 
    121112      IF( id > 0 ) THEN  
    122          CALL iom_get ( inum, jpdom_data, 'closea_mask_rnf', zdata_in ) 
    123          closea_mask_rnf(:,:) = NINT(zdata_in(:,:)) 
    124          WHERE( k_bot(:,:) == 0 ) closea_mask_rnf(:,:) = 0 
    125       ENDIF 
     113         l_sbc_clo = .true. 
     114         ALLOCATE( closea_mask(jpi,jpj) , STAT=ierr ) 
     115         IF( ierr /= 0 )   CALL ctl_stop( 'STOP', 'dom_clo: failed to allocate closea_mask array') 
     116         zdata_in(:,:) = 0.0 
     117         CALL iom_get ( inum, jpdom_data, 'closea_mask', zdata_in ) 
     118         closea_mask(:,:) = NINT(zdata_in(:,:)) * tmask(:,:,1) 
     119         ! number of closed seas = global maximum value in closea_mask field 
     120         jncs = maxval(closea_mask(:,:)) 
     121         IF( lk_mpp ) CALL mpp_max(jncs) 
     122         IF( jncs > 0 ) THEN 
     123            IF( lwp ) WRITE(numout,*) 'Number of closed seas : ',jncs 
     124         ELSE 
     125            CALL ctl_stop( 'Problem with closea_mask field in domain_cfg file. Has no values > 0 so no closed seas defined.') 
     126         ENDIF 
     127      ELSE  
     128         IF( lwp ) WRITE(numout,*) 'closea_mask field not found in domain_cfg file. No closed seas defined.' 
     129         l_sbc_clo = .false. 
     130         jncs = 0  
     131      ENDIF 
     132 
     133      l_clo_rnf = .false. 
     134 
     135      IF( l_sbc_clo ) THEN ! No point reading in closea_mask_rnf or closea_mask_empmr fields if no closed seas defined. 
     136 
     137         id = iom_varid(inum, 'closea_mask_rnf', ldstop = .false.) 
     138         IF( id > 0 ) THEN  
     139            l_clo_rnf = .true.             
     140            ALLOCATE( closea_mask_rnf(jpi,jpj) , STAT=ierr ) 
     141            IF( ierr /= 0 )   CALL ctl_stop( 'STOP', 'dom_clo: failed to allocate closea_mask_rnf array') 
     142            CALL iom_get ( inum, jpdom_data, 'closea_mask_rnf', zdata_in ) 
     143            closea_mask_rnf(:,:) = NINT(zdata_in(:,:)) * tmask(:,:,1) 
     144            ! number of closed seas rnf mappings = global maximum in closea_mask_rnf field 
     145            jncsr = maxval(closea_mask_rnf(:,:)) 
     146            IF( lk_mpp ) CALL mpp_max(jncsr) 
     147            IF( jncsr > 0 ) THEN 
     148               IF( lwp ) WRITE(numout,*) 'Number of closed seas rnf mappings : ',jncsr 
     149            ELSE 
     150               CALL ctl_stop( 'Problem with closea_mask_rnf field in domain_cfg file. Has no values > 0 so no closed seas rnf mappings defined.') 
     151            ENDIF 
     152         ELSE  
     153            IF( lwp ) WRITE(numout,*) 'closea_mask_rnf field not found in domain_cfg file. No closed seas rnf mappings defined.' 
     154            jncsr = 0 
     155         ENDIF 
    126156  
    127       id = iom_varid(inum, 'closea_mask_empmr', ldstop = .false.) 
    128       IF( id > 0 ) THEN  
    129          CALL iom_get ( inum, jpdom_data, 'closea_mask_empmr', zdata_in ) 
    130          closea_mask_empmr(:,:) = NINT(zdata_in(:,:)) 
    131          WHERE( k_bot(:,:) == 0 ) closea_mask_empmr(:,:) = 0 
    132       ENDIF 
    133  
     157         id = iom_varid(inum, 'closea_mask_empmr', ldstop = .false.) 
     158         IF( id > 0 ) THEN  
     159            l_clo_rnf = .true.             
     160            ALLOCATE( closea_mask_empmr(jpi,jpj) , STAT=ierr ) 
     161            IF( ierr /= 0 )   CALL ctl_stop( 'STOP', 'dom_clo: failed to allocate closea_mask_empmr array') 
     162            CALL iom_get ( inum, jpdom_data, 'closea_mask_empmr', zdata_in ) 
     163            closea_mask_empmr(:,:) = NINT(zdata_in(:,:)) * tmask(:,:,1) 
     164            ! number of closed seas empmr mappings = global maximum value in closea_mask_empmr field 
     165            jncse = maxval(closea_mask_empmr(:,:)) 
     166            IF( lk_mpp ) CALL mpp_max(jncse) 
     167            IF( jncse > 0 ) THEN  
     168               IF( lwp ) WRITE(numout,*) 'Number of closed seas empmr mappings : ',jncse 
     169            ELSE 
     170               CALL ctl_stop( 'Problem with closea_mask_empmr field in domain_cfg file. Has no values > 0 so no closed seas empmr mappings defined.') 
     171            ENDIF 
     172         ELSE  
     173            IF( lwp ) WRITE(numout,*) 'closea_mask_empmr field not found in domain_cfg file. No closed seas empmr mappings defined.' 
     174            jncse = 0 
     175         ENDIF 
     176 
     177      ENDIF ! l_sbc_clo 
     178      ! 
    134179      CALL iom_close( inum ) 
    135  
    136       ! number of closed seas = global maximum value in closea_mask field 
    137       jncs = maxval(closea_mask(:,:)) 
    138       IF( lk_mpp ) CALL mpp_max(jncs) 
    139       IF( lwp ) WRITE(numout,*) 'Number of closed seas : ',jncs 
    140       ! 
    141       ! number of closed seas rnf mappings = global maximum in closea_mask_rnf field 
    142       jncsr = maxval(closea_mask_rnf(:,:)) 
    143       IF( lk_mpp ) CALL mpp_max(jncsr) 
    144       IF( lwp ) WRITE(numout,*) 'Number of closed seas rnf mappings : ',jncsr 
    145       ! 
    146       ! number of closed seas empmr mappings = global maximum value in closea_mask_empmr field 
    147       jncse = maxval(closea_mask_empmr(:,:)) 
    148       IF( lk_mpp ) CALL mpp_max(jncse) 
    149       IF( lwp ) WRITE(numout,*) 'Number of closed seas empmr mappings : ',jncse 
    150180      ! 
    151181      CALL wrk_dealloc(jpi,jpj, zdata_in) 
    … …  
    217247         ! 
    218248         !                                        ! surface areas of rnf target areas 
    219          DO jcr = 1, jncsr 
    220             ztmp2d(:,:) = 0.e0_wp 
    221             WHERE( closea_mask_rnf(:,:) == jcr .and. closea_mask(:,:) == 0 ) ztmp2d(:,:) = e1e2t(:,:) * tmask_i(:,:) 
    222             surfr(jcr) = glob_sum( ztmp2d(:,:) ) 
    223          END DO 
     249         IF( jncsr > 0 ) THEN 
     250            DO jcr = 1, jncsr 
     251               ztmp2d(:,:) = 0.e0_wp 
     252               WHERE( closea_mask_rnf(:,:) == jcr .and. closea_mask(:,:) == 0 ) ztmp2d(:,:) = e1e2t(:,:) * tmask_i(:,:) 
     253               surfr(jcr) = glob_sum( ztmp2d(:,:) ) 
     254            END DO 
     255         ENDIF 
    224256         ! 
    225257         !                                        ! surface areas of empmr target areas 
    226          DO jce = 1, jncse 
    227             ztmp2d(:,:) = 0.e0_wp 
    228             WHERE( closea_mask_empmr(:,:) == jcr .and. closea_mask(:,:) == 0 ) ztmp2d(:,:) = e1e2t(:,:) * tmask_i(:,:) 
    229             surfe(jcr) = glob_sum( ztmp2d(:,:) ) 
    230          END DO 
     258         IF( jncse > 0 ) THEN 
     259            DO jce = 1, jncse 
     260               ztmp2d(:,:) = 0.e0_wp 
     261               WHERE( closea_mask_empmr(:,:) == jce .and. closea_mask(:,:) == 0 ) ztmp2d(:,:) = e1e2t(:,:) * tmask_i(:,:) 
     262               surfe(jce) = glob_sum( ztmp2d(:,:) ) 
     263            END DO 
     264         ENDIF 
    231265         ! 
    232266         IF(lwp) WRITE(numout,*)'     Closed sea surface areas (km2)' 
    … …  
    268302      zfwf_total = SUM(zfwf) 
    269303 
    270       ! 
    271       ! 2. Work out total FW fluxes over rnf source areas and add to rnf target areas. If jncsr is zero does not loop. 
    272       !    Where zfwf is negative add flux at specified runoff points and subtract from fluxes for global redistribution. 
    273       !    Where positive leave in global redistribution total. 
    274       ! 
    275304      zfwfr(:) = 0.e0_wp            
    276       DO jcr = 1, jncsr 
    277          ! 
    278          ztmp2d(:,:) = 0.e0_wp 
    279          WHERE( closea_mask_rnf(:,:) == jcr .and. closea_mask(:,:) > 0 ) ztmp2d(:,:) = e1e2t(:,:) * ( emp(:,:)-rnf(:,:) ) * tmask_i(:,:) 
    280          zfwfr(jcr) = glob_sum( ztmp2d(:,:) ) 
    281          IF(lwp) WRITE(numout,*) 'closea runoff output: jcr, zfwfr(jcr), ABS(zfwfr(jcr)/surf(jncs+1) : ',jcr,zfwfr(jcr),ABS(zfwfr(jcr)/surf(jncs+1)) 
    282          ! 
    283          ! The following if avoids the redistribution of the round off 
    284          IF ( ABS(zfwfr(jcr) / surf(jncs+1) ) > rsmall) THEN 
    285             ! 
    286             ! Add residuals to target runoff points if negative and subtract from total to be added globally 
    287             IF(lwp) WRITE(numout,*) 'closea runoff output: passed roundoff test!' 
    288             IF( zfwfr(jcr) < 0.0 ) THEN  
    289                zfwf_total = zfwf_total - zfwfr(jcr) 
    290                zcoef    = zfwfr(jcr) / surfr(jcr) 
     305      IF( jncsr > 0 ) THEN 
     306         ! 
     307         ! 2. Work out total FW fluxes over rnf source areas and add to rnf target areas.  
     308         !    Where zfwf is negative add flux at specified runoff points and subtract from fluxes for global redistribution. 
     309         !    Where positive leave in global redistribution total. 
     310         ! 
     311         DO jcr = 1, jncsr 
     312            ! 
     313            ztmp2d(:,:) = 0.e0_wp 
     314            WHERE( closea_mask_rnf(:,:) == jcr .and. closea_mask(:,:) > 0 ) ztmp2d(:,:) = e1e2t(:,:) * ( emp(:,:)-rnf(:,:) ) * tmask_i(:,:) 
     315            zfwfr(jcr) = glob_sum( ztmp2d(:,:) ) 
     316            ! 
     317            ! The following if avoids the redistribution of the round off 
     318            IF ( ABS(zfwfr(jcr) / surf(jncs+1) ) > rsmall) THEN 
     319               ! 
     320               ! Add residuals to target runoff points if negative and subtract from total to be added globally 
     321               IF( zfwfr(jcr) < 0.0 ) THEN  
     322                  zfwf_total = zfwf_total - zfwfr(jcr) 
     323                  zcoef    = zfwfr(jcr) / surfr(jcr) 
     324                  zcoef1   = rcp * zcoef 
     325                  IF(lwp) WRITE(numout,*) 'NEGATIVE FLUX zfwfr(jcr), zcoef : ',zfwfr(jcr), zcoef 
     326                  WHERE( closea_mask_rnf(:,:) == jcr .and. closea_mask(:,:) == 0.0) 
     327                     emp(:,:) = emp(:,:) + zcoef 
     328                     qns(:,:) = qns(:,:) - zcoef1 * sst_m(:,:) 
     329                  ENDWHERE 
     330               ELSE 
     331                  IF(lwp) WRITE(numout,*) 'POSITIVE FLUX zfwfr(jcr) : ',zfwfr(jcr) 
     332               ENDIF 
     333               ! 
     334            ENDIF 
     335         END DO 
     336      ENDIF  ! jncsr > 0     
     337      ! 
     338      zfwfe(:) = 0.e0_wp            
     339      IF( jncse > 0 ) THEN 
     340         ! 
     341         ! 3. Work out total fluxes over empmr source areas and add to empmr target areas. If jncse is zero does not loop.  
     342         ! 
     343         DO jce = 1, jncse 
     344            ! 
     345            ztmp2d(:,:) = 0.e0_wp 
     346            WHERE( closea_mask_empmr(:,:) == jce .and. closea_mask(:,:) > 0 ) ztmp2d(:,:) = e1e2t(:,:) * ( emp(:,:)-rnf(:,:) ) * tmask_i(:,:) 
     347            zfwfe(jce) = glob_sum( ztmp2d(:,:) ) 
     348            ! 
     349            ! The following if avoids the redistribution of the round off 
     350            IF ( ABS( zfwfe(jce) / surf(jncs+1) ) > rsmall ) THEN 
     351               ! 
     352               ! Add residuals to runoff points and subtract from total to be added globally 
     353               zfwf_total = zfwf_total - zfwfe(jce) 
     354               zcoef    = zfwfe(jce) / surfe(jce) 
    291355               zcoef1   = rcp * zcoef 
    292                IF(lwp) WRITE(numout,*) 'closea runoff output: jcr, zcoef: ',jcr, zcoef 
    293                WHERE( closea_mask_rnf(:,:) == jcr .and. closea_mask(:,:) == 0.0) 
     356               WHERE( closea_mask_empmr(:,:) == jce .and. closea_mask(:,:) == 0.0) 
    294357                  emp(:,:) = emp(:,:) + zcoef 
    295358                  qns(:,:) = qns(:,:) - zcoef1 * sst_m(:,:) 
    296359               ENDWHERE 
     360               ! 
    297361            ENDIF 
    298             ! 
    299          ENDIF 
    300       END DO 
    301      
    302       ! 
    303       ! 3. Work out total fluxes over empmr source areas and add to empmr target areas. If jncse is zero does not loop.  
    304       ! 
    305       zfwfe(:) = 0.e0_wp            
    306       DO jce = 1, jncse 
    307          ! 
    308          ztmp2d(:,:) = 0.e0_wp 
    309          WHERE( closea_mask_empmr(:,:) == jce .and. closea_mask(:,:) > 0 ) ztmp2d(:,:) = e1e2t(:,:) * ( emp(:,:)-rnf(:,:) ) * tmask_i(:,:) 
    310          zfwfe(jce) = glob_sum( ztmp2d(:,:) ) 
    311          ! 
    312          ! The following if avoids the redistribution of the round off 
    313          IF ( ABS(zfwfe(jce) / surf(jncs+1) ) > rsmall) THEN 
    314             ! 
    315             ! Add residuals to runoff points and subtract from total to be added globally 
    316             zfwf_total = zfwf_total - zfwfe(jce) 
    317             zcoef    = zfwfe(jce) / surfe(jce) 
    318             zcoef1   = rcp * zcoef 
    319             WHERE( closea_mask_empmr(:,:) == jce .and. closea_mask(:,:) == 0.0) 
    320                emp(:,:) = emp(:,:) + zcoef 
    321                qns(:,:) = qns(:,:) - zcoef1 * sst_m(:,:) 
    322             ENDWHERE 
    323             ! 
    324          ENDIF 
    325       END DO 
    326      
     362         END DO 
     363      ENDIF ! jncse > 0     
     364 
    327365      ! 
    328366      ! 4. Spread residual flux over global ocean.  
    … …  
    332370         zcoef    = zfwf_total / surf(jncs+1) 
    333371         zcoef1   = rcp * zcoef 
    334          IF(lwp) WRITE(numout,*) 'closea global addition: zfwf_total, zcoef: ', zfwf_total, zcoef 
    335372         WHERE( closea_mask(:,:) == 0 ) 
    336373            emp(:,:) = emp(:,:) + zcoef 
    … …  
    349386            zcoef    = zfwf(jc) / surf(jc) 
    350387            zcoef1   = rcp * zcoef 
    351             IF(lwp) WRITE(numout,*) 'closea subtract mean: jc, zfwf(jc), zcoef: ',jc, zfwf(jc), zcoef 
    352388            WHERE( closea_mask(:,:) == jc ) 
    353389               emp(:,:) = emp(:,:) - zcoef 
    … …  
    386422      !!---------------------------------------------------------------------- 
    387423      ! 
    388       WHERE( ( closea_mask_rnf(:,:) > 0 .or. closea_mask_empmr(:,:) > 0 ) .and. closea_mask(:,:) == 0 ) 
    389          p_rnfmsk(:,:) = MAX( p_rnfmsk(:,:), 1.0_wp ) 
    390       ENDWHERE 
     424      IF( jncsr > 0 ) THEN 
     425         WHERE( closea_mask_rnf(:,:) > 0 .and. closea_mask(:,:) == 0 ) 
     426            p_rnfmsk(:,:) = MAX( p_rnfmsk(:,:), 1.0_wp ) 
     427         ENDWHERE 
     428      ENDIF 
     429      ! 
     430      IF( jncse > 0 ) THEN 
     431         WHERE( closea_mask_empmr(:,:) > 0 .and. closea_mask(:,:) == 0 ) 
     432            p_rnfmsk(:,:) = MAX( p_rnfmsk(:,:), 1.0_wp ) 
     433         ENDWHERE 
     434      ENDIF 
    391435      ! 
    392436   END SUBROUTINE clo_rnf 
    … …  
    397441      !!                  ***  ROUTINE clo_bat  *** 
    398442      !!                     
    399       !! ** Purpose :   suppress closed sea from the domain 
    400       !! 
    401       !! ** Method  :   set first and last ocean level to 0 over the closed seas. 
    402       !! 
    403       !! ** Action  :   set pbat=0 and kbat=0 over closed seas 
     443      !! ** Purpose :   Suppress closed sea from the domain 
     444      !! 
     445      !! ** Method  :   Read in closea_mask field (if it exists) from domain_cfg file. 
     446      !!                Where closea_mask > 0 set first and last ocean level to 0 
     447      !!                (As currently coded you can't define a closea_mask field in  
     448      !!                usr_def_zgr). 
     449      !! 
     450      !! ** Action  :   set k_top=0 and k_bot=0 over closed seas 
    404451      !!---------------------------------------------------------------------- 
    405452      INTEGER, DIMENSION(:,:), INTENT(inout) ::   k_top, k_bot   ! ocean first and last level indices 
    406       !!---------------------------------------------------------------------- 
    407       ! 
    408       WHERE( closea_mask(:,:) > 0 ) 
    409          k_top(:,:) = 0    
    410          k_bot(:,:) = 0    
    411       ENDWHERE 
     453      INTEGER                           :: inum, id 
     454      INTEGER,  POINTER, DIMENSION(:,:) :: closea_mask ! closea_mask field 
     455      REAL(wp), POINTER, DIMENSION(:,:) :: zdata_in ! temporary real array for input 
     456      !!---------------------------------------------------------------------- 
     457      ! 
     458      IF(lwp) THEN                     ! Control print 
     459         WRITE(numout,*) 
     460         WRITE(numout,*) 'clo_bat : suppression of closed seas' 
     461         WRITE(numout,*) '~~~~~~~' 
     462      ENDIF 
     463      ! 
     464      CALL iom_open( cn_domcfg, inum ) 
     465      ! 
     466      id = iom_varid(inum, 'closea_mask', ldstop = .false.)       
     467      IF( id > 0 ) THEN 
     468         IF( lwp ) WRITE(numout,*) 'Suppressing closed seas in bathymetry based on closea_mask field,' 
     469         CALL wrk_alloc( jpi,jpj, zdata_in) 
     470         CALL wrk_alloc( jpi,jpj, closea_mask) 
     471         CALL iom_get ( inum, jpdom_data, 'closea_mask', zdata_in ) 
     472         closea_mask(:,:) = NINT(zdata_in(:,:)) 
     473         WHERE( closea_mask(:,:) > 0 ) 
     474            k_top(:,:) = 0    
     475            k_bot(:,:) = 0    
     476         ENDWHERE 
     477         CALL wrk_dealloc( jpi,jpj, zdata_in) 
     478         CALL wrk_dealloc( jpi,jpj, closea_mask) 
     479      ELSE 
     480         IF( lwp ) WRITE(numout,*) 'No closea_mask field found in domain_cfg file. No suppression of closed seas.' 
     481      ENDIF 
     482      ! 
     483      ! Initialise l_sbc_clo and l_clo_rnf for this case (ln_closea=.false.) 
     484      l_sbc_clo = .false. 
     485      l_clo_rnf = .false. 
     486      ! 
     487      CALL iom_close(inum) 
    412488      ! 
    413489   END SUBROUTINE clo_bat 

  • branches/UKMO/dev_r8600_closea_rewrite/NEMOGCM/NEMO/OPA_SRC/DOM/dom_oce.F90

    r9017 r9078  
    3131   !                                   !!* Namelist namdom : time & space domain * 
    3232   LOGICAL , PUBLIC ::   ln_linssh      !: =T  linear free surface ==>> model level are fixed in time 
    33    LOGICAL , PUBLIC ::   ln_closea      !: =T apply special treatment to closed seas depending on value of nn_closea. 
    34    INTEGER , PUBLIC ::   nn_closea      !: =0 suppress closed sea/lake from the ORCA domain or not (=1) 
    3533   INTEGER , PUBLIC ::   nn_msh         !: >0  create a mesh-mask file (mesh_mask.nc) 
    3634   REAL(wp), PUBLIC ::   rn_isfhmin     !: threshold to discriminate grounded ice to floating ice 

  • branches/UKMO/dev_r8600_closea_rewrite/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90

    r9017 r9078  
    130130      CALL dom_hgr                     ! Horizontal mesh 
    131131      CALL dom_zgr( ik_top, ik_bot )   ! Vertical mesh and bathymetry 
    132       IF( ln_closea ) THEN     
    133           CALL dom_clo( ik_bot )  ! Closed seas and lake  
    134           IF( nn_closea == 0 )   CALL clo_bat( ik_top, ik_bot )    !==  remove closed seas or lakes  ==! 
    135       ENDIF 
    136132      CALL dom_msk( ik_top, ik_bot )   ! Masks 
     133      IF( ln_closea ) CALL dom_clo     ! ln_closea=T : closed seas included in the simulation 
     134                                       ! Read in masks to define closed seas and lakes  
    137135      ! 
    138136      DO jj = 1, jpj                   ! depth of the iceshelves 
    … …  
    288286         &             nn_stock, nn_write , ln_mskland  , ln_clobber   , nn_chunksz, nn_euler  ,     & 
    289287         &             ln_cfmeta, ln_iscpl 
    290       NAMELIST/namdom/ ln_linssh, ln_closea, nn_closea, nn_msh, rn_isfhmin, rn_rdt, rn_atfp, ln_crs 
     288      NAMELIST/namdom/ ln_linssh, nn_msh, rn_isfhmin, rn_rdt, rn_atfp, ln_crs 
    291289#if defined key_netcdf4 
    292290      NAMELIST/namnc4/ nn_nchunks_i, nn_nchunks_j, nn_nchunks_k, ln_nc4zip 
    … …  
    397395         WRITE(numout,*) '   Namelist namdom : space & time domain' 
    398396         WRITE(numout,*) '      linear free surface (=T)              ln_linssh  = ', ln_linssh 
    399          WRITE(numout,*) '      special treatment of closed seas      ln_closea  = ', ln_closea 
    400          IF(ln_closea) WRITE(numout,*) '      suppression of closed seas (=0)       nn_closea  = ', nn_closea 
    401397         WRITE(numout,*) '      create mesh/mask file(s)              nn_msh     = ', nn_msh 
    402398         WRITE(numout,*) '           = 0   no file created           ' 

  • branches/UKMO/dev_r8600_closea_rewrite/NEMOGCM/NEMO/OPA_SRC/DOM/domzgr.F90

    r8995 r9078  
    2929   USE dom_oce        ! ocean domain 
    3030   USE usrdef_zgr     ! user defined vertical coordinate system 
     31   USE closea         ! closed seas 
    3132   USE depth_e3       ! depth <=> e3 
    3233   USE wet_dry, ONLY: ln_wd, ht_wd 
    … …  
    118119         gde3w_0(:,:,jk) = gde3w_0(:,:,jk-1) + e3w_0(:,:,jk) 
    119120      END DO 
     121      ! 
     122      ! Any closed seas (defined by closea_mask > 0 in domain_cfg file) to be filled  
     123      ! in at runtime if ln_closea=.false. 
     124      IF( .not. ln_closea ) CALL clo_bat( k_top, k_bot ) 
    120125      ! 
    121126      IF(lwp) THEN                     ! Control print 

  • branches/UKMO/dev_r8600_closea_rewrite/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90

    r9017 r9078  
    153153         WRITE(numout,*) '         runoff / runoff mouths                     ln_rnf        = ', ln_rnf 
    154154         WRITE(numout,*) '         iceshelf formulation                       ln_isf        = ', ln_isf 
    155          WRITE(numout,*) '         special treatment of closed seas (set in namdom)  ln_closea     = ', ln_closea 
    156          IF(ln_closea) WRITE(numout,*) '         closed sea treatment (=0/1) (set in namdom) nn_closea     = ', nn_closea 
    157155         WRITE(numout,*) '         nb of iterations if land-sea-mask applied  nn_lsm        = ', nn_lsm 
    158156         WRITE(numout,*) '         surface wave                               ln_wave       = ', ln_wave 
    … …  
    445443      IF( nn_fwb    /= 0 )   CALL sbc_fwb( kt, nn_fwb, nn_fsbc )  ! control the freshwater budget 
    446444 
    447       ! treatment of closed sea in the model domain   (update freshwater fluxes) 
    448       ! Should not be ran if ln_diurnal_only 
    449       IF( .NOT.ln_diurnal_only .AND. ln_closea .AND. nn_closea == 1 )   CALL sbc_clo( kt )    
     445      ! Special treatment of freshwater fluxes over closed seas in the model domain 
     446      ! Should not be run if ln_diurnal_only 
     447      IF( l_sbc_clo .AND. (.NOT. ln_diurnal_only) )   CALL sbc_clo( kt )    
    450448 
    451449!RBbug do not understand why see ticket 667 

  • branches/UKMO/dev_r8600_closea_rewrite/NEMOGCM/NEMO/OPA_SRC/SBC/sbcrnf.F90

    r9017 r9078  
    513513      CALL iom_close( inum )                                      ! close file 
    514514      ! 
    515       IF( ln_closea .and. nn_closea == 1 )   CALL clo_rnf( rnfmsk )   ! closed sea inflow set as river mouth 
     515      IF( l_clo_rnf )   CALL clo_rnf( rnfmsk )   ! closed sea inflow set as river mouth 
    516516      ! 
    517517      rnfmsk_z(:)   = 0._wp                                       ! vertical structure 

  • branches/UKMO/dev_r8600_closea_rewrite/NEMOGCM/NEMO/OPA_SRC/nemogcm.F90

    r8995 r9078  
    4848   USE phycst         ! physical constant                  (par_cst routine) 
    4949   USE domain         ! domain initialization   (dom_init & dom_cfg routines) 
     50   USE closea         ! treatment of closed seas (for ln_closea) 
    5051   USE usrdef_nam     ! user defined configuration 
    5152   USE tideini        ! tidal components initialization   (tide_ini routine) 
    … …  
    242243         &             nn_isplt , nn_jsplt, nn_jctls, nn_jctle,   & 
    243244         &             nn_timing, nn_diacfl 
    244       NAMELIST/namcfg/ ln_read_cfg, cn_domcfg, ln_write_cfg, cn_domcfg_out, ln_use_jattr 
     245      NAMELIST/namcfg/ ln_read_cfg, cn_domcfg, ln_write_cfg, cn_domcfg_out, ln_use_jattr, ln_closea 
    245246      !!---------------------------------------------------------------------- 
    246247      ! 

  • branches/UKMO/dev_r8600_closea_rewrite/NEMOGCM/TOOLS/DOMAINcfg/README

    r7828 r9078  
    1818::::::::::::::::::::::::::::::::::::::::::::::::: 
    1919NOTA: it can be run in multiproc mode, but in output there will be domain_cfg_00xx.nc files 
     20 
     21 === Closed seas (closea module) === 
     22 
     23If you want to define closed seas in the bathymetry either to suppress them at runtime or  
     24redistribute freshwater fluxes, then you need to run make_closea_masks.py after you have 
     25created the basic domain_cfg file. This utility will add "closea_mask*" fields to the 
     26domain_cfg file to define the closed seas in the configuration. (If you have closed seas 
     27but don't want to treat them in a special way then you can ignore this step).  
    2028 
    2129================================ 

  • branches/UKMO/dev_r8600_closea_rewrite/NEMOGCM/TOOLS/DOMAINcfg/make_closea_masks.py

    r9017 r9078  
    22 
    33''' 
    4 Routine to create closea mask file for new proposed closea code 
    5 at NEMO 4.0 
     4Routine to create closea mask fields based on old NEMO closea index definitions. 
     5Details of the grid and the bathymetry are read in from the domain_cfg.nc file and 
     6the closea_mask* fields are appended to the same domain_cfg.nc file.  
    67 
    78To use this routine: 
    … …  
    910  1. Provide domain_cfg.nc file for your configuration. 
    1011 
    11   2. Define closed seas with indices in the old style in Section 2.  
    12      (Read the comments on indexing!). 
    13  
    14   3. Choose whether to mask the closea_masks.nc file. Not required 
     12  2. Define closed seas for your configuration in Section 2  
     13     using indices in the old NEMO style. (Read the comments on  
     14     indexing in Section 2!). Examples are given for eORCA025  
     15     (UK version) for the three different options: 
     16        - just defining closed seas (and distribute fluxes over global ocean) 
     17        - defining closed seas with a RNF mapping for the American Great Lakes to the St Laurence Seaway 
     18        - defining closed seas with an EMPMR mapping for the American Great Lakes to the St Laurence Seaway 
     19 
     20  3. Choose whether to mask the closea_mask* fields. Not required 
    1521     but makes the fields easier to check. 
    1622 
    … …  
    2733import numpy.ma as ma 
    2834 
    29 def make_closea_masks(config=None,domcfg_file=None,outfile=None,mask=None): 
     35def make_closea_masks(config=None,domcfg_file=None,mask=None): 
    3036 
    3137#========================= 
    … …  
    3642        raise Exception('configuration must be specified') 
    3743 
    38     if outfile is None: 
    39         outfile='closea_masks_'+config+'.nc' 
     44    if domcfg_file is None: 
     45        raise Exception('domain_cfg file must be specified') 
    4046 
    4147    if mask is None: 
    4248        mask=False 
    4349 
    44     with nc.Dataset(domcfg_file,'r') as domcfg_in: 
    45         lon = domcfg_in.variables['nav_lon'][:] 
    46         lat = domcfg_in.variables['nav_lat'][:] 
    47         top_level = domcfg_in.variables['top_level'][0][:] 
     50    domcfg = nc.Dataset(domcfg_file,'r+') 
     51    lon = domcfg.variables['nav_lon'][:] 
     52    lat = domcfg.variables['nav_lat'][:] 
     53    top_level = domcfg.variables['top_level'][0][:] 
    4854 
    4955    nx = top_level.shape[1] 
    … …  
    225231        ncsi2[4]   = 815  ; ncsj2[4]   = 926  
    226232        # runff points the St Laurence Seaway for all Great Lakes 
    227         ncsir1[4:11]   = 873 ; ncsjr1[4:11]   = 909  
    228         ncsir2[4:11]   = 884 ; ncsjr2[4:11]   = 920  
     233        ncsir1[4:10]   = 873 ; ncsjr1[4:10]   = 909  
     234        ncsir2[4:10]   = 884 ; ncsjr2[4:10]   = 920  
    229235 
    230236        # Lake Michigan 
    … …  
    254260 
    255261        # Lake Victoria 
    256         ncsnr[10]   = 1    ; ncstt[10]   = 1    
     262        ncsnr[10]   = 1    ; ncstt[10]   = 0    
     263        ncsi1[10]   = 1274 ; ncsj1[10]   = 672  
     264        ncsi2[10]   = 1289 ; ncsj2[10]   = 687  
     265 
     266    #================================================================ 
     267    elif config == 'eORCA025_UK_empmr': 
     268 
     269        num_closea = 10 
     270        max_runoff_points = 1 
     271        use_runoff_box = True 
     272 
     273        ncsnr = np.zeros(num_closea+1,dtype=np.int)                     ; ncstt = np.zeros(num_closea+1,dtype=np.int) 
     274        ncsi1 = np.zeros(num_closea+1,dtype=np.int)                     ; ncsj1 = np.zeros(num_closea+1,dtype=np.int) 
     275        ncsi2 = np.zeros(num_closea+1,dtype=np.int)                     ; ncsj2 = np.zeros(num_closea+1,dtype=np.int) 
     276        ncsir1 = np.zeros(num_closea+1,dtype=np.int)                    ; ncsjr1 = np.zeros(num_closea+1,dtype=np.int) 
     277        ncsir2 = np.zeros(num_closea+1,dtype=np.int)                    ; ncsjr2 = np.zeros(num_closea+1,dtype=np.int) 
     278        ncsir = np.zeros((num_closea+1,max_runoff_points+1),dtype=np.int) ; ncsjr = np.zeros((num_closea+1,max_runoff_points+1),dtype=np.int) 
     279 
     280        # Caspian Sea 
     281        ncsnr[1]   = 1    ; ncstt[1]   = 0      
     282        ncsi1[1]   = 1330 ; ncsj1[1]   = 831 
     283        ncsi2[1]   = 1375 ; ncsj2[1]   = 981 
     284 
     285        # Aral Sea 
     286        ncsnr[2]   = 1    ; ncstt[2]   = 0      
     287        ncsi1[2]   = 1376 ; ncsj1[2]   = 900 
     288        ncsi2[2]   = 1400 ; ncsj2[2]   = 981 
     289 
     290        # Azov Sea 
     291        ncsnr[3]   = 1    ; ncstt[3]   = 0      
     292        ncsi1[3]   = 1284 ; ncsj1[3]   = 908 
     293        ncsi2[3]   = 1304 ; ncsj2[3]   = 933 
     294 
     295        # Lake Superior 
     296        ncsnr[4]   = 1    ; ncstt[4]   = 2      
     297        ncsi1[4]   = 781  ; ncsj1[4]   = 904  
     298        ncsi2[4]   = 815  ; ncsj2[4]   = 926  
     299        # runff points the St Laurence Seaway for all Great Lakes 
     300        ncsir1[4:10]   = 873 ; ncsjr1[4:10]   = 909  
     301        ncsir2[4:10]   = 884 ; ncsjr2[4:10]   = 920  
     302 
     303        # Lake Michigan 
     304        ncsnr[5]   = 1    ; ncstt[5]   = 2      
     305        ncsi1[5]   = 795  ; ncsj1[5]   = 871              
     306        ncsi2[5]   = 813  ; ncsj2[5]   = 905  
     307 
     308        # Lake Huron part 1 
     309        ncsnr[6]   = 1    ; ncstt[6]   = 2      
     310        ncsi1[6]   = 814  ; ncsj1[6]   = 882              
     311        ncsi2[6]   = 825  ; ncsj2[6]   = 905  
     312 
     313        # Lake Huron part 2 
     314        ncsnr[7]   = 1    ; ncstt[7]   = 2      
     315        ncsi1[7]   = 826  ; ncsj1[7]   = 889              
     316        ncsi2[7]   = 833  ; ncsj2[7]   = 905  
     317 
     318        # Lake Erie 
     319        ncsnr[8]   = 1    ; ncstt[8]   = 2      
     320        ncsi1[8]   = 816  ; ncsj1[8]   = 871              
     321        ncsi2[8]   = 837  ; ncsj2[8]   = 881  
     322 
     323        # Lake Ontario 
     324        ncsnr[9]   = 1    ; ncstt[9]   = 2      
     325        ncsi1[9]   = 831  ; ncsj1[9]   = 882              
     326        ncsi2[9]   = 847  ; ncsj2[9]   = 889  
     327 
     328        # Lake Victoria 
     329        ncsnr[10]   = 1    ; ncstt[10]   = 0    
    257330        ncsi1[10]   = 1274 ; ncsj1[10]   = 672  
    258331        ncsi2[10]   = 1289 ; ncsj2[10]   = 687  
    … …  
    331404 
    332405#===================================== 
    333 # 4. Generate mask netcdf file 
    334 #===================================== 
    335  
    336     outfile=nc.Dataset(outfile,'w') 
    337     outfile.createDimension('x',size=nx) 
    338     outfile.createDimension('y',size=ny) 
    339     outfile.createVariable('nav_lon',datatype='f',dimensions=('y','x')) 
    340     outfile.createVariable('nav_lat',datatype='f',dimensions=('y','x')) 
    341     outfile.variables['nav_lon'][:]=lon[:] 
    342     outfile.variables['nav_lat'][:]=lat[:] 
    343     outfile.createVariable('closea_mask',datatype='i',dimensions=('y','x'),fill_value=closea_mask.fill_value) 
    344     outfile.variables['closea_mask'][:]=closea_mask 
     406# 4. Append masks to domain_cfg file. 
     407#===================================== 
     408 
     409    domcfg.createVariable('closea_mask',datatype='i',dimensions=('y','x'),fill_value=closea_mask.fill_value,chunksizes=(1000,1000)) 
     410    domcfg.variables['closea_mask'][:]=closea_mask 
    345411    if rnf_count > 0: 
    346         outfile.createVariable('closea_mask_rnf',datatype='i',dimensions=('y','x'),fill_value=closea_mask_rnf.fill_value) 
    347         outfile.variables['closea_mask_rnf'][:]=closea_mask_rnf 
     412        domcfg.createVariable('closea_mask_rnf',datatype='i',dimensions=('y','x'),fill_value=closea_mask_rnf.fill_value,chunksizes=(1000,1000)) 
     413        domcfg.variables['closea_mask_rnf'][:]=closea_mask_rnf 
    348414    if empmr_count > 0: 
    349         outfile.createVariable('closea_mask_empmr',datatype='i',dimensions=('y','x'),fill_value=closea_mask_empmr.fill_value) 
    350         outfile.variables['closea_mask_empmr'][:]=closea_mask_empmr 
    351  
    352     outfile.close() 
     415        domcfg.createVariable('closea_mask_empmr',datatype='i',dimensions=('y','x'),fill_value=closea_mask_empmr.fill_value,chunksizes=(1000,1000)) 
     416        domcfg.variables['closea_mask_empmr'][:]=closea_mask_empmr 
     417 
     418    domcfg.close() 
    353419 
    354420 
    … …  
    360426    parser.add_argument("-d", "--domcfg", action="store",dest="domcfg_file",default=None, 
    361427                    help="domcfg file (input)") 
    362     parser.add_argument("-o", "--out", action="store",dest="outfile",default=None, 
    363                     help="output closea masks file") 
    364428    parser.add_argument("-m", "--mask", action="store_true",dest="mask",default=False, 
    365429                    help="mask output file based on top_level in domcfg file") 
    … …  
    367431    args = parser.parse_args() 
    368432 
    369     make_closea_masks(config=args.config,domcfg_file=args.domcfg_file,outfile=args.outfile,mask=args.mask) 
     433    make_closea_masks(config=args.config,domcfg_file=args.domcfg_file,mask=args.mask) 
Note: See TracChangeset for help on using the changeset viewer.