- Timestamp:
- 2016-01-08T10:35:19+01:00 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/NEMOGCM/NEMO/OFF_SRC/domrea.F90
r4569 r6225 1 1 MODULE domrea 2 !!====================================================================== 3 !! *** MODULE domrea *** 4 !! Ocean initialization : read the ocean domain meshmask file(s) 5 !!====================================================================== 6 !! History : 3.3 ! 2010-05 (C. Ethe) Full reorganization of the off-line 2 !!============================================================================== 3 !! *** MODULE domrea *** 4 !! Ocean initialization : domain initialization 5 !!============================================================================== 6 !! History : OPA ! 1990-10 (C. Levy - G. Madec) Original code 7 !! ! 1992-01 (M. Imbard) insert time step initialization 8 !! ! 1996-06 (G. Madec) generalized vertical coordinate 9 !! ! 1997-02 (G. Madec) creation of domwri.F 10 !! ! 2001-05 (E.Durand - G. Madec) insert closed sea 11 !! NEMO 1.0 ! 2002-08 (G. Madec) F90: Free form and module 7 12 !!---------------------------------------------------------------------- 8 13 9 14 !!---------------------------------------------------------------------- 10 !! dom_rea : read mesh and mask file(s) 11 !! nmsh = 1 : mesh_mask file 12 !! = 2 : mesh and mask file 13 !! = 3 : mesh_hgr, mesh_zgr and mask 15 !! dom_init : initialize the space and time domain 16 !! dom_nam : read and contral domain namelists 17 !! dom_ctl : control print for the ocean domain 14 18 !!---------------------------------------------------------------------- 19 USE oce ! 20 USE trc_oce ! shared ocean/biogeochemical variables 15 21 USE dom_oce ! ocean space and time domain 16 USE dommsk ! domain: masks 22 USE phycst ! physical constants 23 USE domstp ! domain: set the time-step 24 ! 25 USE in_out_manager ! I/O manager 26 USE lib_mpp ! distributed memory computing library 17 27 USE lbclnk ! lateral boundary condition - MPP exchanges 18 USE trc_oce ! shared ocean/biogeochemical variables19 USE lib_mpp20 USE in_out_manager21 28 USE wrk_nemo 22 29 23 30 IMPLICIT NONE 24 31 PRIVATE 25 32 26 PUBLIC dom_rea ! routine called by inidom.F90 27 !! * Substitutions 28 # include "domzgr_substitute.h90" 33 PUBLIC dom_rea ! called by nemogcm.F90 34 35 !! * Substitutions 36 # include "vectopt_loop_substitute.h90" 29 37 !!---------------------------------------------------------------------- 30 !! NEMO/OFF 3. 3 , NEMO Consortium (2010)38 !! NEMO/OFF 3.7 , NEMO Consortium (2015) 31 39 !! $Id$ 32 !! Software governed by the CeCILL licence 40 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 33 41 !!---------------------------------------------------------------------- 34 42 CONTAINS … … 37 45 !!---------------------------------------------------------------------- 38 46 !! *** ROUTINE dom_rea *** 47 !! 48 !! ** Purpose : Domain initialization. Call the routines that are 49 !! required to create the arrays which define the space and time 50 !! domain of the ocean model. 51 !! 52 !! ** Method : 53 !! - dom_stp: defined the model time step 54 !! - dom_rea: read the meshmask file if nmsh=1 55 !!---------------------------------------------------------------------- 56 INTEGER :: jk ! dummy loop index 57 INTEGER :: iconf = 0 ! local integers 58 !!---------------------------------------------------------------------- 59 ! 60 IF(lwp) THEN 61 WRITE(numout,*) 62 WRITE(numout,*) 'dom_init : domain initialization' 63 WRITE(numout,*) '~~~~~~~~' 64 ENDIF 65 ! 66 CALL dom_nam ! read namelist ( namrun, namdom ) 67 CALL dom_zgr ! Vertical mesh and bathymetry option 68 CALL dom_grd ! Create a domain file 69 ! 70 ! ! associated horizontal metrics 71 ! 72 r1_e1t(:,:) = 1._wp / e1t(:,:) ; r1_e2t (:,:) = 1._wp / e2t(:,:) 73 r1_e1u(:,:) = 1._wp / e1u(:,:) ; r1_e2u (:,:) = 1._wp / e2u(:,:) 74 r1_e1v(:,:) = 1._wp / e1v(:,:) ; r1_e2v (:,:) = 1._wp / e2v(:,:) 75 r1_e1f(:,:) = 1._wp / e1f(:,:) ; r1_e2f (:,:) = 1._wp / e2f(:,:) 76 ! 77 !!gm BUG if scale factor reduction !!!! 78 e1e2t (:,:) = e1t(:,:) * e2t(:,:) ; r1_e1e2t(:,:) = 1._wp / e1e2t(:,:) 79 e1e2u (:,:) = e1u(:,:) * e2u(:,:) ; r1_e1e2u(:,:) = 1._wp / e1e2u(:,:) 80 e1e2v (:,:) = e1v(:,:) * e2v(:,:) ; r1_e1e2v(:,:) = 1._wp / e1e2v(:,:) 81 e1e2f (:,:) = e1f(:,:) * e2f(:,:) ; r1_e1e2f(:,:) = 1._wp / e1e2f(:,:) 82 ! 83 e2_e1u(:,:) = e2u(:,:) / e1u(:,:) 84 e1_e2v(:,:) = e1v(:,:) / e2v(:,:) 85 ! 86 hu_n(:,:) = e3u_n(:,:,1) * umask(:,:,1) ! Ocean depth at U- and V-points 87 hv_n(:,:) = e3v_n(:,:,1) * vmask(:,:,1) 88 DO jk = 2, jpk 89 hu_n(:,:) = hu_n(:,:) + e3u_n(:,:,jk) * umask(:,:,jk) 90 hv_n(:,:) = hv_n(:,:) + e3v_n(:,:,jk) * vmask(:,:,jk) 91 END DO 92 ! ! Inverse of the local depth 93 r1_hu_n(:,:) = 1._wp / ( hu_n(:,:) + 1._wp - umask(:,:,1) ) * umask(:,:,1) 94 r1_hv_n(:,:) = 1._wp / ( hv_n(:,:) + 1._wp - vmask(:,:,1) ) * vmask(:,:,1) 95 ! 96 CALL dom_stp ! Time step 97 CALL dom_msk ! Masks 98 CALL dom_ctl ! Domain control 99 ! 100 END SUBROUTINE dom_rea 101 102 103 SUBROUTINE dom_nam 104 !!---------------------------------------------------------------------- 105 !! *** ROUTINE dom_nam *** 106 !! 107 !! ** Purpose : read domaine namelists and print the variables. 108 !! 109 !! ** input : - namrun namelist 110 !! - namdom namelist 111 !!---------------------------------------------------------------------- 112 USE ioipsl 113 INTEGER :: ios ! Local integer output status for namelist read 114 ! 115 NAMELIST/namrun/ cn_ocerst_indir, cn_ocerst_outdir, nn_stocklist, ln_rst_list, & 116 & nn_no , cn_exp , cn_ocerst_in, cn_ocerst_out, ln_rstart , nn_rstctl, & 117 & nn_it000, nn_itend , nn_date0 , nn_time0, nn_leapy , nn_istate , nn_stock , & 118 & nn_write, ln_iscpl , ln_mskland , ln_cfmeta , ln_clobber, nn_chunksz, nn_euler 119 NAMELIST/namdom/ nn_bathy , rn_bathy , rn_e3zps_min, rn_e3zps_rat , nn_msh , rn_hmin , rn_isfhmin,& 120 & rn_atfp , rn_rdt , nn_baro , nn_closea , ln_crs , jphgr_msh, & 121 & ppglam0, ppgphi0, ppe1_deg, ppe2_deg, ppe1_m, ppe2_m, & 122 & ppsur, ppa0, ppa1, ppkth, ppacr, ppdzmin, pphmax, ldbletanh, & 123 & ppa2, ppkth2, ppacr2 124 #if defined key_netcdf4 125 NAMELIST/namnc4/ nn_nchunks_i, nn_nchunks_j, nn_nchunks_k, ln_nc4zip 126 #endif 127 !!---------------------------------------------------------------------- 128 129 REWIND( numnam_ref ) ! Namelist namrun in reference namelist : Parameters of the run 130 READ ( numnam_ref, namrun, IOSTAT = ios, ERR = 901) 131 901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namrun in reference namelist', lwp ) 132 133 REWIND( numnam_cfg ) ! Namelist namrun in configuration namelist : Parameters of the run 134 READ ( numnam_cfg, namrun, IOSTAT = ios, ERR = 902 ) 135 902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namrun in configuration namelist', lwp ) 136 IF(lwm) WRITE ( numond, namrun ) 137 ! 138 IF(lwp) THEN ! control print 139 WRITE(numout,*) 140 WRITE(numout,*) 'dom_nam : domain initialization through namelist read' 141 WRITE(numout,*) '~~~~~~~ ' 142 WRITE(numout,*) ' Namelist namrun' 143 WRITE(numout,*) ' job number nn_no = ', nn_no 144 WRITE(numout,*) ' experiment name for output cn_exp = ', cn_exp 145 WRITE(numout,*) ' restart logical ln_rstart = ', ln_rstart 146 WRITE(numout,*) ' control of time step nn_rstctl = ', nn_rstctl 147 WRITE(numout,*) ' number of the first time step nn_it000 = ', nn_it000 148 WRITE(numout,*) ' number of the last time step nn_itend = ', nn_itend 149 WRITE(numout,*) ' initial calendar date aammjj nn_date0 = ', nn_date0 150 WRITE(numout,*) ' leap year calendar (0/1) nn_leapy = ', nn_leapy 151 WRITE(numout,*) ' initial state output nn_istate = ', nn_istate 152 WRITE(numout,*) ' frequency of restart file nn_stock = ', nn_stock 153 WRITE(numout,*) ' frequency of output file nn_write = ', nn_write 154 WRITE(numout,*) ' mask land points ln_mskland = ', ln_mskland 155 WRITE(numout,*) ' additional CF standard metadata ln_cfmeta = ', ln_cfmeta 156 WRITE(numout,*) ' overwrite an existing file ln_clobber = ', ln_clobber 157 WRITE(numout,*) ' NetCDF chunksize (bytes) nn_chunksz = ', nn_chunksz 158 ENDIF 159 no = nn_no ! conversion DOCTOR names into model names (this should disappear soon) 160 cexper = cn_exp 161 nrstdt = nn_rstctl 162 nit000 = nn_it000 163 nitend = nn_itend 164 ndate0 = nn_date0 165 nleapy = nn_leapy 166 ninist = nn_istate 167 nstock = nn_stock 168 nstocklist = nn_stocklist 169 nwrite = nn_write 170 ! 171 ! ! control of output frequency 172 IF ( nstock == 0 .OR. nstock > nitend ) THEN 173 WRITE(ctmp1,*) 'nstock = ', nstock, ' it is forced to ', nitend 174 CALL ctl_warn( ctmp1 ) 175 nstock = nitend 176 ENDIF 177 IF ( nwrite == 0 ) THEN 178 WRITE(ctmp1,*) 'nwrite = ', nwrite, ' it is forced to ', nitend 179 CALL ctl_warn( ctmp1 ) 180 nwrite = nitend 181 ENDIF 182 183 ! parameters correspondting to nit000 - 1 (as we start the step loop with a call to day) 184 ndastp = ndate0 - 1 ! ndate0 read in the namelist in dom_nam, we assume that we start run at 00:00 185 adatrj = ( REAL( nit000-1, wp ) * rdt ) / rday 186 187 #if defined key_agrif 188 IF( Agrif_Root() ) THEN 189 #endif 190 SELECT CASE ( nleapy ) ! Choose calendar for IOIPSL 191 CASE ( 1 ) 192 CALL ioconf_calendar('gregorian') 193 IF(lwp) WRITE(numout,*) ' The IOIPSL calendar is "gregorian", i.e. leap year' 194 CASE ( 0 ) 195 CALL ioconf_calendar('noleap') 196 IF(lwp) WRITE(numout,*) ' The IOIPSL calendar is "noleap", i.e. no leap year' 197 CASE ( 30 ) 198 CALL ioconf_calendar('360d') 199 IF(lwp) WRITE(numout,*) ' The IOIPSL calendar is "360d", i.e. 360 days in a year' 200 END SELECT 201 #if defined key_agrif 202 ENDIF 203 #endif 204 205 REWIND( numnam_ref ) ! Namelist namdom in reference namelist : space & time domain (bathymetry, mesh, timestep) 206 READ ( numnam_ref, namdom, IOSTAT = ios, ERR = 903) 207 903 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namdom in reference namelist', lwp ) 208 209 REWIND( numnam_cfg ) ! Namelist namdom in configuration namelist : space & time domain (bathymetry, mesh, timestep) 210 READ ( numnam_cfg, namdom, IOSTAT = ios, ERR = 904 ) 211 904 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namdom in configuration namelist', lwp ) 212 IF(lwm) WRITE ( numond, namdom ) 213 214 IF(lwp) THEN 215 WRITE(numout,*) 216 WRITE(numout,*) ' Namelist namdom : space & time domain' 217 WRITE(numout,*) ' flag read/compute bathymetry nn_bathy = ', nn_bathy 218 WRITE(numout,*) ' Depth (if =0 bathy=jpkm1) rn_bathy = ', rn_bathy 219 WRITE(numout,*) ' min depth of the ocean (>0) or rn_hmin = ', rn_hmin 220 WRITE(numout,*) ' minimum thickness of partial rn_e3zps_min = ', rn_e3zps_min, ' (m)' 221 WRITE(numout,*) ' step level rn_e3zps_rat = ', rn_e3zps_rat 222 WRITE(numout,*) ' create mesh/mask file(s) nn_msh = ', nn_msh 223 WRITE(numout,*) ' = 0 no file created ' 224 WRITE(numout,*) ' = 1 mesh_mask ' 225 WRITE(numout,*) ' = 2 mesh and mask ' 226 WRITE(numout,*) ' = 3 mesh_hgr, msh_zgr and mask ' 227 WRITE(numout,*) ' ocean time step rn_rdt = ', rn_rdt 228 WRITE(numout,*) ' asselin time filter parameter rn_atfp = ', rn_atfp 229 WRITE(numout,*) ' time-splitting: nb of sub time-step nn_baro = ', nn_baro 230 WRITE(numout,*) ' suppression of closed seas (=0) nn_closea = ', nn_closea 231 WRITE(numout,*) ' type of horizontal mesh jphgr_msh = ', jphgr_msh 232 WRITE(numout,*) ' longitude of first raw and column T-point ppglam0 = ', ppglam0 233 WRITE(numout,*) ' latitude of first raw and column T-point ppgphi0 = ', ppgphi0 234 WRITE(numout,*) ' zonal grid-spacing (degrees) ppe1_deg = ', ppe1_deg 235 WRITE(numout,*) ' meridional grid-spacing (degrees) ppe2_deg = ', ppe2_deg 236 WRITE(numout,*) ' zonal grid-spacing (degrees) ppe1_m = ', ppe1_m 237 WRITE(numout,*) ' meridional grid-spacing (degrees) ppe2_m = ', ppe2_m 238 WRITE(numout,*) ' ORCA r4, r2 and r05 coefficients ppsur = ', ppsur 239 WRITE(numout,*) ' ppa0 = ', ppa0 240 WRITE(numout,*) ' ppa1 = ', ppa1 241 WRITE(numout,*) ' ppkth = ', ppkth 242 WRITE(numout,*) ' ppacr = ', ppacr 243 WRITE(numout,*) ' Minimum vertical spacing ppdzmin = ', ppdzmin 244 WRITE(numout,*) ' Maximum depth pphmax = ', pphmax 245 WRITE(numout,*) ' Use double tanf function for vertical coordinates ldbletanh = ', ldbletanh 246 WRITE(numout,*) ' Double tanh function parameters ppa2 = ', ppa2 247 WRITE(numout,*) ' ppkth2 = ', ppkth2 248 WRITE(numout,*) ' ppacr2 = ', ppacr2 249 ENDIF 250 251 ntopo = nn_bathy ! conversion DOCTOR names into model names (this should disappear soon) 252 e3zps_min = rn_e3zps_min 253 e3zps_rat = rn_e3zps_rat 254 nmsh = nn_msh 255 atfp = rn_atfp 256 rdt = rn_rdt 257 #if defined key_netcdf4 258 ! ! NetCDF 4 case ("key_netcdf4" defined) 259 REWIND( numnam_ref ) ! Namelist namnc4 in reference namelist : NETCDF 260 READ ( numnam_ref, namnc4, IOSTAT = ios, ERR = 907) 261 907 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namnc4 in reference namelist', lwp ) 262 263 REWIND( numnam_cfg ) ! Namelist namnc4 in configuration namelist : NETCDF 264 READ ( numnam_cfg, namnc4, IOSTAT = ios, ERR = 908 ) 265 908 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namnc4 in configuration namelist', lwp ) 266 IF(lwm) WRITE( numond, namnc4 ) 267 IF(lwp) THEN ! control print 268 WRITE(numout,*) 269 WRITE(numout,*) ' Namelist namnc4 - Netcdf4 chunking parameters' 270 WRITE(numout,*) ' number of chunks in i-dimension nn_nchunks_i = ', nn_nchunks_i 271 WRITE(numout,*) ' number of chunks in j-dimension nn_nchunks_j = ', nn_nchunks_j 272 WRITE(numout,*) ' number of chunks in k-dimension nn_nchunks_k = ', nn_nchunks_k 273 WRITE(numout,*) ' apply netcdf4/hdf5 chunking & compression ln_nc4zip = ', ln_nc4zip 274 ENDIF 275 276 ! Put the netcdf4 settings into a simple structure (snc4set, defined in in_out_manager module) 277 ! Note the chunk size in the unlimited (time) dimension will be fixed at 1 278 snc4set%ni = nn_nchunks_i 279 snc4set%nj = nn_nchunks_j 280 snc4set%nk = nn_nchunks_k 281 snc4set%luse = ln_nc4zip 282 #else 283 snc4set%luse = .FALSE. ! No NetCDF 4 case 284 #endif 285 ! 286 END SUBROUTINE dom_nam 287 288 289 SUBROUTINE dom_zgr 290 !!---------------------------------------------------------------------- 291 !! *** ROUTINE dom_zgr *** 292 !! 293 !! ** Purpose : set the depth of model levels and the resulting 294 !! vertical scale factors. 295 !! 296 !! ** Method : - reference 1D vertical coordinate (gdep._1d, e3._1d) 297 !! - read/set ocean depth and ocean levels (bathy, mbathy) 298 !! - vertical coordinate (gdep., e3.) depending on the 299 !! coordinate chosen : 300 !! ln_zco=T z-coordinate 301 !! ln_zps=T z-coordinate with partial steps 302 !! ln_zco=T s-coordinate 303 !! 304 !! ** Action : define gdep., e3., mbathy and bathy 305 !!---------------------------------------------------------------------- 306 INTEGER :: ioptio = 0 ! temporary integer 307 INTEGER :: ios 308 !! 309 NAMELIST/namzgr/ ln_zco, ln_zps, ln_sco, ln_isfcav, ln_linssh 310 !!---------------------------------------------------------------------- 311 312 REWIND( numnam_ref ) ! Namelist namzgr in reference namelist : Vertical coordinate 313 READ ( numnam_ref, namzgr, IOSTAT = ios, ERR = 901 ) 314 901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namzgr in reference namelist', lwp ) 315 316 REWIND( numnam_cfg ) ! Namelist namzgr in configuration namelist : Vertical coordinate 317 READ ( numnam_cfg, namzgr, IOSTAT = ios, ERR = 902 ) 318 902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namzgr in configuration namelist', lwp ) 319 IF(lwm) WRITE ( numond, namzgr ) 320 321 IF(lwp) THEN ! Control print 322 WRITE(numout,*) 323 WRITE(numout,*) 'dom_zgr : vertical coordinate' 324 WRITE(numout,*) '~~~~~~~' 325 WRITE(numout,*) ' Namelist namzgr : set vertical coordinate' 326 WRITE(numout,*) ' z-coordinate - full steps ln_zco = ', ln_zco 327 WRITE(numout,*) ' z-coordinate - partial steps ln_zps = ', ln_zps 328 WRITE(numout,*) ' s- or hybrid z-s-coordinate ln_sco = ', ln_sco 329 WRITE(numout,*) ' ice shelf cavity ln_isfcav = ', ln_isfcav 330 WRITE(numout,*) ' Linear free surface ln_linssh = ', ln_linssh 331 ENDIF 332 333 ioptio = 0 ! Check Vertical coordinate options 334 IF( ln_zco ) ioptio = ioptio + 1 335 IF( ln_zps ) ioptio = ioptio + 1 336 IF( ln_sco ) ioptio = ioptio + 1 337 IF( ln_isfcav ) ioptio = 33 338 IF ( ioptio /= 1 ) CALL ctl_stop( ' none or several vertical coordinate options used' ) 339 IF ( ioptio == 33 ) CALL ctl_stop( ' isf cavity with off line module not yet done ' ) 340 341 END SUBROUTINE dom_zgr 342 343 344 SUBROUTINE dom_ctl 345 !!---------------------------------------------------------------------- 346 !! *** ROUTINE dom_ctl *** 347 !! 348 !! ** Purpose : Domain control. 349 !! 350 !! ** Method : compute and print extrema of masked scale factors 351 !! 352 !!---------------------------------------------------------------------- 353 INTEGER :: iimi1, ijmi1, iimi2, ijmi2, iima1, ijma1, iima2, ijma2 354 INTEGER, DIMENSION(2) :: iloc ! 355 REAL(wp) :: ze1min, ze1max, ze2min, ze2max 356 !!---------------------------------------------------------------------- 357 358 ! Extrema of the scale factors 359 360 IF(lwp)WRITE(numout,*) 361 IF(lwp)WRITE(numout,*) 'dom_ctl : extrema of the masked scale factors' 362 IF(lwp)WRITE(numout,*) '~~~~~~~' 363 364 IF (lk_mpp) THEN 365 CALL mpp_minloc( e1t(:,:), tmask(:,:,1), ze1min, iimi1,ijmi1 ) 366 CALL mpp_minloc( e2t(:,:), tmask(:,:,1), ze2min, iimi2,ijmi2 ) 367 CALL mpp_maxloc( e1t(:,:), tmask(:,:,1), ze1max, iima1,ijma1 ) 368 CALL mpp_maxloc( e2t(:,:), tmask(:,:,1), ze2max, iima2,ijma2 ) 369 ELSE 370 ze1min = MINVAL( e1t(:,:), mask = tmask(:,:,1) == 1.e0 ) 371 ze2min = MINVAL( e2t(:,:), mask = tmask(:,:,1) == 1.e0 ) 372 ze1max = MAXVAL( e1t(:,:), mask = tmask(:,:,1) == 1.e0 ) 373 ze2max = MAXVAL( e2t(:,:), mask = tmask(:,:,1) == 1.e0 ) 374 375 iloc = MINLOC( e1t(:,:), mask = tmask(:,:,1) == 1.e0 ) 376 iimi1 = iloc(1) + nimpp - 1 377 ijmi1 = iloc(2) + njmpp - 1 378 iloc = MINLOC( e2t(:,:), mask = tmask(:,:,1) == 1.e0 ) 379 iimi2 = iloc(1) + nimpp - 1 380 ijmi2 = iloc(2) + njmpp - 1 381 iloc = MAXLOC( e1t(:,:), mask = tmask(:,:,1) == 1.e0 ) 382 iima1 = iloc(1) + nimpp - 1 383 ijma1 = iloc(2) + njmpp - 1 384 iloc = MAXLOC( e2t(:,:), mask = tmask(:,:,1) == 1.e0 ) 385 iima2 = iloc(1) + nimpp - 1 386 ijma2 = iloc(2) + njmpp - 1 387 ENDIF 388 ! 389 IF(lwp) THEN 390 WRITE(numout,"(14x,'e1t maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") ze1max, iima1, ijma1 391 WRITE(numout,"(14x,'e1t mini: ',1f10.2,' at i = ',i5,' j= ',i5)") ze1min, iimi1, ijmi1 392 WRITE(numout,"(14x,'e2t maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") ze2max, iima2, ijma2 393 WRITE(numout,"(14x,'e2t mini: ',1f10.2,' at i = ',i5,' j= ',i5)") ze2min, iimi2, ijmi2 394 ENDIF 395 ! 396 END SUBROUTINE dom_ctl 397 398 399 SUBROUTINE dom_grd 400 !!---------------------------------------------------------------------- 401 !! *** ROUTINE dom_grd *** 39 402 !! 40 403 !! ** Purpose : Read the NetCDF file(s) which contain(s) all the … … 141 504 CALL iom_get( inum2, jpdom_data, 'facvolt', facvol ) 142 505 #endif 143 144 506 ! ! horizontal mesh (inum3) 145 507 CALL iom_get( inum3, jpdom_data, 'glamt', glamt ) … … 164 526 165 527 CALL iom_get( inum4, jpdom_data, 'mbathy', zmbk ) ! number of ocean t-points 166 mbathy(:,:) = INT( zmbk(:,:) ) 528 mbathy (:,:) = INT( zmbk(:,:) ) 529 misfdep(:,:) = 1 ! ice shelf case not yet done 167 530 168 531 CALL zgr_bot_level ! mbk. arrays (deepest ocean t-, u- & v-points … … 180 543 CALL iom_get( inum4, jpdom_unknown, 'esigw', esigw ) 181 544 182 CALL iom_get( inum4, jpdom_data, 'e3t_0', fse3t_n(:,:,:) ) ! scale factors183 CALL iom_get( inum4, jpdom_data, 'e3u_0', fse3u_n(:,:,:) )184 CALL iom_get( inum4, jpdom_data, 'e3v_0', fse3v_n(:,:,:) )185 CALL iom_get( inum4, jpdom_data, 'e3w_0', fse3w_n(:,:,:) )545 CALL iom_get( inum4, jpdom_data, 'e3t_0', e3t_0(:,:,:) ) ! scale factors 546 CALL iom_get( inum4, jpdom_data, 'e3u_0', e3u_0(:,:,:) ) 547 CALL iom_get( inum4, jpdom_data, 'e3v_0', e3v_0(:,:,:) ) 548 CALL iom_get( inum4, jpdom_data, 'e3w_0', e3w_0(:,:,:) ) 186 549 187 550 CALL iom_get( inum4, jpdom_unknown, 'gdept_1d', gdept_1d ) ! depth … … 197 560 ! 198 561 IF( nmsh <= 6 ) THEN ! 3D vertical scale factors 199 CALL iom_get( inum4, jpdom_data, 'e3t_0', fse3t_n(:,:,:) )200 CALL iom_get( inum4, jpdom_data, 'e3u_0', fse3u_n(:,:,:) )201 CALL iom_get( inum4, jpdom_data, 'e3v_0', fse3v_n(:,:,:) )202 CALL iom_get( inum4, jpdom_data, 'e3w_0', fse3w_n(:,:,:) )562 CALL iom_get( inum4, jpdom_data, 'e3t_0', e3t_0(:,:,:) ) 563 CALL iom_get( inum4, jpdom_data, 'e3u_0', e3u_0(:,:,:) ) 564 CALL iom_get( inum4, jpdom_data, 'e3v_0', e3v_0(:,:,:) ) 565 CALL iom_get( inum4, jpdom_data, 'e3w_0', e3w_0(:,:,:) ) 203 566 ELSE ! 2D bottom scale factors 204 567 CALL iom_get( inum4, jpdom_data, 'e3t_ps', e3tp ) … … 206 569 ! ! deduces the 3D scale factors 207 570 DO jk = 1, jpk 208 fse3t_n(:,:,jk) = e3t_1d(jk) ! set to the ref. factors209 fse3u_n(:,:,jk) = e3t_1d(jk)210 fse3v_n(:,:,jk) = e3t_1d(jk)211 fse3w_n(:,:,jk) = e3w_1d(jk)571 e3t_0(:,:,jk) = e3t_1d(jk) ! set to the ref. factors 572 e3u_0(:,:,jk) = e3t_1d(jk) 573 e3v_0(:,:,jk) = e3t_1d(jk) 574 e3w_0(:,:,jk) = e3w_1d(jk) 212 575 END DO 213 576 DO jj = 1,jpj ! adjust the deepest values 214 577 DO ji = 1,jpi 215 578 ik = mbkt(ji,jj) 216 fse3t_n(ji,jj,ik) = e3tp(ji,jj) * tmask(ji,jj,1) + e3t_1d(1) * ( 1._wp - tmask(ji,jj,1) )217 fse3w_n(ji,jj,ik) = e3wp(ji,jj) * tmask(ji,jj,1) + e3w_1d(1) * ( 1._wp - tmask(ji,jj,1) )579 e3t_0(ji,jj,ik) = e3tp(ji,jj) * tmask(ji,jj,1) + e3t_1d(1) * ( 1._wp - tmask(ji,jj,1) ) 580 e3w_0(ji,jj,ik) = e3wp(ji,jj) * tmask(ji,jj,1) + e3w_1d(1) * ( 1._wp - tmask(ji,jj,1) ) 218 581 END DO 219 582 END DO … … 221 584 DO jj = 1, jpjm1 222 585 DO ji = 1, jpim1 223 fse3u_n(ji,jj,jk) = MIN( fse3t_n(ji,jj,jk), fse3t_n(ji+1,jj,jk) )224 fse3v_n(ji,jj,jk) = MIN( fse3t_n(ji,jj,jk), fse3t_n(ji,jj+1,jk) )586 e3u_0(ji,jj,jk) = MIN( e3t_0(ji,jj,jk), e3t_0(ji+1,jj,jk) ) 587 e3v_0(ji,jj,jk) = MIN( e3t_0(ji,jj,jk), e3t_0(ji,jj+1,jk) ) 225 588 END DO 226 589 END DO 227 590 END DO 228 CALL lbc_lnk( fse3u_n(:,:,:) , 'U', 1._wp ) ; CALL lbc_lnk( fse3uw_n(:,:,:), 'U', 1._wp ) ! lateral boundary conditions229 CALL lbc_lnk( fse3v_n(:,:,:) , 'V', 1._wp ) ; CALL lbc_lnk( fse3vw_n(:,:,:), 'V', 1._wp )591 CALL lbc_lnk( e3u_0(:,:,:) , 'U', 1._wp ) ; CALL lbc_lnk( e3uw_0(:,:,:), 'U', 1._wp ) ! lateral boundary conditions 592 CALL lbc_lnk( e3v_0(:,:,:) , 'V', 1._wp ) ; CALL lbc_lnk( e3vw_0(:,:,:), 'V', 1._wp ) 230 593 ! 231 594 DO jk = 1, jpk ! set to z-scale factor if zero (i.e. along closed boundaries) 232 WHERE( fse3u_n(:,:,jk) == 0._wp ) fse3u_n(:,:,jk) = e3t_1d(jk)233 WHERE( fse3v_n(:,:,jk) == 0._wp ) fse3v_n(:,:,jk) = e3t_1d(jk)595 WHERE( e3u_0(:,:,jk) == 0._wp ) e3u_0(:,:,jk) = e3t_1d(jk) 596 WHERE( e3v_0(:,:,jk) == 0._wp ) e3v_0(:,:,jk) = e3t_1d(jk) 234 597 END DO 235 598 END IF 236 599 237 600 IF( iom_varid( inum4, 'gdept_0', ldstop = .FALSE. ) > 0 ) THEN ! 3D depth of t- and w-level 238 CALL iom_get( inum4, jpdom_data, 'gdept_0', fsdept_n(:,:,:) )239 CALL iom_get( inum4, jpdom_data, 'gdepw_0', fsdepw_n(:,:,:) )601 CALL iom_get( inum4, jpdom_data, 'gdept_0', gdept_0(:,:,:) ) 602 CALL iom_get( inum4, jpdom_data, 'gdepw_0', gdepw_0(:,:,:) ) 240 603 ELSE ! 2D bottom depth 241 604 CALL iom_get( inum4, jpdom_data, 'hdept', zprt ) … … 243 606 ! 244 607 DO jk = 1, jpk ! deduces the 3D depth 245 fsdept_n(:,:,jk) = gdept_1d(jk)246 fsdepw_n(:,:,jk) = gdepw_1d(jk)608 gdept_0(:,:,jk) = gdept_1d(jk) 609 gdepw_0(:,:,jk) = gdepw_1d(jk) 247 610 END DO 248 611 DO jj = 1, jpj … … 250 613 ik = mbkt(ji,jj) 251 614 IF( ik > 0 ) THEN 252 fsdepw_n(ji,jj,ik+1) = zprw(ji,jj)253 fsdept_n(ji,jj,ik ) = zprt(ji,jj)254 fsdept_n(ji,jj,ik+1) = fsdept_n(ji,jj,ik) + fse3t_n(ji,jj,ik)615 gdepw_0(ji,jj,ik+1) = zprw(ji,jj) 616 gdept_0(ji,jj,ik ) = zprt(ji,jj) 617 gdept_0(ji,jj,ik+1) = gdept_0(ji,jj,ik) + e3t_0(ji,jj,ik) 255 618 ENDIF 256 619 END DO … … 266 629 CALL iom_get( inum4, jpdom_unknown, 'e3w_1d' , e3w_1d ) 267 630 DO jk = 1, jpk 268 fse3t_n(:,:,jk) = e3t_1d(jk) ! set to the ref. factors269 fse3u_n(:,:,jk) = e3t_1d(jk)270 fse3v_n(:,:,jk) = e3t_1d(jk)271 fse3w_n(:,:,jk) = e3w_1d(jk)272 fsdept_n(:,:,jk) = gdept_1d(jk)273 fsdepw_n(:,:,jk) = gdepw_1d(jk)631 e3t_0(:,:,jk) = e3t_1d(jk) ! set to the ref. factors 632 e3u_0(:,:,jk) = e3t_1d(jk) 633 e3v_0(:,:,jk) = e3t_1d(jk) 634 e3w_0(:,:,jk) = e3w_1d(jk) 635 gdept_0(:,:,jk) = gdept_1d(jk) 636 gdepw_0(:,:,jk) = gdepw_1d(jk) 274 637 END DO 275 638 ENDIF 639 640 ! 641 ! !== time varying part of coordinate system ==! 642 ! 643 ! before ! now ! after ! 644 ; gdept_b = gdept_0 ; gdept_n = gdept_0 ! --- ! depth of grid-points 645 ; gdepw_b = gdepw_0 ; gdepw_n = gdepw_0 ! --- ! 646 ; ; gde3w_n = gde3w_0 ! --- ! 647 ! 648 ; e3t_b = e3t_0 ; e3t_n = e3t_0 ; e3t_a = e3t_0 ! scale factors 649 ; e3u_b = e3u_0 ; e3u_n = e3u_0 ; e3u_a = e3u_0 ! 650 ; e3v_b = e3v_0 ; e3v_n = e3v_0 ; e3v_a = e3v_0 ! 651 ; ; e3f_n = e3f_0 ! --- ! 652 ; e3w_b = e3w_0 ; e3w_n = e3w_0 ! --- ! 653 ; e3uw_b = e3uw_0 ; e3uw_n = e3uw_0 ! --- ! 654 ; e3vw_b = e3vw_0 ; e3vw_n = e3vw_0 ! --- ! 655 ! 276 656 277 657 !!gm BUG in s-coordinate this does not work! … … 303 683 & e2t (1,jj), e2u (1,jj), & 304 684 & e2v (1,jj), jj = 1, jpj, 10 ) 305 ENDIF306 307 308 IF( nprint == 1 .AND. lwp ) THEN309 WRITE(numout,*) ' e1u e2u '310 CALL prihre( e1u,jpi,jpj,jpi-5,jpi,1,jpj-5,jpj,1,0.,numout )311 CALL prihre( e2u,jpi,jpj,jpi-5,jpi,1,jpj-5,jpj,1,0.,numout )312 WRITE(numout,*) ' e1v e2v '313 CALL prihre( e1v,jpi,jpj,jpi-5,jpi,1,jpj-5,jpj,1,0.,numout )314 CALL prihre( e2v,jpi,jpj,jpi-5,jpi,1,jpj-5,jpj,1,0.,numout )315 685 ENDIF 316 686 … … 343 713 CALL wrk_dealloc( jpi, jpj, zmbk, zprt, zprw ) 344 714 ! 345 END SUBROUTINE dom_ rea715 END SUBROUTINE dom_grd 346 716 347 717 … … 358 728 !! (min value = 1 over land) 359 729 !!---------------------------------------------------------------------- 360 !361 730 INTEGER :: ji, jj ! dummy loop indices 362 731 REAL(wp), POINTER, DIMENSION(:,:) :: zmbk … … 371 740 ! 372 741 mbkt(:,:) = MAX( mbathy(:,:) , 1 ) ! bottom k-index of T-level (=1 over land) 742 mikt(:,:) = 1 ; miku(:,:) = 1; mikv(:,:) = 1; ! top k-index of T-level (=1 over open ocean; >1 beneath ice shelf) 373 743 ! ! bottom k-index of W-level = mbkt+1 374 744 DO jj = 1, jpjm1 ! bottom k-index of u- (v-) level … … 386 756 END SUBROUTINE zgr_bot_level 387 757 758 759 SUBROUTINE dom_msk 760 !!--------------------------------------------------------------------- 761 !! *** ROUTINE dom_msk *** 762 !! 763 !! ** Purpose : Off-line case: defines the interior domain T-mask. 764 !! 765 !! ** Method : The interior ocean/land mask is computed from tmask 766 !! setting to zero the duplicated row and lines due to 767 !! MPP exchange halos, est-west cyclic and north fold 768 !! boundary conditions. 769 !! 770 !! ** Action : tmask_i : interiorland/ocean mask at t-point 771 !! tpol : ??? 772 !!---------------------------------------------------------------------- 773 INTEGER :: ji, jj, jk ! dummy loop indices 774 INTEGER :: iif, iil, ijf, ijl ! local integers 775 INTEGER, POINTER, DIMENSION(:,:) :: imsk 776 !!--------------------------------------------------------------------- 777 778 CALL wrk_alloc( jpi, jpj, imsk ) 779 ! 780 ! Interior domain mask (used for global sum) 781 ! -------------------- 782 ssmask(:,:) = tmask(:,:,1) 783 tmask_i(:,:) = tmask(:,:,1) 784 iif = jpreci ! thickness of exchange halos in i-axis 785 iil = nlci - jpreci + 1 786 ijf = jprecj ! thickness of exchange halos in j-axis 787 ijl = nlcj - jprecj + 1 788 ! 789 tmask_i( 1 :iif, : ) = 0._wp ! first columns 790 tmask_i(iil:jpi, : ) = 0._wp ! last columns (including mpp extra columns) 791 tmask_i( : , 1 :ijf) = 0._wp ! first rows 792 tmask_i( : ,ijl:jpj) = 0._wp ! last rows (including mpp extra rows) 793 ! 794 ! ! north fold mask 795 tpol(1:jpiglo) = 1._wp 796 ! 797 IF( jperio == 3 .OR. jperio == 4 ) tpol(jpiglo/2+1:jpiglo) = 0._wp ! T-point pivot 798 IF( jperio == 5 .OR. jperio == 6 ) tpol( 1 :jpiglo) = 0._wp ! F-point pivot 799 IF( jperio == 3 .OR. jperio == 4 ) THEN ! T-point pivot: only half of the nlcj-1 row 800 IF( mjg(ijl-1) == jpjglo-1 ) THEN 801 DO ji = iif+1, iil-1 802 tmask_i(ji,ijl-1) = tmask_i(ji,ijl-1) * tpol(mig(ji)) 803 END DO 804 ENDIF 805 ENDIF 806 ! 807 ! (ISF) MIN(1,SUM(umask)) is here to check if you have effectively at 808 ! least 1 wet u point 809 DO jj = 1, jpjm1 810 DO ji = 1, fs_jpim1 ! vector loop 811 ssumask(ji,jj) = ssmask(ji,jj) * ssmask(ji+1,jj ) * MIN(1._wp,SUM(umask(ji,jj,:))) 812 ssvmask(ji,jj) = ssmask(ji,jj) * ssmask(ji ,jj+1) * MIN(1._wp,SUM(vmask(ji,jj,:))) 813 END DO 814 DO ji = 1, jpim1 ! NO vector opt. 815 ssfmask(ji,jj) = ssmask(ji,jj ) * ssmask(ji+1,jj ) & 816 & * ssmask(ji,jj+1) * ssmask(ji+1,jj+1) * MIN(1._wp,SUM(fmask(ji,jj,:))) 817 END DO 818 END DO 819 CALL lbc_lnk( ssumask, 'U', 1._wp ) ! Lateral boundary conditions 820 CALL lbc_lnk( ssvmask, 'V', 1._wp ) 821 CALL lbc_lnk( ssfmask, 'F', 1._wp ) 822 823 ! 3. Ocean/land mask at wu-, wv- and w points 824 !---------------------------------------------- 825 wmask (:,:,1) = tmask(:,:,1) ! surface value 826 wumask(:,:,1) = umask(:,:,1) 827 wvmask(:,:,1) = vmask(:,:,1) 828 DO jk = 2, jpk ! deeper value 829 wmask (:,:,jk) = tmask(:,:,jk) * tmask(:,:,jk-1) 830 wumask(:,:,jk) = umask(:,:,jk) * umask(:,:,jk-1) 831 wvmask(:,:,jk) = vmask(:,:,jk) * vmask(:,:,jk-1) 832 END DO 833 ! 834 CALL wrk_dealloc( jpi, jpj, imsk ) 835 ! 836 END SUBROUTINE dom_msk 837 388 838 !!====================================================================== 389 839 END MODULE domrea 840
Note: See TracChangeset
for help on using the changeset viewer.
