Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Normal
Revision Log

dommsk.F90 @ 6808

Last change on this file since 6808 was 6808, checked in by jamesharle, 8 years ago
merge with trunk@6232 for consistency with SSB code
Property svn:keywords set to `Id`
File size: 17.4 KB

Rev	Line
[3]	1	MODULE dommsk
[1566]	2	!!======================================================================
[3]	3	!! * MODULE dommsk *
	4	!! Ocean initialization : domain land/sea mask
[1566]	5	!!======================================================================
	6	!! History : OPA ! 1987-07 (G. Madec) Original code
[2528]	7	!! 6.0 ! 1993-03 (M. Guyon) symetrical conditions (M. Guyon)
	8	!! 7.0 ! 1996-01 (G. Madec) suppression of common work arrays
[6808]	9	!! - ! 1996-05 (G. Madec) mask computed from tmask
[2528]	10	!! 8.0 ! 1997-02 (G. Madec) mesh information put in domhgr.F
	11	!! 8.1 ! 1997-07 (G. Madec) modification of mbathy and fmask
[1566]	12	!! - ! 1998-05 (G. Roullet) free surface
[2528]	13	!! 8.2 ! 2000-03 (G. Madec) no slip accurate
[1566]	14	!! - ! 2001-09 (J.-M. Molines) Open boundaries
	15	!! NEMO 1.0 ! 2002-08 (G. Madec) F90: Free form and module
	16	!! - ! 2005-11 (V. Garnier) Surface pressure gradient organization
	17	!! 3.2 ! 2009-07 (R. Benshila) Suppression of rigid-lid option
[6808]	18	!! 3.6 ! 2015-05 (P. Mathiot) ISF: add wmask,wumask and wvmask
[1566]	19	!!----------------------------------------------------------------------
[3]	20
	21	!!----------------------------------------------------------------------
	22	!! dom_msk : compute land/ocean mask
	23	!!----------------------------------------------------------------------
	24	USE oce ! ocean dynamics and tracers
	25	USE dom_oce ! ocean space and time domain
[6808]	26	!
[3]	27	USE in_out_manager ! I/O manager
	28	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
[6808]	29	USE lib_mpp !
[3294]	30	USE wrk_nemo ! Memory allocation
	31	USE timing ! Timing
[3]	32
	33	IMPLICIT NONE
	34	PRIVATE
	35
[6808]	36	PUBLIC dom_msk ! routine called by inidom.F90
[3]	37
[1601]	38	! !!* Namelist namlbc : lateral boundary condition *
[4147]	39	REAL(wp) :: rn_shlat ! type of lateral boundary condition on velocity
	40	LOGICAL, PUBLIC :: ln_vorlat ! consistency of vorticity boundary condition
[3294]	41	! with analytical eqs.
[2715]	42
[3]	43	!! * Substitutions
	44	# include "vectopt_loop_substitute.h90"
[1566]	45	!!----------------------------------------------------------------------
	46	!! NEMO/OPA 3.2 , LODYC-IPSL (2009)
[1152]	47	!! $Id$
[2528]	48	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
[1566]	49	!!----------------------------------------------------------------------
[3]	50	CONTAINS
[2715]	51
[3]	52	SUBROUTINE dom_msk
	53	!!---------------------------------------------------------------------
	54	!! * ROUTINE dom_msk *
	55	!!
	56	!! ** Purpose : Compute land/ocean mask arrays at tracer points, hori-
[6808]	57	!! zontal velocity points (u & v), vorticity points (f) points.
[3]	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.
[1528]	61	!! mbathy equals 0 over continental T-point
	62	!! and the number of ocean level over the ocean.
[3]	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)
[2528]	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.
[3]	78	!!
	79	!! The lateral friction is set through the value of fmask along
[1601]	80	!! the coast and topography. This value is defined by rn_shlat, a
[3]	81	!! namelist parameter:
[1601]	82	!! rn_shlat = 0, free slip (no shear along the coast)
	83	!! rn_shlat = 2, no slip (specified zero velocity at the coast)
	84	!! 0 < rn_shlat < 2, partial slip \| non-linear velocity profile
	85	!! 2 < rn_shlat, strong slip \| in the lateral boundary layer
[3]	86	!!
	87	!! N.B. If nperio not equal to 0, the land/ocean mask arrays
[6808]	88	!! are defined with the proper value at lateral domain boundaries.
[3]	89	!!
[4328]	90	!! In case of open boundaries (lk_bdy=T):
[3]	91	!! - tmask is set to 1 on the points to be computed bay the open
	92	!! boundaries routines.
	93	!!
[1566]	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.)
[1601]	98	!! =rn_shlat along lateral boundaries
[2528]	99	!! tmask_i : interior ocean mask
[3]	100	!!----------------------------------------------------------------------
[6808]	101	INTEGER :: ji, jj, jk ! dummy loop indices
[2715]	102	INTEGER :: iif, iil, ii0, ii1, ii ! local integers
	103	INTEGER :: ijf, ijl, ij0, ij1 ! - -
[4147]	104	INTEGER :: ios
[5385]	105	INTEGER :: isrow ! index for ORCA1 starting row
[3294]	106	INTEGER , POINTER, DIMENSION(:,:) :: imsk
	107	REAL(wp), POINTER, DIMENSION(:,:) :: zwf
[1601]	108	!!
[3294]	109	NAMELIST/namlbc/ rn_shlat, ln_vorlat
[3]	110	!!---------------------------------------------------------------------
[3294]	111	!
	112	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	!
[4147]	117	REWIND( numnam_ref ) ! Namelist namlbc in reference namelist : Lateral momentum boundary condition
	118	READ ( numnam_ref, namlbc, IOSTAT = ios, ERR = 901 )
	119	901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namlbc in reference namelist', lwp )
	120
	121	REWIND( numnam_cfg ) ! Namelist namlbc in configuration namelist : Lateral momentum boundary condition
	122	READ ( numnam_cfg, namlbc, IOSTAT = ios, ERR = 902 )
	123	902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namlbc in configuration namelist', lwp )
[4624]	124	IF(lwm) WRITE ( numond, namlbc )
[1566]	125
	126	IF(lwp) THEN ! control print
[3]	127	WRITE(numout,*)
	128	WRITE(numout,*) 'dommsk : ocean mask '
	129	WRITE(numout,*) '~~~~~~'
[1566]	130	WRITE(numout,*) ' Namelist namlbc'
[3294]	131	WRITE(numout,*) ' lateral momentum boundary cond. rn_shlat = ',rn_shlat
	132	WRITE(numout,*) ' consistency with analytical form ln_vorlat = ',ln_vorlat
[3]	133	ENDIF
	134
[2528]	135	IF ( rn_shlat == 0. ) THEN ; IF(lwp) WRITE(numout,*) ' ocean lateral free-slip '
[1601]	136	ELSEIF ( rn_shlat == 2. ) THEN ; IF(lwp) WRITE(numout,*) ' ocean lateral no-slip '
	137	ELSEIF ( 0. < rn_shlat .AND. rn_shlat < 2. ) THEN ; IF(lwp) WRITE(numout,*) ' ocean lateral partial-slip '
	138	ELSEIF ( 2. < rn_shlat ) THEN ; IF(lwp) WRITE(numout,*) ' ocean lateral strong-slip '
	139	ELSE
	140	WRITE(ctmp1,*) ' rn_shlat is negative = ', rn_shlat
	141	CALL ctl_stop( ctmp1 )
[3]	142	ENDIF
	143
	144	! 1. Ocean/land mask at t-point (computed from mbathy)
	145	! -----------------------------
[1566]	146	! N.B. tmask has already the right boundary conditions since mbathy is ok
	147	!
[2528]	148	tmask(:,:,:) = 0._wp
[3]	149	DO jk = 1, jpk
	150	DO jj = 1, jpj
	151	DO ji = 1, jpi
[2528]	152	IF( REAL( mbathy(ji,jj) - jk, wp ) + 0.1_wp >= 0._wp ) tmask(ji,jj,jk) = 1._wp
[3]	153	END DO
	154	END DO
	155	END DO
[4990]	156
	157	! (ISF) define barotropic mask and mask the ice shelf point
	158	ssmask(:,:)=tmask(:,:,1) ! at this stage ice shelf is not masked
	159
	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
[3]	169
	170	! Interior domain mask (used for global sum)
	171	! --------------------
[4990]	172	tmask_i(:,:) = ssmask(:,:) ! (ISH) tmask_i = 1 even on the ice shelf
[6808]	173
	174	tmask_h(:,:) = 1._wp ! 0 on the halo and 1 elsewhere
[3]	175	iif = jpreci ! ???
	176	iil = nlci - jpreci + 1
	177	ijf = jprecj ! ???
	178	ijl = nlcj - jprecj + 1
	179
[6808]	180	tmask_h( 1 :iif, : ) = 0._wp ! first columns
	181	tmask_h(iil:jpi, : ) = 0._wp ! last columns (including mpp extra columns)
	182	tmask_h( : , 1 :ijf) = 0._wp ! first rows
	183	tmask_h( : ,ijl:jpj) = 0._wp ! last rows (including mpp extra rows)
[3]	184
	185	! north fold mask
[1566]	186	! ---------------
[2528]	187	tpol(1:jpiglo) = 1._wp
	188	fpol(1:jpiglo) = 1._wp
[3]	189	IF( jperio == 3 .OR. jperio == 4 ) THEN ! T-point pivot
[2528]	190	tpol(jpiglo/2+1:jpiglo) = 0._wp
	191	fpol( 1 :jpiglo) = 0._wp
[1566]	192	IF( mjg(nlej) == jpjglo ) THEN ! only half of the nlcj-1 row
[291]	193	DO ji = iif+1, iil-1
[6808]	194	tmask_h(ji,nlej-1) = tmask_h(ji,nlej-1) * tpol(mig(ji))
[291]	195	END DO
	196	ENDIF
[3]	197	ENDIF
[6808]	198
	199	tmask_i(:,:) = tmask_i(:,:) * tmask_h(:,:)
	200
[3]	201	IF( jperio == 5 .OR. jperio == 6 ) THEN ! F-point pivot
[2528]	202	tpol( 1 :jpiglo) = 0._wp
	203	fpol(jpiglo/2+1:jpiglo) = 0._wp
[3]	204	ENDIF
	205
	206	! 2. Ocean/land mask at u-, v-, and z-points (computed from tmask)
	207	! -------------------------------------------
	208	DO jk = 1, jpk
	209	DO jj = 1, jpjm1
	210	DO ji = 1, fs_jpim1 ! vector loop
	211	umask(ji,jj,jk) = tmask(ji,jj ,jk) * tmask(ji+1,jj ,jk)
	212	vmask(ji,jj,jk) = tmask(ji,jj ,jk) * tmask(ji ,jj+1,jk)
[1271]	213	END DO
[1694]	214	DO ji = 1, jpim1 ! NO vector opt.
[3]	215	fmask(ji,jj,jk) = tmask(ji,jj ,jk) * tmask(ji+1,jj ,jk) &
[62]	216	& * tmask(ji,jj+1,jk) * tmask(ji+1,jj+1,jk)
[3]	217	END DO
	218	END DO
	219	END DO
[6808]	220	! (ISF) MIN(1,SUM(umask)) is here to check if you have effectively at least 1 wet cell at u point
[4990]	221	DO jj = 1, jpjm1
	222	DO ji = 1, fs_jpim1 ! vector loop
[6808]	223	ssumask(ji,jj) = ssmask(ji,jj) * ssmask(ji+1,jj ) * MIN(1._wp,SUM(umask(ji,jj,:)))
	224	ssvmask(ji,jj) = ssmask(ji,jj) * ssmask(ji ,jj+1) * MIN(1._wp,SUM(vmask(ji,jj,:)))
[4990]	225	END DO
	226	DO ji = 1, jpim1 ! NO vector opt.
[6808]	227	ssfmask(ji,jj) = ssmask(ji,jj ) * ssmask(ji+1,jj ) &
[4990]	228	& * ssmask(ji,jj+1) * ssmask(ji+1,jj+1) * MIN(1._wp,SUM(fmask(ji,jj,:)))
	229	END DO
	230	END DO
[6808]	231	CALL lbc_lnk( umask , 'U', 1._wp ) ! Lateral boundary conditions
	232	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 conditions
	235	CALL lbc_lnk( ssvmask, 'V', 1._wp )
	236	CALL lbc_lnk( ssfmask, 'F', 1._wp )
[3]	237
[5120]	238	! 3. Ocean/land mask at wu-, wv- and w points
	239	!----------------------------------------------
[6808]	240	wmask (:,:,1) = tmask(:,:,1) ! surface
	241	wumask(:,:,1) = umask(:,:,1)
	242	wvmask(:,:,1) = vmask(:,:,1)
	243	DO jk = 2, jpk ! interior values
	244	wmask (:,:,jk) = tmask(:,:,jk) * tmask(:,:,jk-1)
	245	wumask(:,:,jk) = umask(:,:,jk) * umask(:,:,jk-1)
	246	wvmask(:,:,jk) = vmask(:,:,jk) * vmask(:,:,jk-1)
[5120]	247	END DO
[3]	248
	249	! Lateral boundary conditions on velocity (modify fmask)
[1566]	250	! ---------------------------------------
[3]	251	DO jk = 1, jpk
[1566]	252	zwf(:,:) = fmask(:,:,jk)
[3]	253	DO jj = 2, jpjm1
	254	DO ji = fs_2, fs_jpim1 ! vector opt.
[2528]	255	IF( fmask(ji,jj,jk) == 0._wp ) THEN
	256	fmask(ji,jj,jk) = rn_shlat * MIN( 1._wp , MAX( zwf(ji+1,jj), zwf(ji,jj+1), &
	257	& zwf(ji-1,jj), zwf(ji,jj-1) ) )
[3]	258	ENDIF
	259	END DO
	260	END DO
	261	DO jj = 2, jpjm1
[2528]	262	IF( fmask(1,jj,jk) == 0._wp ) THEN
	263	fmask(1 ,jj,jk) = rn_shlat * MIN( 1._wp , MAX( zwf(2,jj), zwf(1,jj+1), zwf(1,jj-1) ) )
[3]	264	ENDIF
[2528]	265	IF( fmask(jpi,jj,jk) == 0._wp ) THEN
	266	fmask(jpi,jj,jk) = rn_shlat * MIN( 1._wp , MAX( zwf(jpi,jj+1), zwf(jpim1,jj), zwf(jpi,jj-1) ) )
[3]	267	ENDIF
[1566]	268	END DO
[3]	269	DO ji = 2, jpim1
[2528]	270	IF( fmask(ji,1,jk) == 0._wp ) THEN
	271	fmask(ji, 1 ,jk) = rn_shlat * MIN( 1._wp , MAX( zwf(ji+1,1), zwf(ji,2), zwf(ji-1,1) ) )
[3]	272	ENDIF
[2528]	273	IF( fmask(ji,jpj,jk) == 0._wp ) THEN
	274	fmask(ji,jpj,jk) = rn_shlat * MIN( 1._wp , MAX( zwf(ji+1,jpj), zwf(ji-1,jpj), zwf(ji,jpjm1) ) )
[3]	275	ENDIF
	276	END DO
	277	END DO
[1566]	278	!
	279	IF( cp_cfg == "orca" .AND. jp_cfg == 2 ) THEN ! ORCA_R2 configuration
	280	! ! Increased lateral friction near of some straits
[6808]	281	! ! Gibraltar strait : partial slip (fmask=0.5)
	282	ij0 = 101 ; ij1 = 101
	283	ii0 = 139 ; ii1 = 140 ; fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) = 0.5_wp
	284	ij0 = 102 ; ij1 = 102
	285	ii0 = 139 ; ii1 = 140 ; fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) = 0.5_wp
	286	!
	287	! ! Bab el Mandeb : partial slip (fmask=1)
	288	ij0 = 87 ; ij1 = 88
	289	ii0 = 160 ; ii1 = 160 ; fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) = 1._wp
	290	ij0 = 88 ; ij1 = 88
	291	ii0 = 159 ; ii1 = 159 ; fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) = 1._wp
	292	!
[1707]	293	! ! Danish straits : strong slip (fmask > 2)
	294	! We keep this as an example but it is instable in this case
	295	! ij0 = 115 ; ij1 = 115
[2528]	296	! ii0 = 145 ; ii1 = 146 ; fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) = 4._wp
[1707]	297	! ij0 = 116 ; ij1 = 116
[2528]	298	! ii0 = 145 ; ii1 = 146 ; fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) = 4._wp
[3]	299	!
	300	ENDIF
[1566]	301	!
[2528]	302	IF( cp_cfg == "orca" .AND. jp_cfg == 1 ) THEN ! ORCA R1 configuration
	303	! ! Increased lateral friction near of some straits
[5506]	304	! This dirty section will be suppressed by simplification process:
	305	! all this will come back in input files
	306	! Currently these hard-wired indices relate to configuration with
	307	! extend grid (jpjglo=332)
	308	!
	309	isrow = 332 - jpjglo
	310	!
[2528]	311	IF(lwp) WRITE(numout,*)
	312	IF(lwp) WRITE(numout,*) ' orca_r1: increase friction near the following straits : '
	313	IF(lwp) WRITE(numout,*) ' Gibraltar '
[5385]	314	ii0 = 282 ; ii1 = 283 ! Gibraltar Strait
[6808]	315	ij0 = 241 - isrow ; ij1 = 241 - isrow ; fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp
[1566]	316
[2528]	317	IF(lwp) WRITE(numout,*) ' Bhosporus '
[5385]	318	ii0 = 314 ; ii1 = 315 ! Bhosporus Strait
[6808]	319	ij0 = 248 - isrow ; ij1 = 248 - isrow ; fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp
[2528]	320
	321	IF(lwp) WRITE(numout,*) ' Makassar (Top) '
[5385]	322	ii0 = 48 ; ii1 = 48 ! Makassar Strait (Top)
[6808]	323	ij0 = 189 - isrow ; ij1 = 190 - isrow ; fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 3._wp
[2528]	324
	325	IF(lwp) WRITE(numout,*) ' Lombok '
[5385]	326	ii0 = 44 ; ii1 = 44 ! Lombok Strait
[6808]	327	ij0 = 164 - isrow ; ij1 = 165 - isrow ; fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp
[2528]	328
	329	IF(lwp) WRITE(numout,*) ' Ombai '
[5385]	330	ii0 = 53 ; ii1 = 53 ! Ombai Strait
[6808]	331	ij0 = 164 - isrow ; ij1 = 165 - isrow ; fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp
[2528]	332
	333	IF(lwp) WRITE(numout,*) ' Timor Passage '
[5385]	334	ii0 = 56 ; ii1 = 56 ! Timor Passage
[6808]	335	ij0 = 164 - isrow ; ij1 = 165 - isrow ; fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp
[2528]	336
	337	IF(lwp) WRITE(numout,*) ' West Halmahera '
[5385]	338	ii0 = 58 ; ii1 = 58 ! West Halmahera Strait
[6808]	339	ij0 = 181 - isrow ; ij1 = 182 - isrow ; fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 3._wp
[2528]	340
	341	IF(lwp) WRITE(numout,*) ' East Halmahera '
[5385]	342	ii0 = 55 ; ii1 = 55 ! East Halmahera Strait
[6808]	343	ij0 = 181 - isrow ; ij1 = 182 - isrow ; fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 3._wp
[2528]	344	!
	345	ENDIF
	346	!
	347	CALL lbc_lnk( fmask, 'F', 1._wp ) ! Lateral boundary conditions on fmask
[6808]	348	!
[3294]	349	! CAUTION : The fmask may be further modified in dyn_vor_init ( dynvor.F90 )
[1566]	350	!
[3294]	351	CALL wrk_dealloc( jpi, jpj, imsk )
	352	CALL wrk_dealloc( jpi, jpj, zwf )
[2715]	353	!
[3294]	354	IF( nn_timing == 1 ) CALL timing_stop('dom_msk')
	355	!
[3]	356	END SUBROUTINE dom_msk
	357
	358	!!======================================================================
	359	END MODULE dommsk

Note: See TracBrowser for help on using the repository browser.

New URL for NEMO forge! http://forge.nemo-ocean.eu

Context Navigation

source: branches/NERC/dev_r5549_BDY_ZEROGRAD/NEMOGCM/NEMO/OPA_SRC/DOM/dommsk.F90 @ 6808

Download in other formats: