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.

dommsk.F90 in branches/UKMO/AMM15_v3_6_STABLE_package/NEMOGCM/NEMO/OPA_SRC/DOM – NEMO

Context Navigation

source: branches/UKMO/AMM15_v3_6_STABLE_package/NEMOGCM/NEMO/OPA_SRC/DOM/dommsk.F90 @ 10064

Last change on this file since 10064 was 8059, checked in by jgraham, 7 years ago
Merged branches required for AMM15 simulations, see ticket #1904. Merged branches include: branches/UKMO/CO6_KD490 branches/UKMO/CO6_Restartable_Tidal_Analysis branches/UKMO/AMM15_v3_6_STABLE
File size: 30.7 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
[1566]	9	!! - ! 1996-05 (G. Madec) mask computed from tmask and sup-
	10	!! ! pression of the double computation of bmask
[2528]	11	!! 8.0 ! 1997-02 (G. Madec) mesh information put in domhgr.F
	12	!! 8.1 ! 1997-07 (G. Madec) modification of mbathy and fmask
[1566]	13	!! - ! 1998-05 (G. Roullet) free surface
[2528]	14	!! 8.2 ! 2000-03 (G. Madec) no slip accurate
[1566]	15	!! - ! 2001-09 (J.-M. Molines) Open boundaries
	16	!! NEMO 1.0 ! 2002-08 (G. Madec) F90: Free form and module
	17	!! - ! 2005-11 (V. Garnier) Surface pressure gradient organization
	18	!! 3.2 ! 2009-07 (R. Benshila) Suppression of rigid-lid option
	19	!!----------------------------------------------------------------------
[3]	20
	21	!!----------------------------------------------------------------------
	22	!! dom_msk : compute land/ocean mask
[1566]	23	!! dom_msk_nsa : update land/ocean mask when no-slip accurate option is used.
[3]	24	!!----------------------------------------------------------------------
	25	USE oce ! ocean dynamics and tracers
	26	USE dom_oce ! ocean space and time domain
	27	USE in_out_manager ! I/O manager
	28	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
[32]	29	USE lib_mpp
[367]	30	USE dynspg_oce ! choice/control of key cpp for surface pressure gradient
[3294]	31	USE wrk_nemo ! Memory allocation
	32	USE timing ! Timing
[8059]	33	USE iom ! slwa
[3]	34
	35	IMPLICIT NONE
	36	PRIVATE
	37
[2715]	38	PUBLIC dom_msk ! routine called by inidom.F90
	39	PUBLIC dom_msk_alloc ! routine called by nemogcm.F90
[3]	40
[1601]	41	! !!* Namelist namlbc : lateral boundary condition *
[4147]	42	REAL(wp) :: rn_shlat ! type of lateral boundary condition on velocity
	43	LOGICAL, PUBLIC :: ln_vorlat ! consistency of vorticity boundary condition
[3294]	44	! with analytical eqs.
[2715]	45
[3294]	46
[2715]	47	INTEGER, ALLOCATABLE, SAVE, DIMENSION(:,:) :: icoord ! Workspace for dom_msk_nsa()
	48
[3]	49	!! * Substitutions
	50	# include "vectopt_loop_substitute.h90"
[1566]	51	!!----------------------------------------------------------------------
	52	!! NEMO/OPA 3.2 , LODYC-IPSL (2009)
[8058]	53	!! $Id$
[2528]	54	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
[1566]	55	!!----------------------------------------------------------------------
[3]	56	CONTAINS
	57
[2715]	58	INTEGER FUNCTION dom_msk_alloc()
	59	!!---------------------------------------------------------------------
	60	!! * FUNCTION dom_msk_alloc *
	61	!!---------------------------------------------------------------------
	62	dom_msk_alloc = 0
	63	#if defined key_noslip_accurate
	64	ALLOCATE(icoord(jpijpjjpk,3), STAT=dom_msk_alloc)
	65	#endif
	66	IF( dom_msk_alloc /= 0 ) CALL ctl_warn('dom_msk_alloc: failed to allocate icoord array')
	67	!
	68	END FUNCTION dom_msk_alloc
	69
	70
[3]	71	SUBROUTINE dom_msk
	72	!!---------------------------------------------------------------------
	73	!! * ROUTINE dom_msk *
	74	!!
	75	!! ** Purpose : Compute land/ocean mask arrays at tracer points, hori-
	76	!! zontal velocity points (u & v), vorticity points (f) and baro-
	77	!! tropic stream function points (b).
	78	!!
	79	!! ** Method : The ocean/land mask is computed from the basin bathy-
	80	!! metry in level (mbathy) which is defined or read in dommba.
[1528]	81	!! mbathy equals 0 over continental T-point
	82	!! and the number of ocean level over the ocean.
[3]	83	!!
	84	!! At a given position (ji,jj,jk) the ocean/land mask is given by:
	85	!! t-point : 0. IF mbathy( ji ,jj) =< 0
	86	!! 1. IF mbathy( ji ,jj) >= jk
	87	!! u-point : 0. IF mbathy( ji ,jj) or mbathy(ji+1, jj ) =< 0
	88	!! 1. IF mbathy( ji ,jj) and mbathy(ji+1, jj ) >= jk.
	89	!! v-point : 0. IF mbathy( ji ,jj) or mbathy( ji ,jj+1) =< 0
	90	!! 1. IF mbathy( ji ,jj) and mbathy( ji ,jj+1) >= jk.
	91	!! f-point : 0. IF mbathy( ji ,jj) or mbathy( ji ,jj+1)
	92	!! or mbathy(ji+1,jj) or mbathy(ji+1,jj+1) =< 0
	93	!! 1. IF mbathy( ji ,jj) and mbathy( ji ,jj+1)
[2528]	94	!! and mbathy(ji+1,jj) and mbathy(ji+1,jj+1) >= jk.
[3]	95	!! b-point : the same definition as for f-point of the first ocean
	96	!! level (surface level) but with 0 along coastlines.
[2528]	97	!! tmask_i : interior ocean mask at t-point, i.e. excluding duplicated
	98	!! rows/lines due to cyclic or North Fold boundaries as well
	99	!! as MPP halos.
[3]	100	!!
	101	!! The lateral friction is set through the value of fmask along
[1601]	102	!! the coast and topography. This value is defined by rn_shlat, a
[3]	103	!! namelist parameter:
[1601]	104	!! rn_shlat = 0, free slip (no shear along the coast)
	105	!! rn_shlat = 2, no slip (specified zero velocity at the coast)
	106	!! 0 < rn_shlat < 2, partial slip \| non-linear velocity profile
	107	!! 2 < rn_shlat, strong slip \| in the lateral boundary layer
[3]	108	!!
	109	!! N.B. If nperio not equal to 0, the land/ocean mask arrays
	110	!! are defined with the proper value at lateral domain boundaries,
	111	!! but bmask. indeed, bmask defined the domain over which the
	112	!! barotropic stream function is computed. this domain cannot
	113	!! contain identical columns because the matrix associated with
	114	!! the barotropic stream function equation is then no more inverti-
	115	!! ble. therefore bmask is set to 0 along lateral domain boundaries
	116	!! even IF nperio is not zero.
	117	!!
[4328]	118	!! In case of open boundaries (lk_bdy=T):
[3]	119	!! - tmask is set to 1 on the points to be computed bay the open
	120	!! boundaries routines.
	121	!! - bmask is set to 0 on the open boundaries.
	122	!!
[1566]	123	!! ** Action : tmask : land/ocean mask at t-point (=0. or 1.)
	124	!! umask : land/ocean mask at u-point (=0. or 1.)
	125	!! vmask : land/ocean mask at v-point (=0. or 1.)
	126	!! fmask : land/ocean mask at f-point (=0. or 1.)
[1601]	127	!! =rn_shlat along lateral boundaries
[1566]	128	!! bmask : land/ocean mask at barotropic stream
	129	!! function point (=0. or 1.) and set to 0 along lateral boundaries
[2528]	130	!! tmask_i : interior ocean mask
[3]	131	!!----------------------------------------------------------------------
[2715]	132	!
[454]	133	INTEGER :: ji, jj, jk ! dummy loop indices
[2715]	134	INTEGER :: iif, iil, ii0, ii1, ii ! local integers
	135	INTEGER :: ijf, ijl, ij0, ij1 ! - -
[4147]	136	INTEGER :: ios
[5385]	137	INTEGER :: isrow ! index for ORCA1 starting row
[3294]	138	INTEGER , POINTER, DIMENSION(:,:) :: imsk
	139	REAL(wp), POINTER, DIMENSION(:,:) :: zwf
[1601]	140	!!
[3294]	141	NAMELIST/namlbc/ rn_shlat, ln_vorlat
[3]	142	!!---------------------------------------------------------------------
[3294]	143	!
	144	IF( nn_timing == 1 ) CALL timing_start('dom_msk')
	145	!
	146	CALL wrk_alloc( jpi, jpj, imsk )
	147	CALL wrk_alloc( jpi, jpj, zwf )
	148	!
[4147]	149	REWIND( numnam_ref ) ! Namelist namlbc in reference namelist : Lateral momentum boundary condition
	150	READ ( numnam_ref, namlbc, IOSTAT = ios, ERR = 901 )
	151	901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namlbc in reference namelist', lwp )
	152
	153	REWIND( numnam_cfg ) ! Namelist namlbc in configuration namelist : Lateral momentum boundary condition
	154	READ ( numnam_cfg, namlbc, IOSTAT = ios, ERR = 902 )
	155	902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namlbc in configuration namelist', lwp )
[4624]	156	IF(lwm) WRITE ( numond, namlbc )
[1566]	157
	158	IF(lwp) THEN ! control print
[3]	159	WRITE(numout,*)
	160	WRITE(numout,*) 'dommsk : ocean mask '
	161	WRITE(numout,*) '~~~~~~'
[1566]	162	WRITE(numout,*) ' Namelist namlbc'
[3294]	163	WRITE(numout,*) ' lateral momentum boundary cond. rn_shlat = ',rn_shlat
	164	WRITE(numout,*) ' consistency with analytical form ln_vorlat = ',ln_vorlat
[3]	165	ENDIF
	166
[2528]	167	IF ( rn_shlat == 0. ) THEN ; IF(lwp) WRITE(numout,*) ' ocean lateral free-slip '
[1601]	168	ELSEIF ( rn_shlat == 2. ) THEN ; IF(lwp) WRITE(numout,*) ' ocean lateral no-slip '
	169	ELSEIF ( 0. < rn_shlat .AND. rn_shlat < 2. ) THEN ; IF(lwp) WRITE(numout,*) ' ocean lateral partial-slip '
	170	ELSEIF ( 2. < rn_shlat ) THEN ; IF(lwp) WRITE(numout,*) ' ocean lateral strong-slip '
	171	ELSE
	172	WRITE(ctmp1,*) ' rn_shlat is negative = ', rn_shlat
	173	CALL ctl_stop( ctmp1 )
[3]	174	ENDIF
	175
	176	! 1. Ocean/land mask at t-point (computed from mbathy)
	177	! -----------------------------
[1566]	178	! N.B. tmask has already the right boundary conditions since mbathy is ok
	179	!
[2528]	180	tmask(:,:,:) = 0._wp
[3]	181	DO jk = 1, jpk
	182	DO jj = 1, jpj
	183	DO ji = 1, jpi
[2528]	184	IF( REAL( mbathy(ji,jj) - jk, wp ) + 0.1_wp >= 0._wp ) tmask(ji,jj,jk) = 1._wp
[3]	185	END DO
	186	END DO
	187	END DO
[4990]	188
	189	! (ISF) define barotropic mask and mask the ice shelf point
	190	ssmask(:,:)=tmask(:,:,1) ! at this stage ice shelf is not masked
	191
	192	DO jk = 1, jpk
	193	DO jj = 1, jpj
	194	DO ji = 1, jpi
	195	IF( REAL( misfdep(ji,jj) - jk, wp ) - 0.1_wp >= 0._wp ) THEN
	196	tmask(ji,jj,jk) = 0._wp
	197	END IF
	198	END DO
	199	END DO
	200	END DO
[3]	201
[1566]	202	!!gm ????
[255]	203	#if defined key_zdfkpp
[1601]	204	IF( cp_cfg == 'orca' ) THEN
	205	IF( jp_cfg == 2 ) THEN ! land point on Bab el Mandeb zonal section
[291]	206	ij0 = 87 ; ij1 = 88
	207	ii0 = 160 ; ii1 = 161
[2528]	208	tmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) = 0._wp
[291]	209	ELSE
	210	IF(lwp) WRITE(numout,*)
	211	IF(lwp) WRITE(numout,cform_war)
	212	IF(lwp) WRITE(numout,*)
	213	IF(lwp) WRITE(numout,*)' A mask must be applied on Bab el Mandeb strait'
	214	IF(lwp) WRITE(numout,*)' in case of ORCAs configurations'
	215	IF(lwp) WRITE(numout,*)' This is a problem which is not yet solved'
	216	IF(lwp) WRITE(numout,*)
	217	ENDIF
	218	ENDIF
[255]	219	#endif
[1566]	220	!!gm end
[291]	221
[3]	222	! Interior domain mask (used for global sum)
	223	! --------------------
[4990]	224	tmask_i(:,:) = ssmask(:,:) ! (ISH) tmask_i = 1 even on the ice shelf
[3]	225	iif = jpreci ! ???
	226	iil = nlci - jpreci + 1
	227	ijf = jprecj ! ???
	228	ijl = nlcj - jprecj + 1
	229
[2528]	230	tmask_i( 1 :iif, : ) = 0._wp ! first columns
	231	tmask_i(iil:jpi, : ) = 0._wp ! last columns (including mpp extra columns)
	232	tmask_i( : , 1 :ijf) = 0._wp ! first rows
	233	tmask_i( : ,ijl:jpj) = 0._wp ! last rows (including mpp extra rows)
[3]	234
	235	! north fold mask
[1566]	236	! ---------------
[2528]	237	tpol(1:jpiglo) = 1._wp
	238	fpol(1:jpiglo) = 1._wp
[3]	239	IF( jperio == 3 .OR. jperio == 4 ) THEN ! T-point pivot
[2528]	240	tpol(jpiglo/2+1:jpiglo) = 0._wp
	241	fpol( 1 :jpiglo) = 0._wp
[1566]	242	IF( mjg(nlej) == jpjglo ) THEN ! only half of the nlcj-1 row
[291]	243	DO ji = iif+1, iil-1
	244	tmask_i(ji,nlej-1) = tmask_i(ji,nlej-1) * tpol(mig(ji))
	245	END DO
	246	ENDIF
[3]	247	ENDIF
	248	IF( jperio == 5 .OR. jperio == 6 ) THEN ! F-point pivot
[2528]	249	tpol( 1 :jpiglo) = 0._wp
	250	fpol(jpiglo/2+1:jpiglo) = 0._wp
[3]	251	ENDIF
	252
	253	! 2. Ocean/land mask at u-, v-, and z-points (computed from tmask)
	254	! -------------------------------------------
	255	DO jk = 1, jpk
	256	DO jj = 1, jpjm1
	257	DO ji = 1, fs_jpim1 ! vector loop
	258	umask(ji,jj,jk) = tmask(ji,jj ,jk) * tmask(ji+1,jj ,jk)
	259	vmask(ji,jj,jk) = tmask(ji,jj ,jk) * tmask(ji ,jj+1,jk)
[1271]	260	END DO
[1694]	261	DO ji = 1, jpim1 ! NO vector opt.
[3]	262	fmask(ji,jj,jk) = tmask(ji,jj ,jk) * tmask(ji+1,jj ,jk) &
[62]	263	& * tmask(ji,jj+1,jk) * tmask(ji+1,jj+1,jk)
[3]	264	END DO
	265	END DO
	266	END DO
[4990]	267	! (ISF) MIN(1,SUM(umask)) is here to check if you have effectively at least 1 wet u point
	268	DO jj = 1, jpjm1
	269	DO ji = 1, fs_jpim1 ! vector loop
	270	umask_i(ji,jj) = ssmask(ji,jj) * ssmask(ji+1,jj ) * MIN(1._wp,SUM(umask(ji,jj,:)))
	271	vmask_i(ji,jj) = ssmask(ji,jj) * ssmask(ji ,jj+1) * MIN(1._wp,SUM(vmask(ji,jj,:)))
	272	END DO
	273	DO ji = 1, jpim1 ! NO vector opt.
	274	fmask_i(ji,jj) = ssmask(ji,jj ) * ssmask(ji+1,jj ) &
	275	& * ssmask(ji,jj+1) * ssmask(ji+1,jj+1) * MIN(1._wp,SUM(fmask(ji,jj,:)))
	276	END DO
	277	END DO
[2528]	278	CALL lbc_lnk( umask, 'U', 1._wp ) ! Lateral boundary conditions
	279	CALL lbc_lnk( vmask, 'V', 1._wp )
	280	CALL lbc_lnk( fmask, 'F', 1._wp )
[4990]	281	CALL lbc_lnk( umask_i, 'U', 1._wp ) ! Lateral boundary conditions
	282	CALL lbc_lnk( vmask_i, 'V', 1._wp )
	283	CALL lbc_lnk( fmask_i, 'F', 1._wp )
[3]	284
[5120]	285	! 3. Ocean/land mask at wu-, wv- and w points
	286	!----------------------------------------------
	287	wmask (:,:,1) = tmask(:,:,1) ! ????????
	288	wumask(:,:,1) = umask(:,:,1) ! ????????
	289	wvmask(:,:,1) = vmask(:,:,1) ! ????????
	290	DO jk=2,jpk
	291	wmask (:,:,jk)=tmask(:,:,jk) * tmask(:,:,jk-1)
	292	wumask(:,:,jk)=umask(:,:,jk) * umask(:,:,jk-1)
	293	wvmask(:,:,jk)=vmask(:,:,jk) * vmask(:,:,jk-1)
	294	END DO
[3]	295
	296	! 4. ocean/land mask for the elliptic equation
	297	! --------------------------------------------
[4990]	298	bmask(:,:) = ssmask(:,:) ! elliptic equation is written at t-point
[1566]	299	!
	300	! ! Boundary conditions
	301	! ! cyclic east-west : bmask must be set to 0. on rows 1 and jpi
[3]	302	IF( nperio == 1 .OR. nperio == 4 .OR. nperio == 6 ) THEN
[2528]	303	bmask( 1 ,:) = 0._wp
	304	bmask(jpi,:) = 0._wp
[3]	305	ENDIF
[1566]	306	IF( nperio == 2 ) THEN ! south symmetric : bmask must be set to 0. on row 1
[2528]	307	bmask(:, 1 ) = 0._wp
[3]	308	ENDIF
[1566]	309	! ! north fold :
	310	IF( nperio == 3 .OR. nperio == 4 ) THEN ! T-pt pivot : bmask set to 0. on row jpj and on half jpjglo-1 row
	311	DO ji = 1, jpi
[1528]	312	ii = ji + nimpp - 1
	313	bmask(ji,jpj-1) = bmask(ji,jpj-1) * tpol(ii)
[2528]	314	bmask(ji,jpj ) = 0._wp
[1528]	315	END DO
[3]	316	ENDIF
[1566]	317	IF( nperio == 5 .OR. nperio == 6 ) THEN ! F-pt pivot and T-pt elliptic eq. : bmask set to 0. on row jpj
[2528]	318	bmask(:,jpj) = 0._wp
[3]	319	ENDIF
[1566]	320	!
	321	IF( lk_mpp ) THEN ! mpp specificities
	322	! ! bmask is set to zero on the overlap region
[2528]	323	IF( nbondi /= -1 .AND. nbondi /= 2 ) bmask( 1 :jpreci,:) = 0._wp
	324	IF( nbondi /= 1 .AND. nbondi /= 2 ) bmask(nlci:jpi ,:) = 0._wp
	325	IF( nbondj /= -1 .AND. nbondj /= 2 ) bmask(:, 1 :jprecj) = 0._wp
	326	IF( nbondj /= 1 .AND. nbondj /= 2 ) bmask(:,nlcj:jpj ) = 0._wp
[1566]	327	!
	328	IF( npolj == 3 .OR. npolj == 4 ) THEN ! north fold : bmask must be set to 0. on rows jpj-1 and jpj
[1528]	329	DO ji = 1, nlci
	330	ii = ji + nimpp - 1
	331	bmask(ji,nlcj-1) = bmask(ji,nlcj-1) * tpol(ii)
[2528]	332	bmask(ji,nlcj ) = 0._wp
[1528]	333	END DO
[32]	334	ENDIF
[1566]	335	IF( npolj == 5 .OR. npolj == 6 ) THEN ! F-pt pivot and T-pt elliptic eq. : bmask set to 0. on row jpj
[1528]	336	DO ji = 1, nlci
[2528]	337	bmask(ji,nlcj ) = 0._wp
[1528]	338	END DO
[32]	339	ENDIF
[3]	340	ENDIF
	341
	342
	343	! mask for second order calculation of vorticity
	344	! ----------------------------------------------
	345	CALL dom_msk_nsa
	346
	347
	348	! Lateral boundary conditions on velocity (modify fmask)
[1566]	349	! ---------------------------------------
[3]	350	DO jk = 1, jpk
[1566]	351	zwf(:,:) = fmask(:,:,jk)
[3]	352	DO jj = 2, jpjm1
	353	DO ji = fs_2, fs_jpim1 ! vector opt.
[2528]	354	IF( fmask(ji,jj,jk) == 0._wp ) THEN
	355	fmask(ji,jj,jk) = rn_shlat * MIN( 1._wp , MAX( zwf(ji+1,jj), zwf(ji,jj+1), &
	356	& zwf(ji-1,jj), zwf(ji,jj-1) ) )
[3]	357	ENDIF
	358	END DO
	359	END DO
	360	DO jj = 2, jpjm1
[2528]	361	IF( fmask(1,jj,jk) == 0._wp ) THEN
	362	fmask(1 ,jj,jk) = rn_shlat * MIN( 1._wp , MAX( zwf(2,jj), zwf(1,jj+1), zwf(1,jj-1) ) )
[3]	363	ENDIF
[2528]	364	IF( fmask(jpi,jj,jk) == 0._wp ) THEN
	365	fmask(jpi,jj,jk) = rn_shlat * MIN( 1._wp , MAX( zwf(jpi,jj+1), zwf(jpim1,jj), zwf(jpi,jj-1) ) )
[3]	366	ENDIF
[1566]	367	END DO
[3]	368	DO ji = 2, jpim1
[2528]	369	IF( fmask(ji,1,jk) == 0._wp ) THEN
	370	fmask(ji, 1 ,jk) = rn_shlat * MIN( 1._wp , MAX( zwf(ji+1,1), zwf(ji,2), zwf(ji-1,1) ) )
[3]	371	ENDIF
[2528]	372	IF( fmask(ji,jpj,jk) == 0._wp ) THEN
	373	fmask(ji,jpj,jk) = rn_shlat * MIN( 1._wp , MAX( zwf(ji+1,jpj), zwf(ji-1,jpj), zwf(ji,jpjm1) ) )
[3]	374	ENDIF
	375	END DO
	376	END DO
[1566]	377	!
	378	IF( cp_cfg == "orca" .AND. jp_cfg == 2 ) THEN ! ORCA_R2 configuration
	379	! ! Increased lateral friction near of some straits
[2528]	380	IF( nn_cla == 0 ) THEN
[1273]	381	! ! Gibraltar strait : partial slip (fmask=0.5)
	382	ij0 = 101 ; ij1 = 101
[2528]	383	ii0 = 139 ; ii1 = 140 ; fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) = 0.5_wp
[32]	384	ij0 = 102 ; ij1 = 102
[2528]	385	ii0 = 139 ; ii1 = 140 ; fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) = 0.5_wp
[1273]	386	!
	387	! ! Bab el Mandeb : partial slip (fmask=1)
	388	ij0 = 87 ; ij1 = 88
[2528]	389	ii0 = 160 ; ii1 = 160 ; fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) = 1._wp
[1273]	390	ij0 = 88 ; ij1 = 88
[2528]	391	ii0 = 159 ; ii1 = 159 ; fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) = 1._wp
[1273]	392	!
[3]	393	ENDIF
[1707]	394	! ! Danish straits : strong slip (fmask > 2)
	395	! We keep this as an example but it is instable in this case
	396	! ij0 = 115 ; ij1 = 115
[2528]	397	! ii0 = 145 ; ii1 = 146 ; fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) = 4._wp
[1707]	398	! ij0 = 116 ; ij1 = 116
[2528]	399	! ii0 = 145 ; ii1 = 146 ; fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) = 4._wp
[3]	400	!
	401	ENDIF
[1566]	402	!
[2528]	403	IF( cp_cfg == "orca" .AND. jp_cfg == 1 ) THEN ! ORCA R1 configuration
	404	! ! Increased lateral friction near of some straits
[5506]	405	! This dirty section will be suppressed by simplification process:
	406	! all this will come back in input files
	407	! Currently these hard-wired indices relate to configuration with
	408	! extend grid (jpjglo=332)
	409	!
	410	isrow = 332 - jpjglo
	411	!
[2528]	412	IF(lwp) WRITE(numout,*)
	413	IF(lwp) WRITE(numout,*) ' orca_r1: increase friction near the following straits : '
	414	IF(lwp) WRITE(numout,*) ' Gibraltar '
[5385]	415	ii0 = 282 ; ii1 = 283 ! Gibraltar Strait
[5551]	416	ij0 = 241 - isrow ; ij1 = 241 - isrow ; fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp
[1566]	417
[2528]	418	IF(lwp) WRITE(numout,*) ' Bhosporus '
[5385]	419	ii0 = 314 ; ii1 = 315 ! Bhosporus Strait
[5551]	420	ij0 = 248 - isrow ; ij1 = 248 - isrow ; fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp
[2528]	421
	422	IF(lwp) WRITE(numout,*) ' Makassar (Top) '
[5385]	423	ii0 = 48 ; ii1 = 48 ! Makassar Strait (Top)
[5551]	424	ij0 = 189 - isrow ; ij1 = 190 - isrow ; fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 3._wp
[2528]	425
	426	IF(lwp) WRITE(numout,*) ' Lombok '
[5385]	427	ii0 = 44 ; ii1 = 44 ! Lombok Strait
[5551]	428	ij0 = 164 - isrow ; ij1 = 165 - isrow ; fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp
[2528]	429
	430	IF(lwp) WRITE(numout,*) ' Ombai '
[5385]	431	ii0 = 53 ; ii1 = 53 ! Ombai Strait
[5551]	432	ij0 = 164 - isrow ; ij1 = 165 - isrow ; fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp
[2528]	433
	434	IF(lwp) WRITE(numout,*) ' Timor Passage '
[5385]	435	ii0 = 56 ; ii1 = 56 ! Timor Passage
[5551]	436	ij0 = 164 - isrow ; ij1 = 165 - isrow ; fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp
[2528]	437
	438	IF(lwp) WRITE(numout,*) ' West Halmahera '
[5385]	439	ii0 = 58 ; ii1 = 58 ! West Halmahera Strait
[5551]	440	ij0 = 181 - isrow ; ij1 = 182 - isrow ; fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 3._wp
[2528]	441
	442	IF(lwp) WRITE(numout,*) ' East Halmahera '
[5385]	443	ii0 = 55 ; ii1 = 55 ! East Halmahera Strait
[5551]	444	ij0 = 181 - isrow ; ij1 = 182 - isrow ; fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 3._wp
[2528]	445	!
	446	ENDIF
	447	!
	448	CALL lbc_lnk( fmask, 'F', 1._wp ) ! Lateral boundary conditions on fmask
	449
[3294]	450	! CAUTION : The fmask may be further modified in dyn_vor_init ( dynvor.F90 )
[2528]	451
[1566]	452	IF( nprint == 1 .AND. lwp ) THEN ! Control print
[3]	453	imsk(:,:) = INT( tmask_i(:,:) )
	454	WRITE(numout,*) ' tmask_i : '
	455	CALL prihin( imsk(:,:), jpi, jpj, 1, jpi, 1, &
	456	& 1, jpj, 1, 1, numout)
	457	WRITE (numout,*)
	458	WRITE (numout,*) ' dommsk: tmask for each level'
	459	WRITE (numout,*) ' ----------------------------'
	460	DO jk = 1, jpk
	461	imsk(:,:) = INT( tmask(:,:,jk) )
	462
	463	WRITE(numout,*)
	464	WRITE(numout,*) ' level = ',jk
	465	CALL prihin( imsk(:,:), jpi, jpj, 1, jpi, 1, &
	466	& 1, jpj, 1, 1, numout)
	467	END DO
	468	WRITE(numout,*)
	469	WRITE(numout,*) ' dom_msk: vmask for each level'
	470	WRITE(numout,*) ' -----------------------------'
	471	DO jk = 1, jpk
	472	imsk(:,:) = INT( vmask(:,:,jk) )
	473	WRITE(numout,*)
	474	WRITE(numout,*) ' level = ',jk
	475	CALL prihin( imsk(:,:), jpi, jpj, 1, jpi, 1, &
	476	& 1, jpj, 1, 1, numout)
	477	END DO
	478	WRITE(numout,*)
	479	WRITE(numout,*) ' dom_msk: fmask for each level'
	480	WRITE(numout,*) ' -----------------------------'
	481	DO jk = 1, jpk
	482	imsk(:,:) = INT( fmask(:,:,jk) )
	483	WRITE(numout,*)
	484	WRITE(numout,*) ' level = ',jk
	485	CALL prihin( imsk(:,:), jpi, jpj, 1, jpi, 1, &
	486	& 1, jpj, 1, 1, numout )
	487	END DO
	488	WRITE(numout,*)
	489	WRITE(numout,*) ' dom_msk: bmask '
	490	WRITE(numout,*) ' ---------------'
	491	WRITE(numout,*)
	492	imsk(:,:) = INT( bmask(:,:) )
	493	CALL prihin( imsk(:,:), jpi, jpj, 1, jpi, 1, &
[2528]	494	& 1, jpj, 1, 1, numout )
[3]	495	ENDIF
[1566]	496	!
[3294]	497	CALL wrk_dealloc( jpi, jpj, imsk )
	498	CALL wrk_dealloc( jpi, jpj, zwf )
[2715]	499	!
[3294]	500	IF( nn_timing == 1 ) CALL timing_stop('dom_msk')
	501	!
[3]	502	END SUBROUTINE dom_msk
	503
	504	#if defined key_noslip_accurate
	505	!!----------------------------------------------------------------------
	506	!! 'key_noslip_accurate' : accurate no-slip boundary condition
	507	!!----------------------------------------------------------------------
	508
	509	SUBROUTINE dom_msk_nsa
	510	!!---------------------------------------------------------------------
	511	!! * ROUTINE dom_msk_nsa *
	512	!!
	513	!! ** Purpose :
	514	!!
	515	!! ** Method :
	516	!!
	517	!! ** Action :
	518	!!----------------------------------------------------------------------
[2715]	519	INTEGER :: ji, jj, jk, jl ! dummy loop indices
[1566]	520	INTEGER :: ine, inw, ins, inn, itest, ierror, iind, ijnd
	521	REAL(wp) :: zaa
[3]	522	!!---------------------------------------------------------------------
[3294]	523	!
	524	IF( nn_timing == 1 ) CALL timing_start('dom_msk_nsa')
	525	!
[2715]	526	IF(lwp) WRITE(numout,*)
	527	IF(lwp) WRITE(numout,*) 'dom_msk_nsa : noslip accurate boundary condition'
	528	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ using Schchepetkin and O Brian scheme'
	529	IF( lk_mpp ) CALL ctl_stop( ' mpp version is not yet implemented' )
[3]	530
	531	! mask for second order calculation of vorticity
	532	! ----------------------------------------------
	533	! noslip boundary condition: fmask=1 at convex corner, store
	534	! index of straight coast meshes ( 'west', refering to a coast,
	535	! means west of the ocean, aso)
	536
	537	DO jk = 1, jpk
	538	DO jl = 1, 4
	539	npcoa(jl,jk) = 0
	540	DO ji = 1, 2*(jpi+jpj)
	541	nicoa(ji,jl,jk) = 0
	542	njcoa(ji,jl,jk) = 0
	543	END DO
	544	END DO
	545	END DO
	546
	547	IF( jperio == 2 ) THEN
	548	WRITE(numout,*) ' '
	549	WRITE(numout,*) ' symetric boundary conditions need special'
	550	WRITE(numout,*) ' treatment not implemented. we stop.'
	551	STOP
	552	ENDIF
	553
	554	! convex corners
	555
	556	DO jk = 1, jpkm1
	557	DO jj = 1, jpjm1
	558	DO ji = 1, jpim1
	559	zaa = tmask(ji ,jj,jk) + tmask(ji ,jj+1,jk) &
[32]	560	&+ tmask(ji+1,jj,jk) + tmask(ji+1,jj+1,jk)
[2528]	561	IF( ABS(zaa-3._wp) <= 0.1_wp ) fmask(ji,jj,jk) = 1._wp
[3]	562	END DO
	563	END DO
	564	END DO
	565
	566	! north-south straight coast
	567
	568	DO jk = 1, jpkm1
	569	inw = 0
	570	ine = 0
	571	DO jj = 2, jpjm1
	572	DO ji = 2, jpim1
	573	zaa = tmask(ji+1,jj,jk) + tmask(ji+1,jj+1,jk)
[2528]	574	IF( ABS(zaa-2._wp) <= 0.1_wp .AND. fmask(ji,jj,jk) == 0._wp ) THEN
[3]	575	inw = inw + 1
	576	nicoa(inw,1,jk) = ji
	577	njcoa(inw,1,jk) = jj
	578	IF( nprint == 1 ) WRITE(numout,*) ' west : ', jk, inw, ji, jj
	579	ENDIF
	580	zaa = tmask(ji,jj,jk) + tmask(ji,jj+1,jk)
[2528]	581	IF( ABS(zaa-2._wp) <= 0.1_wp .AND. fmask(ji,jj,jk) == 0._wp ) THEN
[3]	582	ine = ine + 1
	583	nicoa(ine,2,jk) = ji
	584	njcoa(ine,2,jk) = jj
	585	IF( nprint == 1 ) WRITE(numout,*) ' east : ', jk, ine, ji, jj
	586	ENDIF
	587	END DO
	588	END DO
	589	npcoa(1,jk) = inw
	590	npcoa(2,jk) = ine
	591	END DO
	592
	593	! west-east straight coast
	594
	595	DO jk = 1, jpkm1
	596	ins = 0
	597	inn = 0
	598	DO jj = 2, jpjm1
	599	DO ji =2, jpim1
	600	zaa = tmask(ji,jj+1,jk) + tmask(ji+1,jj+1,jk)
[2528]	601	IF( ABS(zaa-2._wp) <= 0.1_wp .AND. fmask(ji,jj,jk) == 0._wp ) THEN
[3]	602	ins = ins + 1
	603	nicoa(ins,3,jk) = ji
	604	njcoa(ins,3,jk) = jj
	605	IF( nprint == 1 ) WRITE(numout,*) ' south : ', jk, ins, ji, jj
	606	ENDIF
	607	zaa = tmask(ji+1,jj,jk) + tmask(ji,jj,jk)
[2528]	608	IF( ABS(zaa-2._wp) <= 0.1_wp .AND. fmask(ji,jj,jk) == 0._wp ) THEN
[3]	609	inn = inn + 1
	610	nicoa(inn,4,jk) = ji
	611	njcoa(inn,4,jk) = jj
	612	IF( nprint == 1 ) WRITE(numout,*) ' north : ', jk, inn, ji, jj
	613	ENDIF
	614	END DO
	615	END DO
	616	npcoa(3,jk) = ins
	617	npcoa(4,jk) = inn
	618	END DO
	619
	620	itest = 2 * ( jpi + jpj )
	621	DO jk = 1, jpk
	622	IF( npcoa(1,jk) > itest .OR. npcoa(2,jk) > itest .OR. &
	623	npcoa(3,jk) > itest .OR. npcoa(4,jk) > itest ) THEN
[474]	624
	625	WRITE(ctmp1,*) ' level jk = ',jk
	626	WRITE(ctmp2,*) ' straight coast index arraies are too small.:'
	627	WRITE(ctmp3,*) ' npe, npw, nps, npn = ', npcoa(1,jk), npcoa(2,jk), &
[32]	628	& npcoa(3,jk), npcoa(4,jk)
[474]	629	WRITE(ctmp4,) ' 2(jpi+jpj) = ',itest,'. we stop.'
	630	CALL ctl_stop( ctmp1, ctmp2, ctmp3, ctmp4 )
[3]	631	ENDIF
	632	END DO
	633
	634	ierror = 0
	635	iind = 0
	636	ijnd = 0
[2528]	637	IF( nperio == 1 .OR. nperio == 4 .OR. nperio == 6 ) iind = 2
	638	IF( nperio == 3 .OR. nperio == 4 .OR. nperio == 5 .OR. nperio == 6 ) ijnd = 2
[3]	639	DO jk = 1, jpk
	640	DO jl = 1, npcoa(1,jk)
	641	IF( nicoa(jl,1,jk)+3 > jpi+iind ) THEN
	642	ierror = ierror+1
	643	icoord(ierror,1) = nicoa(jl,1,jk)
	644	icoord(ierror,2) = njcoa(jl,1,jk)
	645	icoord(ierror,3) = jk
	646	ENDIF
	647	END DO
	648	DO jl = 1, npcoa(2,jk)
	649	IF(nicoa(jl,2,jk)-2 < 1-iind ) THEN
	650	ierror = ierror + 1
	651	icoord(ierror,1) = nicoa(jl,2,jk)
	652	icoord(ierror,2) = njcoa(jl,2,jk)
	653	icoord(ierror,3) = jk
	654	ENDIF
	655	END DO
	656	DO jl = 1, npcoa(3,jk)
	657	IF( njcoa(jl,3,jk)+3 > jpj+ijnd ) THEN
	658	ierror = ierror + 1
	659	icoord(ierror,1) = nicoa(jl,3,jk)
	660	icoord(ierror,2) = njcoa(jl,3,jk)
	661	icoord(ierror,3) = jk
	662	ENDIF
	663	END DO
[2528]	664	DO jl = 1, npcoa(4,jk)
[3]	665	IF( njcoa(jl,4,jk)-2 < 1) THEN
[2528]	666	ierror=ierror + 1
	667	icoord(ierror,1) = nicoa(jl,4,jk)
	668	icoord(ierror,2) = njcoa(jl,4,jk)
	669	icoord(ierror,3) = jk
[3]	670	ENDIF
	671	END DO
	672	END DO
	673
	674	IF( ierror > 0 ) THEN
	675	IF(lwp) WRITE(numout,*)
	676	IF(lwp) WRITE(numout,*) ' Problem on lateral conditions'
	677	IF(lwp) WRITE(numout,*) ' Bad marking off at points:'
	678	DO jl = 1, ierror
	679	IF(lwp) WRITE(numout,*) 'Level:',icoord(jl,3), &
[32]	680	& ' Point(',icoord(jl,1),',',icoord(jl,2),')'
[3]	681	END DO
[474]	682	CALL ctl_stop( 'We stop...' )
[3]	683	ENDIF
[2715]	684	!
[3294]	685	IF( nn_timing == 1 ) CALL timing_stop('dom_msk_nsa')
	686	!
[3]	687	END SUBROUTINE dom_msk_nsa
	688
	689	#else
	690	!!----------------------------------------------------------------------
	691	!! Default option : Empty routine
	692	!!----------------------------------------------------------------------
	693	SUBROUTINE dom_msk_nsa
	694	END SUBROUTINE dom_msk_nsa
	695	#endif
	696
	697	!!======================================================================
	698	END MODULE dommsk

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

Download in other formats: