- Timestamp:
- 2016-06-06T07:57:00+02:00 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2016/dev_r6409_SIMPLIF_2_usrdef/NEMOGCM/NEMO/OPA_SRC/DOM/dommsk.F90
r6140 r6667 9 9 !! - ! 1996-05 (G. Madec) mask computed from tmask 10 10 !! 8.0 ! 1997-02 (G. Madec) mesh information put in domhgr.F 11 !! 8.1 ! 1997-07 (G. Madec) modification of mbathyand fmask11 !! 8.1 ! 1997-07 (G. Madec) modification of kbat and fmask 12 12 !! - ! 1998-05 (G. Roullet) free surface 13 13 !! 8.2 ! 2000-03 (G. Madec) no slip accurate … … 50 50 CONTAINS 51 51 52 SUBROUTINE dom_msk 52 SUBROUTINE dom_msk( k_top, k_bot ) 53 53 !!--------------------------------------------------------------------- 54 54 !! *** ROUTINE dom_msk *** … … 57 57 !! zontal velocity points (u & v), vorticity points (f) points. 58 58 !! 59 !! ** Method : The ocean/land mask is computed from the basin bathy- 60 !! metry in level (mbathy) which is defined or read in dommba. 61 !! mbathy equals 0 over continental T-point 62 !! and the number of ocean level over the ocean. 63 !! 64 !! At a given position (ji,jj,jk) the ocean/land mask is given by: 65 !! t-point : 0. IF mbathy( ji ,jj) =< 0 66 !! 1. IF mbathy( ji ,jj) >= jk 67 !! u-point : 0. IF mbathy( ji ,jj) or mbathy(ji+1, jj ) =< 0 68 !! 1. IF mbathy( ji ,jj) and mbathy(ji+1, jj ) >= jk. 69 !! v-point : 0. IF mbathy( ji ,jj) or mbathy( ji ,jj+1) =< 0 70 !! 1. IF mbathy( ji ,jj) and mbathy( ji ,jj+1) >= jk. 71 !! f-point : 0. IF mbathy( ji ,jj) or mbathy( ji ,jj+1) 72 !! or mbathy(ji+1,jj) or mbathy(ji+1,jj+1) =< 0 73 !! 1. IF mbathy( ji ,jj) and mbathy( ji ,jj+1) 74 !! and mbathy(ji+1,jj) and mbathy(ji+1,jj+1) >= jk. 75 !! tmask_i : interior ocean mask at t-point, i.e. excluding duplicated 76 !! rows/lines due to cyclic or North Fold boundaries as well 77 !! as MPP halos. 78 !! 79 !! The lateral friction is set through the value of fmask along 80 !! the coast and topography. This value is defined by rn_shlat, a 81 !! namelist parameter: 59 !! ** Method : The ocean/land mask at t-point is deduced from ko_top 60 !! and ko_bot, the indices of the fist and last ocean t-levels which 61 !! are either defined in usrdef_zgr or read in zgr_read. 62 !! The velocity masks (umask, vmask, wmask, wumask, wvmask) 63 !! are deduced from a product of the two neighboring tmask. 64 !! The vorticity mask (fmask) is deduced from tmask taking 65 !! into account the choice of lateral boundary condition (rn_shlat) : 82 66 !! rn_shlat = 0, free slip (no shear along the coast) 83 67 !! rn_shlat = 2, no slip (specified zero velocity at the coast) … … 85 69 !! 2 < rn_shlat, strong slip | in the lateral boundary layer 86 70 !! 87 !! N.B. If nperio not equal to 0, the land/ocean mask arrays 88 !! are defined with the proper value at lateral domain boundaries. 71 !! tmask_i : interior ocean mask at t-point, i.e. excluding duplicated 72 !! rows/lines due to cyclic or North Fold boundaries as well 73 !! as MPP halos. 74 !! tmask_h : halo mask at t-point, i.e. excluding duplicated rows/lines 75 !! due to cyclic or North Fold boundaries as well 76 !! as MPP halos. 89 77 !! 90 78 !! In case of open boundaries (lk_bdy=T): 91 !! - tmask is set to 1 on the points to be computed b ay the open79 !! - tmask is set to 1 on the points to be computed by the open 92 80 !! boundaries routines. 93 81 !! 94 !! ** Action : tmask : land/ocean mask at t-point (=0. or 1.) 95 !! umask : land/ocean mask at u-point (=0. or 1.) 96 !! vmask : land/ocean mask at v-point (=0. or 1.) 97 !! fmask : land/ocean mask at f-point (=0. or 1.) 98 !! =rn_shlat along lateral boundaries 99 !! tmask_i : interior ocean mask 82 !! ** Action : tmask, umask, vmask, wmask, wumask, wvmask : land/ocean mask 83 !! at t-, u-, v- w, wu-, and wv-points (=0. or 1.) 84 !! fmask : land/ocean mask at f-point (=0., or =1., or 85 !! =rn_shlat along lateral boundaries) 86 !! tmask_i : interior ocean mask 87 !! tmask_h : halo mask 88 !! ssmask , ssumask, ssvmask, ssfmask : 2D ocean mask 100 89 !!---------------------------------------------------------------------- 90 INTEGER, DIMENSION(:,:), INTENT(in) :: k_top, k_bot ! first and last ocean level 91 ! 101 92 INTEGER :: ji, jj, jk ! dummy loop indices 102 93 INTEGER :: iif, iil, ii0, ii1, ii ! local integers 103 94 INTEGER :: ijf, ijl, ij0, ij1 ! - - 95 INTEGER :: iktop, ikbot ! - - 104 96 INTEGER :: ios 105 97 INTEGER :: isrow ! index for ORCA1 starting row 106 INTEGER , POINTER, DIMENSION(:,:) :: imsk107 98 REAL(wp), POINTER, DIMENSION(:,:) :: zwf 108 99 !! … … 111 102 ! 112 103 IF( nn_timing == 1 ) CALL timing_start('dom_msk') 113 !114 CALL wrk_alloc( jpi, jpj, imsk )115 CALL wrk_alloc( jpi, jpj, zwf )116 104 ! 117 105 REWIND( numnam_ref ) ! Namelist namlbc in reference namelist : Lateral momentum boundary condition … … 142 130 ENDIF 143 131 144 ! 1. Ocean/land mask at t-point (computed from mbathy) 145 ! ----------------------------- 146 ! N.B. tmask has already the right boundary conditions since mbathy is ok 147 ! 148 tmask(:,:,:) = 0._wp 149 DO jk = 1, jpk 150 DO jj = 1, jpj 151 DO ji = 1, jpi 152 IF( REAL( mbathy(ji,jj) - jk, wp ) + 0.1_wp >= 0._wp ) tmask(ji,jj,jk) = 1._wp 153 END DO 132 133 ! Ocean/land mask at t-point (computed from ko_top and ko_bot) 134 ! ---------------------------- 135 ! 136 DO jj = 1, jpj 137 DO ji = 1, jpi 138 iktop = k_top(ji,jj) 139 ikbot = k_bot(ji,jj) 140 tmask(ji,jj, 1:iktop-1) = 0._wp ! mask the iceshelves 141 tmask(ji,jj,iktop :ikbot ) = 1._wp 142 tmask(ji,jj,ikbot+1:jpkglo ) = 0._wp ! mask the ocean topography 154 143 END DO 155 144 END DO 145 146 ! 2D ocean mask (=1 if at least one level of the water column is ocean, =0 otherwise) 147 WHERE( k_bot(:,:) > 0 ) ; ssmask(:,:) = 1._wp 148 ELSEWHERE ; ssmask(:,:) = 0._wp 149 END WHERE 156 150 157 ! (ISF) define barotropic mask and mask the ice shelf point158 ssmask(:,:)=tmask(:,:,1) ! at this stage ice shelf is not masked159 151 160 DO jk = 1, jpk 161 DO jj = 1, jpj 162 DO ji = 1, jpi 163 IF( REAL( misfdep(ji,jj) - jk, wp ) - 0.1_wp >= 0._wp ) THEN 164 tmask(ji,jj,jk) = 0._wp 165 END IF 166 END DO 167 END DO 168 END DO 169 170 ! Interior domain mask (used for global sum) 152 ! Interior domain mask (used for global sum) 171 153 ! -------------------- 172 tmask_i(:,:) = ssmask(:,:) ! (ISH) tmask_i = 1 even on the ice shelf 173 174 tmask_h(:,:) = 1._wp ! 0 on the halo and 1 elsewhere 175 iif = jpreci ! ??? 176 iil = nlci - jpreci + 1 177 ijf = jprecj ! ??? 178 ijl = nlcj - jprecj + 1 179 154 ! 155 iif = jpreci ; iil = nlci - jpreci + 1 156 ijf = jprecj ; ijl = nlcj - jprecj + 1 157 ! 158 ! ! halo mask : 0 on the halo and 1 elsewhere 159 tmask_h(:,:) = 1._wp 180 160 tmask_h( 1 :iif, : ) = 0._wp ! first columns 181 161 tmask_h(iil:jpi, : ) = 0._wp ! last columns (including mpp extra columns) 182 162 tmask_h( : , 1 :ijf) = 0._wp ! first rows 183 163 tmask_h( : ,ijl:jpj) = 0._wp ! last rows (including mpp extra rows) 184 185 ! north fold mask 186 ! --------------- 164 ! 165 ! ! north fold mask 187 166 tpol(1:jpiglo) = 1._wp 188 167 fpol(1:jpiglo) = 1._wp … … 190 169 tpol(jpiglo/2+1:jpiglo) = 0._wp 191 170 fpol( 1 :jpiglo) = 0._wp 192 IF( mjg(nlej) == jpjglo ) THEN ! only half of the nlcj-1 row 171 IF( mjg(nlej) == jpjglo ) THEN ! only half of the nlcj-1 row for tmask_h 193 172 DO ji = iif+1, iil-1 194 173 tmask_h(ji,nlej-1) = tmask_h(ji,nlej-1) * tpol(mig(ji)) … … 196 175 ENDIF 197 176 ENDIF 198 199 tmask_i(:,:) = tmask_i(:,:) * tmask_h(:,:) 200 177 ! 201 178 IF( jperio == 5 .OR. jperio == 6 ) THEN ! F-point pivot 202 179 tpol( 1 :jpiglo) = 0._wp 203 180 fpol(jpiglo/2+1:jpiglo) = 0._wp 204 181 ENDIF 205 206 ! 2. Ocean/land mask at u-, v-, and z-points (computed from tmask) 207 ! ------------------------------------------- 182 ! 183 ! ! interior mask : 2D ocean mask x halo mask 184 tmask_i(:,:) = ssmask(:,:) * tmask_h(:,:) 185 186 187 ! Ocean/land mask at u-, v-, and z-points (computed from tmask) 188 ! ---------------------------------------- 208 189 DO jk = 1, jpk 209 190 DO jj = 1, jpjm1 … … 221 202 DO jj = 1, jpjm1 222 203 DO ji = 1, fs_jpim1 ! vector loop 204 !!gm simpler : 205 ! ssumask(ji,jj) = MIN( 1._wp , SUM( umask(ji,jj,:) ) ) 206 ! ssvmask(ji,jj) = MIN( 1._wp , SUM( vmask(ji,jj,:) ) ) 207 !!gm 208 !!gm faster : 209 ! ssumask(ji,jj) = ssmask(ji,jj) * tmask(ji+1,jj ) 210 ! ssvmask(ji,jj) = ssmask(ji,jj) * tmask(ji ,jj+1) 211 !!gm 223 212 ssumask(ji,jj) = ssmask(ji,jj) * ssmask(ji+1,jj ) * MIN(1._wp,SUM(umask(ji,jj,:))) 224 213 ssvmask(ji,jj) = ssmask(ji,jj) * ssmask(ji ,jj+1) * MIN(1._wp,SUM(vmask(ji,jj,:))) 214 !!end 225 215 END DO 226 216 DO ji = 1, jpim1 ! NO vector opt. 217 !!gm faster 218 ! ssfmask(ji,jj) = ssmask(ji,jj ) * ssmask(ji+1,jj ) & 219 ! & * ssmask(ji,jj+1) * ssmask(ji+1,jj+1) 220 !!gm 227 221 ssfmask(ji,jj) = ssmask(ji,jj ) * ssmask(ji+1,jj ) & 228 222 & * ssmask(ji,jj+1) * ssmask(ji+1,jj+1) * MIN(1._wp,SUM(fmask(ji,jj,:))) 223 !!gm 229 224 END DO 230 225 END DO 231 226 CALL lbc_lnk( umask , 'U', 1._wp ) ! Lateral boundary conditions 232 227 CALL lbc_lnk( vmask , 'V', 1._wp ) 233 CALL lbc_lnk( fmask , 'F', 1._wp ) 234 CALL lbc_lnk( ssumask, 'U', 1._wp ) ! Lateral boundary conditions228 ! CALL lbc_lnk( fmask , 'F', 1._wp ) ! applied after the specification of lateral b.c. 229 CALL lbc_lnk( ssumask, 'U', 1._wp ) 235 230 CALL lbc_lnk( ssvmask, 'V', 1._wp ) 236 231 CALL lbc_lnk( ssfmask, 'F', 1._wp ) 237 232 238 ! 3. Ocean/land mask at wu-, wv- and w points 233 234 ! Ocean/land mask at wu-, wv- and w points 239 235 !---------------------------------------------- 240 236 wmask (:,:,1) = tmask(:,:,1) ! surface … … 247 243 END DO 248 244 245 249 246 ! Lateral boundary conditions on velocity (modify fmask) 250 247 ! --------------------------------------- 248 CALL wrk_alloc( jpi,jpj, zwf ) 249 ! 251 250 DO jk = 1, jpk 252 251 zwf(:,:) = fmask(:,:,jk) … … 277 276 END DO 278 277 ! 278 CALL wrk_dealloc( jpi,jpj, zwf ) 279 ! 279 280 IF( cp_cfg == "orca" .AND. jp_cfg == 2 ) THEN ! ORCA_R2 configuration 280 281 ! ! Increased lateral friction near of some straits … … 349 350 ! CAUTION : The fmask may be further modified in dyn_vor_init ( dynvor.F90 ) 350 351 ! 351 CALL wrk_dealloc( jpi, jpj, imsk )352 CALL wrk_dealloc( jpi, jpj, zwf )353 !354 352 IF( nn_timing == 1 ) CALL timing_stop('dom_msk') 355 353 !
Note: See TracChangeset
for help on using the changeset viewer.