Changeset 11207 for NEMO/branches/2019/ENHANCE-03_closea/src/OCE/DOM
- Timestamp:
- 2019-07-02T19:32:41+02:00 (5 years ago)
- Location:
- NEMO/branches/2019/ENHANCE-03_closea
- Files:
-
- 3 edited
- 1 copied
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2019/ENHANCE-03_closea/src/OCE/DOM/closea.F90
r10425 r11207 11 11 !! 4.0 ! 2016-06 (G. Madec) move to usrdef_closea, remove clo_ups 12 12 !! 4.0 ! 2017-12 (D. Storkey) new formulation based on masks read from file 13 !! 4.1 ! 2019-07 (P. Mathiot) update to the new domcfg.nc input file 13 14 !!---------------------------------------------------------------------- 14 15 … … 35 36 36 37 PUBLIC dom_clo ! called by domain module 38 PUBLIC sbc_clo_init ! called by sbcmod module 37 39 PUBLIC sbc_clo ! called by sbcmod module 38 40 PUBLIC clo_rnf ! called by sbcrnf module 39 41 PUBLIC clo_bat ! called in domzgr module 40 42 41 LOGICAL, PUBLIC :: ln_closea !: T => keep closed seas (defined by closea_mask field) in the domain and apply 42 !: special treatment of freshwater fluxes. 43 !: F => suppress closed seas (defined by closea_mask field) from the bathymetry 44 !: at runtime. 45 !: If there is no closea_mask field in the domain_cfg file or we do not use 46 !: a domain_cfg file then this logical does nothing. 47 !: 48 LOGICAL, PUBLIC :: l_sbc_clo !: T => Closed seas defined, apply special treatment of freshwater fluxes. 49 !: F => No closed seas defined (closea_mask field not found). 50 LOGICAL, PUBLIC :: l_clo_rnf !: T => Some closed seas output freshwater (RNF or EMPMR) to specified runoff points. 51 INTEGER, PUBLIC :: jncs !: number of closed seas (inferred from closea_mask field) 52 INTEGER, PUBLIC :: jncsr !: number of closed seas rnf mappings (inferred from closea_mask_rnf field) 53 INTEGER, PUBLIC :: jncse !: number of closed seas empmr mappings (inferred from closea_mask_empmr field) 54 55 INTEGER, PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: closea_mask !: mask of integers defining closed seas 56 INTEGER, PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: closea_mask_rnf !: mask of integers defining closed seas rnf mappings 57 INTEGER, PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: closea_mask_empmr !: mask of integers defining closed seas empmr mappings 58 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:) :: surf !: closed sea surface areas 59 !: (and residual global surface area) 60 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:) :: surfr !: closed sea target rnf surface areas 61 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:) :: surfe !: closed sea target empmr surface areas 43 LOGICAL, PUBLIC :: ln_maskcs !: mask all closed sea 44 LOGICAL, PUBLIC :: ln_mask_csundef !: mask all closed sea 45 LOGICAL, PUBLIC :: ln_clo_rnf !: closed sea treated as runoff (update rnf mask) 46 47 LOGICAL, PUBLIC :: l_sbc_clo !: T => net evap/precip over closed seas spread outover the globe/river mouth 48 LOGICAL, PUBLIC :: l_clo_rnf !: T => Some closed seas output freshwater (RNF) to specified runoff points. 49 50 INTEGER, PUBLIC :: jncsg !: number of closed seas global mappings (inferred from closea_mask_glo field) 51 INTEGER, PUBLIC :: jncsr !: number of closed seas rnf mappings (inferred from closea_mask_rnf field) 52 INTEGER, PUBLIC :: jncse !: number of closed seas empmr mappings (inferred from closea_mask_emp field) 53 54 INTEGER, PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: jcsgrpg, jcsgrpr, jcsgrpe 55 56 INTEGER, PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:) :: mask_opnsea, mask_csundef 57 58 INTEGER, PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:) :: mask_csglo, mask_csgrpglo !: mask of integers defining closed seas 59 INTEGER, PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:) :: mask_csrnf, mask_csgrprnf !: mask of integers defining closed seas rnf mappings 60 INTEGER, PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:) :: mask_csemp, mask_csgrpemp !: mask of integers defining closed seas empmr mappings 61 62 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: rsurfsrcg, rsurftrgg !: closed sea target glo surface areas 63 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: rsurfsrcr, rsurftrgr !: closed sea target rnf surface areas 64 REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:) :: rsurfsrce, rsurftrge !: closed sea target emp surface areas 62 65 63 66 !! * Substitutions … … 79 82 !! just the thermodynamic processes are applied. 80 83 !! 81 !! ** Action : Read closea_mask* fields (if they exist) from domain_cfg file and infer 82 !! number of closed seas from closea_mask field. 83 !! closea_mask : integer values defining closed seas (or groups of closed seas) 84 !! closea_mask_rnf : integer values defining mappings from closed seas or groups of 85 !! closed seas to a runoff area for downwards flux only. 86 !! closea_mask_empmr : integer values defining mappings from closed seas or groups of 87 !! closed seas to a runoff area for net fluxes. 88 !! 89 !! Python code to generate the closea_masks* fields from the old-style indices 90 !! definitions is available at TOOLS/DOMAINcfg/make_closea_masks.py 91 !!---------------------------------------------------------------------- 92 INTEGER :: inum ! input file identifier 93 INTEGER :: ierr ! error code 94 INTEGER :: id ! netcdf variable ID 95 96 REAL(wp), DIMENSION(jpi,jpj) :: zdata_in ! temporary real array for input 97 !!---------------------------------------------------------------------- 98 ! 84 !! ** Action : Read mask_cs* fields (if needed) from domain_cfg file and infer 85 !! number of closed seas for each case (glo, rnf, emp) from mask_cs* field. 86 !! mask_csglo and mask_csgrpglo : integer values defining mappings from closed seas and associated groups to the open ocean for net fluxes. 87 !! mask_csrnf and mask_csgrprnf : integer values defining mappings from closed seas and associated groups to a runoff area for downwards flux only. 88 !! mask_csemp and mask_csgrpemp : integer values defining mappings from closed seas and associated groups to a runoff area for net fluxes. 89 !!---------------------------------------------------------------------- 90 INTEGER :: ios ! io status 91 !! 92 NAMELIST/namclo/ ln_maskcs, ln_mask_csundef, ln_clo_rnf 93 !!--------------------------------------------------------------------- 94 !! 95 REWIND( numnam_ref ) ! Namelist namclo in reference namelist : Lateral momentum boundary condition 96 READ ( numnam_ref, namclo, IOSTAT = ios, ERR = 901 ) 97 901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namclo in reference namelist', lwp ) 98 REWIND( numnam_cfg ) ! Namelist namclo in configuration namelist : Lateral momentum boundary condition 99 READ ( numnam_cfg, namclo, IOSTAT = ios, ERR = 902 ) 100 902 IF( ios > 0 ) CALL ctl_nam ( ios , 'namclo in configuration namelist', lwp ) 101 IF(lwm) WRITE ( numond, namclo ) 102 !! 99 103 IF(lwp) WRITE(numout,*) 100 104 IF(lwp) WRITE(numout,*)'dom_clo : read in masks to define closed seas ' 101 105 IF(lwp) WRITE(numout,*)'~~~~~~~' 106 !! 107 !! check option compatibility 108 IF( .NOT. ln_read_cfg ) THEN 109 CALL ctl_stop('Suppression of closed seas does not work with ln_read_cfg = .true. . Set ln_closea = .false. .') 110 ENDIF 111 !! 112 IF( (.NOT. ln_maskcs) .AND. ln_diurnal_only ) THEN 113 CALL ctl_stop('Special handling of freshwater fluxes over closed seas not compatible with ln_diurnal_only.') 114 END IF 115 ! 116 l_sbc_clo = .false. ; l_clo_rnf = .false. 117 IF (.NOT. ln_maskcs) l_sbc_clo = .true. 118 IF (.NOT. ln_maskcs .AND. ln_clo_rnf) l_clo_rnf = .true. 102 119 ! 103 120 ! read the closed seas masks (if they exist) from domain_cfg file (if it exists) 104 121 ! ------------------------------------------------------------------------------ 105 122 ! 106 IF( ln_read_cfg) THEN 107 ! 108 CALL iom_open( cn_domcfg, inum ) 109 ! 110 id = iom_varid(inum, 'closea_mask', ldstop = .false.) 111 IF( id > 0 ) THEN 112 l_sbc_clo = .true. 113 ALLOCATE( closea_mask(jpi,jpj) , STAT=ierr ) 114 IF( ierr /= 0 ) CALL ctl_stop( 'STOP', 'dom_clo: failed to allocate closea_mask array') 115 zdata_in(:,:) = 0.0 116 CALL iom_get ( inum, jpdom_data, 'closea_mask', zdata_in ) 117 closea_mask(:,:) = NINT(zdata_in(:,:)) * tmask(:,:,1) 118 ! number of closed seas = global maximum value in closea_mask field 119 jncs = maxval(closea_mask(:,:)) 120 CALL mpp_max('closea', jncs) 121 IF( jncs > 0 ) THEN 122 IF( lwp ) WRITE(numout,*) 'Number of closed seas : ',jncs 123 ELSE 124 CALL ctl_stop( 'Problem with closea_mask field in domain_cfg file. Has no values > 0 so no closed seas defined.') 125 ENDIF 126 ELSE 127 IF( lwp ) WRITE(numout,*) 128 IF( lwp ) WRITE(numout,*) ' ==>>> closea_mask field not found in domain_cfg file.' 129 IF( lwp ) WRITE(numout,*) ' No closed seas defined.' 130 IF( lwp ) WRITE(numout,*) 131 l_sbc_clo = .false. 132 jncs = 0 133 ENDIF 134 135 l_clo_rnf = .false. 136 137 IF( l_sbc_clo ) THEN ! No point reading in closea_mask_rnf or closea_mask_empmr fields if no closed seas defined. 138 139 id = iom_varid(inum, 'closea_mask_rnf', ldstop = .false.) 140 IF( id > 0 ) THEN 141 l_clo_rnf = .true. 142 ALLOCATE( closea_mask_rnf(jpi,jpj) , STAT=ierr ) 143 IF( ierr /= 0 ) CALL ctl_stop( 'STOP', 'dom_clo: failed to allocate closea_mask_rnf array') 144 CALL iom_get ( inum, jpdom_data, 'closea_mask_rnf', zdata_in ) 145 closea_mask_rnf(:,:) = NINT(zdata_in(:,:)) * tmask(:,:,1) 146 ! number of closed seas rnf mappings = global maximum in closea_mask_rnf field 147 jncsr = maxval(closea_mask_rnf(:,:)) 148 CALL mpp_max('closea', jncsr) 149 IF( jncsr > 0 ) THEN 150 IF( lwp ) WRITE(numout,*) 'Number of closed seas rnf mappings : ',jncsr 151 ELSE 152 CALL ctl_stop( 'Problem with closea_mask_rnf field in domain_cfg file. Has no values > 0 so no closed seas rnf mappings defined.') 153 ENDIF 154 ELSE 155 IF( lwp ) WRITE(numout,*) 'closea_mask_rnf field not found in domain_cfg file. No closed seas rnf mappings defined.' 156 jncsr = 0 157 ENDIF 158 159 id = iom_varid(inum, 'closea_mask_empmr', ldstop = .false.) 160 IF( id > 0 ) THEN 161 l_clo_rnf = .true. 162 ALLOCATE( closea_mask_empmr(jpi,jpj) , STAT=ierr ) 163 IF( ierr /= 0 ) CALL ctl_stop( 'STOP', 'dom_clo: failed to allocate closea_mask_empmr array') 164 CALL iom_get ( inum, jpdom_data, 'closea_mask_empmr', zdata_in ) 165 closea_mask_empmr(:,:) = NINT(zdata_in(:,:)) * tmask(:,:,1) 166 ! number of closed seas empmr mappings = global maximum value in closea_mask_empmr field 167 jncse = maxval(closea_mask_empmr(:,:)) 168 CALL mpp_max('closea', jncse) 169 IF( jncse > 0 ) THEN 170 IF( lwp ) WRITE(numout,*) 'Number of closed seas empmr mappings : ',jncse 171 ELSE 172 CALL ctl_stop( 'Problem with closea_mask_empmr field in domain_cfg file. Has no values > 0 so no closed seas empmr mappings defined.') 173 ENDIF 174 ELSE 175 IF( lwp ) WRITE(numout,*) 'closea_mask_empmr field not found in domain_cfg file. No closed seas empmr mappings defined.' 176 jncse = 0 177 ENDIF 178 179 ENDIF ! l_sbc_clo 180 ! 181 CALL iom_close( inum ) 182 ! 183 ELSE ! ln_read_cfg = .false. so no domain_cfg file 184 IF( lwp ) WRITE(numout,*) 'No domain_cfg file so no closed seas defined.' 185 l_sbc_clo = .false. 186 l_clo_rnf = .false. 187 ENDIF 188 ! 123 ! load mask of open sea and undefined closed seas 124 CALL alloc_csmask( mask_opnsea ) 125 CALL read_csmask( cn_domcfg, 'mask_opensea' , mask_opnsea ) 126 ! 127 IF ( ln_maskcs ) THEN 128 ! not special treatment of closed sea 129 l_sbc_clo = .false. ; l_clo_rnf = .false. 130 ELSE 131 ! special treatment of closed seas 132 l_sbc_clo = .true. 133 ! 134 ! river mouth from lakes added to rnf mask for special treatment 135 IF ( ln_clo_rnf) l_clo_rnf = .true. 136 ! 137 IF ( ln_mask_csundef) THEN 138 ! load undef cs mask (1 in undef closed sea) 139 CALL alloc_csmask( mask_csundef ) 140 CALL read_csmask( cn_domcfg, 'mask_csundef', mask_csundef ) 141 ! revert the mask for masking in domzgr 142 mask_csundef = 1 - mask_csundef 143 END IF 144 ! 145 ! allocate source mask 146 CALL alloc_csmask( mask_csglo ) 147 CALL alloc_csmask( mask_csrnf ) 148 CALL alloc_csmask( mask_csemp ) 149 ! 150 ! load source mask of cs for each cases 151 CALL read_csmask( cn_domcfg, 'mask_csglo', mask_csglo ) 152 CALL read_csmask( cn_domcfg, 'mask_csrnf', mask_csrnf ) 153 CALL read_csmask( cn_domcfg, 'mask_csemp', mask_csemp ) 154 ! 155 ! compute number of cs for each cases 156 jncsg = MAXVAL( mask_csglo(:,:) ) ; CALL mpp_max( 'closea', jncsg ) 157 jncsr = MAXVAL( mask_csrnf(:,:) ) ; CALL mpp_max( 'closea', jncsr ) 158 jncse = MAXVAL( mask_csemp(:,:) ) ; CALL mpp_max( 'closea', jncse ) 159 ! 160 ! allocate closed sea group masks 161 CALL alloc_csmask( mask_csgrpglo ) 162 CALL alloc_csmask( mask_csgrprnf ) 163 CALL alloc_csmask( mask_csgrpemp ) 164 165 ! load mask of cs group for each cases 166 CALL read_csmask( cn_domcfg, 'mask_csgrpglo', mask_csgrpglo ) 167 CALL read_csmask( cn_domcfg, 'mask_csgrprnf', mask_csgrprnf ) 168 CALL read_csmask( cn_domcfg, 'mask_csgrpemp', mask_csgrpemp ) 169 ! 170 ! Allocate cs variables (surf) 171 CALL alloc_cssurf( jncsg, rsurfsrcg, rsurftrgg ) 172 CALL alloc_cssurf( jncsr, rsurfsrcr, rsurftrgr ) 173 CALL alloc_cssurf( jncse, rsurfsrce, rsurftrge ) 174 ! 175 ! Allocate cs group variables (jcsgrp) 176 CALL alloc_csgrp( jncsg, jcsgrpg ) 177 CALL alloc_csgrp( jncsr, jcsgrpr ) 178 CALL alloc_csgrp( jncse, jcsgrpe ) 179 ! 180 END IF 189 181 END SUBROUTINE dom_clo 190 182 191 192 SUBROUTINE sbc_clo( kt ) 183 SUBROUTINE sbc_clo_init 184 185 ! compute source surface area 186 CALL get_cssrcsurf( jncsg, mask_csglo, rsurfsrcg ) 187 CALL get_cssrcsurf( jncsr, mask_csrnf, rsurfsrcr ) 188 CALL get_cssrcsurf( jncse, mask_csemp, rsurfsrce ) 189 ! 190 ! compute target surface area and group number (jcsgrp) for all cs and cases 191 ! glo could be simpler but for lisibility, all treated the same way 192 ! It is only done once, so not a big deal 193 CALL get_cstrgsurf( jncsg, mask_csglo, mask_csgrpglo, rsurftrgg, jcsgrpg ) 194 CALL get_cstrgsurf( jncsr, mask_csrnf, mask_csgrprnf, rsurftrgr, jcsgrpr ) 195 CALL get_cstrgsurf( jncse, mask_csemp, mask_csgrpemp, rsurftrge, jcsgrpe ) 196 ! 197 ! print out in ocean.ouput 198 CALL prt_csctl( jncsg, rsurfsrcg, rsurftrgg, jcsgrpg, 'glo' ) 199 CALL prt_csctl( jncsr, rsurfsrcr, rsurftrgr, jcsgrpr, 'rnf' ) 200 CALL prt_csctl( jncse, rsurfsrce, rsurftrge, jcsgrpe, 'emp' ) 201 202 END SUBROUTINE sbc_clo_init 203 204 SUBROUTINE get_cssrcsurf(kncs, pmaskcs, psurfsrc) 205 206 ! subroutine parameters 207 INTEGER, INTENT(in ) :: kncs 208 INTEGER, DIMENSION(:,:), INTENT(in ) :: pmaskcs 209 210 REAL(wp), DIMENSION(:) , INTENT(inout) :: psurfsrc 211 212 ! local variables 213 INTEGER, DIMENSION(jpi,jpj) :: zmsksrc 214 INTEGER :: jcs 215 216 DO jcs = 1,kncs 217 ! 218 ! build river mouth mask for this lake 219 WHERE ( pmaskcs == jcs ) 220 zmsksrc = 1 221 ELSE WHERE 222 zmsksrc = 0 223 END WHERE 224 ! 225 ! compute target area 226 psurfsrc(jcs) = glob_sum('closea', e1e2t(:,:) * zmsksrc(:,:) ) 227 ! 228 END DO 229 230 END SUBROUTINE 231 232 SUBROUTINE get_cstrgsurf(kncs, pmaskcs, pmaskcsgrp, psurftrg, kcsgrp ) 233 234 ! subroutine parameters 235 INTEGER, INTENT(in ) :: kncs 236 INTEGER, DIMENSION(:), INTENT(inout) :: kcsgrp 237 INTEGER, DIMENSION(:,:), INTENT(in ) :: pmaskcs, pmaskcsgrp 238 239 REAL(wp), DIMENSION(:) , INTENT(inout) :: psurftrg 240 241 ! local variables 242 INTEGER, DIMENSION(jpi,jpj) :: zmskgrp, zmsksrc, zmsktrg 243 INTEGER :: jcs, jtmp 244 245 DO jcs = 1,kncs 246 ! 247 ! find group number 248 zmskgrp = pmaskcsgrp 249 zmsksrc = pmaskcs 250 ! 251 ! set value where cs is 252 zmsktrg = HUGE(1) 253 WHERE ( zmsksrc == jcs ) zmsktrg = jcs 254 ! 255 ! zmsk = HUGE outside the cs number jcs 256 ! ktmp = jcs - group number 257 ! jgrp = group corresponding to the cs jcs 258 zmsktrg = zmsktrg - zmskgrp 259 jtmp = MINVAL(zmsktrg) ; CALL mpp_min('closea',jtmp) 260 kcsgrp(jcs) = jcs - jtmp 261 ! 262 ! build river mouth mask for this lake 263 WHERE ( zmskgrp * mask_opnsea == kcsgrp(jcs) ) 264 zmsktrg = 1 265 ELSE WHERE 266 zmsktrg = 0 267 END WHERE 268 ! 269 ! compute target area 270 psurftrg(jcs) = glob_sum('closea', e1e2t(:,:) * zmsktrg(:,:) ) 271 ! 272 END DO 273 274 END SUBROUTINE 275 276 SUBROUTINE prt_csctl(kncs, psurfsrc, psurftrg, kcsgrp, pcstype) 277 ! subroutine parameters 278 INTEGER, INTENT(in ) :: kncs 279 INTEGER, DIMENSION(:), INTENT(in ) :: kcsgrp 280 ! 281 REAL(wp), DIMENSION(:), INTENT(in ) :: psurfsrc, psurftrg 282 ! 283 CHARACTER(256), INTENT(in ) :: pcstype 284 ! 285 ! local variable 286 INTEGER :: jcs 287 288 IF ( lwp .AND. kncs > 0 ) THEN 289 WRITE(numout,*)'' 290 ! 291 WRITE(numout,*)'Closed sea target ',TRIM(pcstype),' : ' 292 ! 293 DO jcs = 1,kncs 294 WRITE(numout,FMT='(3a,i3,a,i3)') ' ',TRIM(pcstype),' closed sea id is ',jcs,' and trg id is : ', kcsgrp(jcs) 295 WRITE(numout,FMT='(a,f12.2)' ) ' src surface areas (km2) : ', psurfsrc(jcs) * 1.0e-6 296 WRITE(numout,FMT='(a,f12.2)' ) ' trg surface areas (km2) : ', psurftrg(jcs) * 1.0e-6 297 END DO 298 ! 299 WRITE(numout,*)'' 300 END IF 301 END SUBROUTINE 302 303 SUBROUTINE sbc_clo( kt ) ! to be move in SBC in a file sbcclo ??? 193 304 !!--------------------------------------------------------------------- 194 305 !! *** ROUTINE sbc_clo *** … … 203 314 !!---------------------------------------------------------------------- 204 315 INTEGER , INTENT(in ) :: kt ! ocean model time step 205 !206 INTEGER :: ierr207 INTEGER :: jc, jcr, jce ! dummy loop indices208 REAL(wp), PARAMETER :: rsmall = 1.e-20_wp ! Closed sea correction epsilon209 REAL(wp) :: zfwf_total, zcoef, zcoef1 !210 REAL(wp), DIMENSION(jncs) :: zfwf !:211 REAL(wp), DIMENSION(jncsr+1) :: zfwfr !: freshwater fluxes over closed seas212 REAL(wp), DIMENSION(jncse+1) :: zfwfe !:213 REAL(wp), DIMENSION(jpi,jpj) :: ztmp2d ! 2D workspace214 316 !!---------------------------------------------------------------------- 215 317 ! 216 318 IF( ln_timing ) CALL timing_start('sbc_clo') 217 319 ! 218 ! !------------------! 219 IF( kt == nit000 ) THEN ! Initialisation ! 220 ! !------------------! 221 IF(lwp) WRITE(numout,*) 222 IF(lwp) WRITE(numout,*)'sbc_clo : closed seas ' 223 IF(lwp) WRITE(numout,*)'~~~~~~~' 224 225 ALLOCATE( surf(jncs+1) , STAT=ierr ) 226 IF( ierr /= 0 ) CALL ctl_stop( 'STOP', 'sbc_clo: failed to allocate surf array') 227 surf(:) = 0.e0_wp 228 ! 229 ! jncsr can be zero so add 1 to avoid allocating zero-length array 230 ALLOCATE( surfr(jncsr+1) , STAT=ierr ) 231 IF( ierr /= 0 ) CALL ctl_stop( 'STOP', 'sbc_clo: failed to allocate surfr array') 232 surfr(:) = 0.e0_wp 233 ! 234 ! jncse can be zero so add 1 to avoid allocating zero-length array 235 ALLOCATE( surfe(jncse+1) , STAT=ierr ) 236 IF( ierr /= 0 ) CALL ctl_stop( 'STOP', 'sbc_clo: failed to allocate surfe array') 237 surfe(:) = 0.e0_wp 238 ! 239 surf(jncs+1) = glob_sum( 'closea', e1e2t(:,:) ) ! surface of the global ocean 240 ! 241 ! ! surface areas of closed seas 242 DO jc = 1, jncs 243 ztmp2d(:,:) = 0.e0_wp 244 WHERE( closea_mask(:,:) == jc ) ztmp2d(:,:) = e1e2t(:,:) * tmask_i(:,:) 245 surf(jc) = glob_sum( 'closea', ztmp2d(:,:) ) 246 END DO 247 ! 248 ! jncs+1 : surface area of global ocean, closed seas excluded 249 surf(jncs+1) = surf(jncs+1) - SUM(surf(1:jncs)) 250 ! 251 ! ! surface areas of rnf target areas 252 IF( jncsr > 0 ) THEN 253 DO jcr = 1, jncsr 254 ztmp2d(:,:) = 0.e0_wp 255 WHERE( closea_mask_rnf(:,:) == jcr .and. closea_mask(:,:) == 0 ) ztmp2d(:,:) = e1e2t(:,:) * tmask_i(:,:) 256 surfr(jcr) = glob_sum( 'closea', ztmp2d(:,:) ) 257 END DO 258 ENDIF 259 ! 260 ! ! surface areas of empmr target areas 261 IF( jncse > 0 ) THEN 262 DO jce = 1, jncse 263 ztmp2d(:,:) = 0.e0_wp 264 WHERE( closea_mask_empmr(:,:) == jce .and. closea_mask(:,:) == 0 ) ztmp2d(:,:) = e1e2t(:,:) * tmask_i(:,:) 265 surfe(jce) = glob_sum( 'closea', ztmp2d(:,:) ) 266 END DO 267 ENDIF 268 ! 269 IF(lwp) WRITE(numout,*)' Closed sea surface areas (km2)' 270 DO jc = 1, jncs 271 IF(lwp) WRITE(numout,FMT='(1I3,5X,ES12.2)') jc, surf(jc) * 1.0e-6 272 END DO 273 IF(lwp) WRITE(numout,FMT='(A,ES12.2)') 'Global surface area excluding closed seas (km2): ', surf(jncs+1) * 1.0e-6 274 ! 275 IF(jncsr > 0) THEN 276 IF(lwp) WRITE(numout,*)' Closed sea target rnf surface areas (km2)' 277 DO jcr = 1, jncsr 278 IF(lwp) WRITE(numout,FMT='(1I3,5X,ES12.2)') jcr, surfr(jcr) * 1.0e-6 279 END DO 280 ENDIF 281 ! 282 IF(jncse > 0) THEN 283 IF(lwp) WRITE(numout,*)' Closed sea target empmr surface areas (km2)' 284 DO jce = 1, jncse 285 IF(lwp) WRITE(numout,FMT='(1I3,5X,ES12.2)') jce, surfe(jce) * 1.0e-6 286 END DO 287 ENDIF 288 ENDIF 289 ! 290 ! !--------------------! 291 ! ! update emp ! 292 ! !--------------------! 293 294 zfwf_total = 0._wp 295 296 ! 297 ! 1. Work out total freshwater fluxes over closed seas from EMP - RNF. 298 ! 299 zfwf(:) = 0.e0_wp 300 DO jc = 1, jncs 301 ztmp2d(:,:) = 0.e0_wp 302 WHERE( closea_mask(:,:) == jc ) ztmp2d(:,:) = e1e2t(:,:) * ( emp(:,:)-rnf(:,:) ) * tmask_i(:,:) 303 zfwf(jc) = glob_sum( 'closea', ztmp2d(:,:) ) 304 END DO 305 zfwf_total = SUM(zfwf) 306 307 zfwfr(:) = 0.e0_wp 308 IF( jncsr > 0 ) THEN 309 ! 310 ! 2. Work out total FW fluxes over rnf source areas and add to rnf target areas. 311 ! Where zfwf is negative add flux at specified runoff points and subtract from fluxes for global redistribution. 312 ! Where positive leave in global redistribution total. 313 ! 314 DO jcr = 1, jncsr 315 ! 316 ztmp2d(:,:) = 0.e0_wp 317 WHERE( closea_mask_rnf(:,:) == jcr .and. closea_mask(:,:) > 0 ) ztmp2d(:,:) = e1e2t(:,:) * ( emp(:,:)-rnf(:,:) ) * tmask_i(:,:) 318 zfwfr(jcr) = glob_sum( 'closea', ztmp2d(:,:) ) 319 ! 320 ! The following if avoids the redistribution of the round off 321 IF ( ABS(zfwfr(jcr) / surf(jncs+1) ) > rsmall) THEN 322 ! 323 ! Add residuals to target runoff points if negative and subtract from total to be added globally 324 IF( zfwfr(jcr) < 0.0 ) THEN 325 zfwf_total = zfwf_total - zfwfr(jcr) 326 zcoef = zfwfr(jcr) / surfr(jcr) 327 zcoef1 = rcp * zcoef 328 WHERE( closea_mask_rnf(:,:) == jcr .and. closea_mask(:,:) == 0.0) 329 emp(:,:) = emp(:,:) + zcoef 330 qns(:,:) = qns(:,:) - zcoef1 * sst_m(:,:) 331 ENDWHERE 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. 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( 'closea', 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) 355 zcoef1 = rcp * zcoef 356 WHERE( closea_mask_empmr(:,:) == jce .and. closea_mask(:,:) == 0.0) 357 emp(:,:) = emp(:,:) + zcoef 358 qns(:,:) = qns(:,:) - zcoef1 * sst_m(:,:) 359 ENDWHERE 360 ! 361 ENDIF 362 END DO 363 ENDIF ! jncse > 0 364 365 ! 366 ! 4. Spread residual flux over global ocean. 367 ! 368 ! The following if avoids the redistribution of the round off 369 IF ( ABS(zfwf_total / surf(jncs+1) ) > rsmall) THEN 370 zcoef = zfwf_total / surf(jncs+1) 320 ! update emp and qns 321 CALL sbc_csupdate( jncsg, jcsgrpg, mask_csglo, mask_csgrpglo, rsurfsrcg, rsurftrgg, 'glo', mask_opnsea, rsurftrgg ) 322 CALL sbc_csupdate( jncsr, jcsgrpr, mask_csrnf, mask_csgrprnf, rsurfsrcr, rsurftrgr, 'rnf', mask_opnsea, rsurftrgg ) 323 CALL sbc_csupdate( jncse, jcsgrpe, mask_csemp, mask_csgrpemp, rsurfsrce, rsurftrge, 'emp', mask_opnsea, rsurftrgg ) 324 ! 325 ! is this really useful ?????? 326 emp(:,:) = emp(:,:) * tmask(:,:,1) 327 qns(:,:) = qns(:,:) * tmask(:,:,1) 328 ! 329 ! is this really useful ?????? 330 CALL lbc_lnk( 'closea', emp , 'T', 1._wp ) 331 CALL lbc_lnk( 'closea', qns , 'T', 1._wp ) 332 ! 333 END SUBROUTINE sbc_clo 334 335 SUBROUTINE sbc_csupdate(kncs, kcsgrp, pmsk_src, pmsk_trg, psurfsrc, psurftrg, pcstype, pmsk_opnsea, psurf_opnsea) 336 337 ! subroutine parameters 338 INTEGER, INTENT(in ) :: kncs 339 INTEGER, DIMENSION(: ), INTENT(in ) :: kcsgrp 340 INTEGER, DIMENSION(:,:), INTENT(in ) :: pmsk_src, pmsk_trg, pmsk_opnsea 341 342 REAL(wp), DIMENSION(:), INTENT(inout) :: psurfsrc, psurftrg, psurf_opnsea 343 344 CHARACTER(256), INTENT(in ) :: pcstype 345 346 ! local variables 347 INTEGER :: jcs 348 INTEGER, DIMENSION(jpi,jpj) :: zmsk_src, zmsk_trg 349 350 REAL(wp) :: zcoef, zcoef1, ztmp 351 REAL(wp) :: zcsfwf 352 REAL(wp) :: zsurftrg 353 354 DO jcs = 1, kncs 355 !! 356 !! 0. get mask of each closed sea 357 zmsk_src(:,:) = 0 358 WHERE ( pmsk_src(:,:) == jcs ) zmsk_src = 1 359 !! 360 !! 1. Work out net freshwater fluxes over each closed seas from EMP - RNF. 361 zcsfwf = glob_sum( 'closea', e1e2t(:,:) * ( emp(:,:)-rnf(:,:) ) * zmsk_src ) 362 !! 363 !! 2. Deal with runoff special case (net evaporation spread globally) 364 IF (pcstype == 'rnf' .AND. zcsfwf > 0) THEN 365 zsurftrg = psurf_opnsea(1) 366 zmsk_trg = pmsk_opnsea 367 ELSE 368 zsurftrg = psurftrg(jcs) 369 zmsk_trg = pmsk_trg 370 END IF 371 zmsk_trg = zmsk_trg * pmsk_opnsea 372 !! 373 !! 3. Add residuals to target points 374 zcoef = zcsfwf / zsurftrg 371 375 zcoef1 = rcp * zcoef 372 WHERE( closea_mask(:,:) == 0)376 WHERE( zmsk_trg(:,:) == kcsgrp(jcs) ) 373 377 emp(:,:) = emp(:,:) + zcoef 374 378 qns(:,:) = qns(:,:) - zcoef1 * sst_m(:,:) 375 379 ENDWHERE 376 ENDIF 377 378 ! 379 ! 5. Subtract area means from emp (and qns) over closed seas to give zero mean FW flux over each sea. 380 ! 381 DO jc = 1, jncs 382 ! The following if avoids the redistribution of the round off 383 IF ( ABS(zfwf(jc) / surf(jncs+1) ) > rsmall) THEN 384 ! 385 ! Subtract residuals from fluxes over closed sea 386 zcoef = zfwf(jc) / surf(jc) 387 zcoef1 = rcp * zcoef 388 WHERE( closea_mask(:,:) == jc ) 389 emp(:,:) = emp(:,:) - zcoef 390 qns(:,:) = qns(:,:) + zcoef1 * sst_m(:,:) 391 ENDWHERE 392 ! 393 ENDIF 394 END DO 395 ! 396 emp (:,:) = emp (:,:) * tmask(:,:,1) 397 ! 398 CALL lbc_lnk( 'closea', emp , 'T', 1._wp ) 399 ! 400 END SUBROUTINE sbc_clo 380 !! 381 !! 4. Subtract residuals from source points 382 zcoef = zcsfwf / psurfsrc(jcs) 383 zcoef1 = rcp * zcoef 384 WHERE( pmsk_src(:,:) == jcs ) 385 emp(:,:) = emp(:,:) - zcoef 386 qns(:,:) = qns(:,:) + zcoef1 * sst_m(:,:) 387 ENDWHERE 388 !! 389 END DO ! jcs 390 391 END SUBROUTINE 392 401 393 402 394 SUBROUTINE clo_rnf( p_rnfmsk ) 403 395 !!--------------------------------------------------------------------- 404 !! *** ROUTINE sbc_rnf ***396 !! *** ROUTINE clo_rnf *** 405 397 !! 406 398 !! ** Purpose : allow the treatment of closed sea outflow grid-points … … 412 404 !! ** Action : update (p_)mskrnf (set 1 at closed sea outflow) 413 405 !!---------------------------------------------------------------------- 406 !! 407 !! subroutine parameter 414 408 REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: p_rnfmsk ! river runoff mask (rnfmsk array) 415 !!---------------------------------------------------------------------- 416 ! 417 IF( jncsr > 0 ) THEN 418 WHERE( closea_mask_rnf(:,:) > 0 .and. closea_mask(:,:) == 0 ) 419 p_rnfmsk(:,:) = MAX( p_rnfmsk(:,:), 1.0_wp ) 420 ENDWHERE 421 ENDIF 422 ! 423 IF( jncse > 0 ) THEN 424 WHERE( closea_mask_empmr(:,:) > 0 .and. closea_mask(:,:) == 0 ) 425 p_rnfmsk(:,:) = MAX( p_rnfmsk(:,:), 1.0_wp ) 426 ENDWHERE 427 ENDIF 409 !! 410 !! local variables 411 REAL(wp), DIMENSION(jpi,jpj) :: zmsk 412 !!---------------------------------------------------------------------- 413 ! 414 ! zmsk > 0 where cs river mouth defined (case rnf and emp) 415 zmsk(:,:) = ( mask_csgrprnf (:,:) + mask_csgrpemp(:,:) ) * mask_opnsea(:,:) 416 WHERE( zmsk(:,:) > 0 ) 417 p_rnfmsk(:,:) = 1.0_wp 418 END WHERE 428 419 ! 429 420 END SUBROUTINE clo_rnf 430 431 421 432 SUBROUTINE clo_bat( k_top, k_bot )422 SUBROUTINE clo_bat( k_top, k_bot, p_mask, p_prt ) 433 423 !!--------------------------------------------------------------------- 434 424 !! *** ROUTINE clo_bat *** … … 443 433 !! ** Action : set k_top=0 and k_bot=0 over closed seas 444 434 !!---------------------------------------------------------------------- 435 !! 436 !! subroutine parameter 445 437 INTEGER, DIMENSION(:,:), INTENT(inout) :: k_top, k_bot ! ocean first and last level indices 446 INTEGER :: inum, id 447 INTEGER, DIMENSION(jpi,jpj) :: closea_mask ! closea_mask field 448 REAL(wp), DIMENSION(jpi,jpj) :: zdata_in ! temporary real array for input 449 !!---------------------------------------------------------------------- 450 ! 451 IF(lwp) THEN ! Control print 438 INTEGER, DIMENSION(:,:), INTENT(in ) :: p_mask ! mask used to mask ktop and k_bot 439 CHARACTER(256), INTENT(in ) :: p_prt ! text for control print 440 !! 441 !! local variables 442 !!---------------------------------------------------------------------- 443 !! 444 IF ( lwp ) THEN 452 445 WRITE(numout,*) 453 WRITE(numout,*) 'clo_bat : suppression of closed seas'446 WRITE(numout,*) 'clo_bat : Suppression closed seas based on ',TRIM(p_prt),' field.' 454 447 WRITE(numout,*) '~~~~~~~' 448 WRITE(numout,*) 455 449 ENDIF 456 ! 457 IF( ln_read_cfg ) THEN 458 ! 459 CALL iom_open( cn_domcfg, inum ) 460 ! 461 id = iom_varid(inum, 'closea_mask', ldstop = .false.) 462 IF( id > 0 ) THEN 463 IF( lwp ) WRITE(numout,*) 'Suppressing closed seas in bathymetry based on closea_mask field,' 464 CALL iom_get ( inum, jpdom_data, 'closea_mask', zdata_in ) 465 closea_mask(:,:) = NINT(zdata_in(:,:)) 466 WHERE( closea_mask(:,:) > 0 ) 467 k_top(:,:) = 0 468 k_bot(:,:) = 0 469 ENDWHERE 470 ELSE 471 IF( lwp ) WRITE(numout,*) 'No closea_mask field found in domain_cfg file. No suppression of closed seas.' 472 ENDIF 473 ! 474 CALL iom_close(inum) 475 ! 476 ELSE 477 IF( lwp ) WRITE(numout,*) 'No domain_cfg file => no suppression of closed seas.' 478 ENDIF 479 ! 480 ! Initialise l_sbc_clo and l_clo_rnf for this case (ln_closea=.false.) 481 l_sbc_clo = .false. 482 l_clo_rnf = .false. 483 ! 450 !! 451 k_top(:,:) = k_top(:,:) * p_mask(:,:) 452 k_bot(:,:) = k_bot(:,:) * p_mask(:,:) 453 !! 484 454 END SUBROUTINE clo_bat 455 456 SUBROUTINE read_csmask(p_file, p_var, p_mskout) 457 ! 458 ! subroutine parameter 459 CHARACTER(256), INTENT(in ) :: p_file, p_var 460 INTEGER, DIMENSION(:,:), INTENT(inout) :: p_mskout 461 ! 462 ! local variables 463 INTEGER :: ics 464 REAL(wp), DIMENSION(jpi,jpj) :: zdta 465 ! 466 CALL iom_open ( p_file, ics ) 467 CALL iom_get ( ics, jpdom_data, TRIM(p_var), zdta ) 468 CALL iom_close( ics ) 469 p_mskout(:,:) = NINT(zdta(:,:)) 470 ! 471 END SUBROUTINE read_csmask 472 473 SUBROUTINE alloc_csmask( pmask ) 474 ! 475 ! subroutine parameter 476 INTEGER, ALLOCATABLE, DIMENSION(:,:), INTENT(inout) :: pmask 477 ! 478 ! local variables 479 INTEGER :: ierr 480 ! 481 ALLOCATE( pmask(jpi,jpj) , STAT=ierr ) 482 IF( ierr /= 0 ) CALL ctl_stop( 'STOP', 'alloc_csmask: failed to allocate surf array') 483 ! 484 END SUBROUTINE 485 486 487 SUBROUTINE alloc_cssurf( klen, pvarsrc, pvartrg ) 488 ! 489 ! subroutine parameter 490 INTEGER, INTENT(in) :: klen 491 REAL(wp), ALLOCATABLE, DIMENSION(:), INTENT(inout) :: pvarsrc, pvartrg 492 ! 493 ! local variables 494 INTEGER :: ierr 495 ! 496 ! klen (number of lake) can be zero so use MAX(klen,1) to avoid 0 length array 497 ALLOCATE( pvarsrc(MAX(klen,1)) , pvartrg(MAX(klen,1)) , STAT=ierr ) 498 IF( ierr /= 0 ) CALL ctl_stop( 'STOP', 'sbc_clo: failed to allocate surf array') 499 ! initialise to 0 500 pvarsrc(:) = 0.e0_wp 501 pvartrg(:) = 0.e0_wp 502 END SUBROUTINE 503 504 SUBROUTINE alloc_csgrp( klen, kvar ) 505 ! 506 ! subroutine parameter 507 INTEGER, INTENT(in) :: klen 508 INTEGER, ALLOCATABLE, DIMENSION(:), INTENT(inout) :: kvar 509 ! 510 ! local variables 511 INTEGER :: ierr 512 ! 513 ! klen (number of lake) can be zero so use MAX(klen,1) to avoid 0 length array 514 ALLOCATE( kvar(MAX(klen,1)) , STAT=ierr ) 515 IF( ierr /= 0 ) CALL ctl_stop( 'STOP', 'sbc_clo: failed to allocate group array') 516 ! initialise to 0 517 kvar(:) = 0 518 END SUBROUTINE 485 519 486 520 !!====================================================================== -
NEMO/branches/2019/ENHANCE-03_closea/src/OCE/DOM/domain.F90
r10425 r11207 134 134 ENDIF 135 135 ! 136 CALL dom_hgr ! Horizontal mesh 137 CALL dom_zgr( ik_top, ik_bot ) ! Vertical mesh and bathymetry 138 CALL dom_msk( ik_top, ik_bot ) ! Masks 139 IF( ln_closea ) CALL dom_clo ! ln_closea=T : closed seas included in the simulation 140 ! Read in masks to define closed seas and lakes 141 ! 142 DO jj = 1, jpj ! depth of the iceshelves 136 CALL dom_hgr ! Horizontal mesh 137 138 IF( ln_closea ) CALL dom_clo ! Read in masks to define closed seas and lakes 139 140 CALL dom_zgr( ik_top, ik_bot ) ! Vertical mesh and bathymetry 141 142 CALL dom_msk( ik_top, ik_bot ) ! Masks 143 ! 144 DO jj = 1, jpj ! depth of the iceshelves 143 145 DO ji = 1, jpi 144 146 ik = mikt(ji,jj) -
NEMO/branches/2019/ENHANCE-03_closea/src/OCE/DOM/domzgr.F90
r10425 r11207 118 118 ! Any closed seas (defined by closea_mask > 0 in domain_cfg file) to be filled 119 119 ! in at runtime if ln_closea=.false. 120 IF( .NOT.ln_closea ) CALL clo_bat( k_top, k_bot ) 120 IF( ln_closea ) THEN 121 IF ( ln_maskcs ) THEN 122 ! mask all the closed sea 123 CALL clo_bat( k_top, k_bot, mask_opnsea, 'mask_opensea' ) 124 ELSE IF ( ln_mask_csundef ) THEN 125 ! defined closed sea are kept 126 ! mask all the undefined closed sea 127 CALL clo_bat( k_top, k_bot, mask_csundef, 'mask_csundef' ) 128 END IF 129 END IF 121 130 ! 122 131 IF(lwp) THEN ! Control print
Note: See TracChangeset
for help on using the changeset viewer.