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.

wet_dry.F90 in NEMO/branches/2020/dev_r14116_HPC-04_mcastril_Mixed_Precision_implementation_final/src/OCE/DYN – NEMO

Context Navigation

source: NEMO/branches/2020/dev_r14116_HPC-04_mcastril_Mixed_Precision_implementation_final/src/OCE/DYN/wet_dry.F90 @ 14644

Last change on this file since 14644 was 14644, checked in by sparonuz, 3 years ago
Merge trunk -r14642:HEAD
Property svn:keywords set to `Id`
File size: 19.4 KB

Rev	Line
[6152]	1	MODULE wet_dry
[9023]	2
	3	!! includes updates to namelist namwad for diagnostic outputs of ROMS wetting and drying
	4
[6152]	5	!!==============================================================================
	6	!! * MODULE wet_dry *
	7	!! Wetting and drying includes initialisation routine and routines to
	8	!! compute and apply flux limiters and preserve water depth positivity
[9023]	9	!! only effects if wetting/drying is on (ln_wd_il == .true. or ln_wd_dl==.true. )
[6152]	10	!!==============================================================================
[7646]	11	!! History : 3.6 ! 2014-09 ((H.Liu) Original code
[6152]	12	!! ! will add the runoff and periodic BC case later
	13	!!----------------------------------------------------------------------
	14
	15	!!----------------------------------------------------------------------
[9019]	16	!! wad_init : initialisation of wetting and drying
	17	!! wad_lmt : horizontal flux limiter and limited velocity when wetting and drying happens
	18	!! wad_lmt_bt : same as wad_lmt for the barotropic stepping (dynspg_ts)
[6152]	19	!!----------------------------------------------------------------------
[9019]	20	USE oce ! ocean dynamics and tracers
	21	USE dom_oce ! ocean space and time domain
	22	USE sbc_oce , ONLY: ln_rnf ! surface boundary condition: ocean
	23	USE sbcrnf ! river runoff
	24	!
	25	USE in_out_manager ! I/O manager
	26	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
	27	USE lib_mpp ! MPP library
[9124]	28	USE timing ! timing of the main modules
[6152]	29
	30	IMPLICIT NONE
	31	PRIVATE
	32
[12377]	33	!! * Substitutions
	34	# include "do_loop_substitute.h90"
[13237]	35	# include "domzgr_substitute.h90"
[6152]	36	!!----------------------------------------------------------------------
	37	!! critical depths,filters, limiters,and masks for Wetting and Drying
	38	!! ---------------------------------------------------------------------
	39
[9019]	40	REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: wdmask !: u- and v- limiter
	41	! ! (can include negative depths)
[9023]	42	REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: wdramp, wdrampu, wdrampv !: for hpg limiting
[6152]	43
[9023]	44	LOGICAL, PUBLIC :: ln_wd_il !: Wetting/drying il activation switch (T:on,F:off)
	45	LOGICAL, PUBLIC :: ln_wd_dl !: Wetting/drying dl activation switch (T:on,F:off)
	46	REAL(wp), PUBLIC :: rn_wdmin0 !: depth at which wetting/drying starts
[6152]	47	REAL(wp), PUBLIC :: rn_wdmin1 !: minimum water depth on dried cells
[9023]	48	REAL(wp), PUBLIC :: r_rn_wdmin1 !: 1/minimum water depth on dried cells
	49	REAL(wp), PUBLIC :: rn_wdmin2 !: tolerance of minimum water depth on dried cells
[10499]	50	REAL(wp), PUBLIC :: rn_wd_sbcdep !: Depth at which to taper sbc fluxes
	51	REAL(wp), PUBLIC :: rn_wd_sbcfra !: Fraction of SBC at taper depth
[9019]	52	REAL(wp), PUBLIC :: rn_wdld !: land elevation below which wetting/drying will be considered
[6152]	53	INTEGER , PUBLIC :: nn_wdit !: maximum number of iteration for W/D limiter
[9023]	54	LOGICAL, PUBLIC :: ln_wd_dl_bc !: DL scheme: True implies 3D velocities are set to the barotropic values at points
	55	!: where the flow is from wet points on less than half the barotropic sub-steps
	56	LOGICAL, PUBLIC :: ln_wd_dl_rmp !: use a ramp for the dl flux limiter between 2 rn_wdmin1 and rn_wdmin1 (rather than a cut-off at rn_wdmin1)
	57	REAL(wp), PUBLIC :: ssh_ref !: height of z=0 with respect to the geoid;
[6152]	58
[13558]	59	LOGICAL, PUBLIC :: ll_wd = .FALSE. !: Wetting/drying activation switch (ln_wd_il or ln_wd_dl) <- default def if wad_init not called
[9023]	60
[6152]	61	PUBLIC wad_init ! initialisation routine called by step.F90
	62	PUBLIC wad_lmt ! routine called by sshwzv.F90
	63	PUBLIC wad_lmt_bt ! routine called by dynspg_ts.F90
	64
	65	!! * Substitutions
[9019]	66	!!----------------------------------------------------------------------
[6152]	67	CONTAINS
	68
	69	SUBROUTINE wad_init
	70	!!----------------------------------------------------------------------
	71	!! * ROUTINE wad_init *
	72	!!
	73	!! ** Purpose : read wetting and drying namelist and print the variables.
	74	!!
	75	!! ** input : - namwad namelist
	76	!!----------------------------------------------------------------------
[9124]	77	INTEGER :: ios, ierr ! Local integer
[9019]	78	!!
[9124]	79	NAMELIST/namwad/ ln_wd_il, ln_wd_dl , rn_wdmin0, rn_wdmin1, rn_wdmin2, rn_wdld, &
[10499]	80	& nn_wdit , ln_wd_dl_bc, ln_wd_dl_rmp, rn_wd_sbcdep,rn_wd_sbcfra
[6152]	81	!!----------------------------------------------------------------------
[9019]	82	!
[6152]	83	READ ( numnam_ref, namwad, IOSTAT = ios, ERR = 905)
[11536]	84	905 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namwad in reference namelist' )
[6152]	85	READ ( numnam_cfg, namwad, IOSTAT = ios, ERR = 906)
[11536]	86	906 IF( ios > 0 ) CALL ctl_nam ( ios , 'namwad in configuration namelist' )
[6152]	87	IF(lwm) WRITE ( numond, namwad )
[9019]	88	!
[10499]	89	IF( rn_wd_sbcfra>=1 ) CALL ctl_stop( 'STOP', 'rn_wd_sbcfra >=1 : must be < 1' )
[6152]	90	IF(lwp) THEN ! control print
	91	WRITE(numout,*)
	92	WRITE(numout,*) 'wad_init : Wetting and drying initialization through namelist read'
[7646]	93	WRITE(numout,*) '~~~~~~~~'
[6152]	94	WRITE(numout,*) ' Namelist namwad'
[9023]	95	WRITE(numout,*) ' Logical for Iter Lim wd option ln_wd_il = ', ln_wd_il
	96	WRITE(numout,*) ' Logical for Dir. Lim wd option ln_wd_dl = ', ln_wd_dl
	97	WRITE(numout,*) ' Depth at which wet/drying starts rn_wdmin0 = ', rn_wdmin0
[6152]	98	WRITE(numout,*) ' Minimum wet depth on dried cells rn_wdmin1 = ', rn_wdmin1
[9023]	99	WRITE(numout,*) ' Tolerance of min wet depth rn_wdmin2 = ', rn_wdmin2
[6152]	100	WRITE(numout,*) ' land elevation threshold rn_wdld = ', rn_wdld
	101	WRITE(numout,*) ' Max iteration for W/D limiter nn_wdit = ', nn_wdit
[9023]	102	WRITE(numout,*) ' T => baroclinic u,v=0 at dry pts: ln_wd_dl_bc = ', ln_wd_dl_bc
	103	WRITE(numout,*) ' use a ramp for rwd limiter: ln_wd_dl_rwd_rmp = ', ln_wd_dl_rmp
[10499]	104	WRITE(numout,*) ' cut off depth sbc for wd rn_wd_sbcdep = ', rn_wd_sbcdep
	105	WRITE(numout,*) ' fraction to start sbc wgt rn_wd_sbcfra = ', rn_wd_sbcfra
[7646]	106	ENDIF
[9023]	107	IF( .NOT. ln_read_cfg ) THEN
[9042]	108	IF(lwp) WRITE(numout,*) ' No configuration file so seting ssh_ref to zero '
[9124]	109	ssh_ref=0._wp
[9023]	110	ENDIF
	111
[9124]	112	r_rn_wdmin1 = 1 / rn_wdmin1
	113	IF( ln_wd_il .OR. ln_wd_dl ) THEN
[9023]	114	ll_wd = .TRUE.
	115	ALLOCATE( wdmask(jpi,jpj), STAT=ierr )
	116	ALLOCATE( wdramp(jpi,jpj), wdrampu(jpi,jpj), wdrampv(jpi,jpj), STAT=ierr )
[6152]	117	IF( ierr /= 0 ) CALL ctl_stop('STOP', 'wad_init : Array allocation error')
	118	ENDIF
[7646]	119	!
[6152]	120	END SUBROUTINE wad_init
	121
[7646]	122
[12377]	123	SUBROUTINE wad_lmt( psshb1, psshemp, z2dt, Kmm, puu, pvv )
[6152]	124	!!----------------------------------------------------------------------
	125	!! * ROUTINE wad_lmt *
	126	!!
	127	!! ** Purpose : generate flux limiters for wetting/drying
	128	!!
	129	!! ** Method : - Prevent negative depth occurring (Not ready for Agrif)
	130	!!
	131	!! ** Action : - calculate flux limiter and W/D flag
	132	!!----------------------------------------------------------------------
[14219]	133	REAL(dp), DIMENSION(:,:) , INTENT(inout) :: psshb1
[12377]	134	REAL(wp), DIMENSION(:,:) , INTENT(in ) :: psshemp
	135	REAL(wp) , INTENT(in ) :: z2dt
	136	INTEGER , INTENT(in ) :: Kmm ! time level index
[14219]	137	REAL(dp), DIMENSION(jpi,jpj,jpk,jpt), INTENT(inout) :: puu, pvv ! velocity arrays
[6152]	138	!
	139	INTEGER :: ji, jj, jk, jk1 ! dummy loop indices
[7646]	140	INTEGER :: jflag ! local scalar
[6152]	141	REAL(wp) :: zcoef, zdep1, zdep2 ! local scalars
	142	REAL(wp) :: zzflxp, zzflxn ! local scalars
	143	REAL(wp) :: zdepwd ! local scalar, always wet cell depth
	144	REAL(wp) :: ztmp ! local scalars
[9019]	145	REAL(wp), DIMENSION(jpi,jpj) :: zwdlmtu, zwdlmtv ! W/D flux limiters
	146	REAL(wp), DIMENSION(jpi,jpj) :: zflxp , zflxn ! local 2D workspace
	147	REAL(wp), DIMENSION(jpi,jpj) :: zflxu , zflxv ! local 2D workspace
	148	REAL(wp), DIMENSION(jpi,jpj) :: zflxu1 , zflxv1 ! local 2D workspace
[6152]	149	!!----------------------------------------------------------------------
[9124]	150	IF( ln_timing ) CALL timing_start('wad_lmt') !
[6152]	151	!
[9023]	152	DO jk = 1, jpkm1
[12377]	153	puu(:,:,jk,Kmm) = puu(:,:,jk,Kmm) * zwdlmtu(:,:)
	154	pvv(:,:,jk,Kmm) = pvv(:,:,jk,Kmm) * zwdlmtv(:,:)
[9023]	155	END DO
[9019]	156	jflag = 0
[9124]	157	zdepwd = 50._wp ! maximum depth on which that W/D could possibly happen
	158	!
[9019]	159	zflxp(:,:) = 0._wp
	160	zflxn(:,:) = 0._wp
	161	zflxu(:,:) = 0._wp
	162	zflxv(:,:) = 0._wp
[9124]	163	!
	164	zwdlmtu(:,:) = 1._wp
	165	zwdlmtv(:,:) = 1._wp
	166	!
	167	DO jk = 1, jpkm1 ! Horizontal Flux in u and v direction
[12377]	168	zflxu(:,:) = zflxu(:,:) + e3u(:,:,jk,Kmm) * puu(:,:,jk,Kmm) * umask(:,:,jk)
	169	zflxv(:,:) = zflxv(:,:) + e3v(:,:,jk,Kmm) * pvv(:,:,jk,Kmm) * vmask(:,:,jk)
[9019]	170	END DO
	171	zflxu(:,:) = zflxu(:,:) * e2u(:,:)
	172	zflxv(:,:) = zflxv(:,:) * e1v(:,:)
[9124]	173	!
	174	wdmask(:,:) = 1._wp
[13295]	175	DO_2D( 0, 1, 0, 1 )
[12377]	176	!
	177	IF( tmask(ji,jj,1) < 0.5_wp ) CYCLE ! we don't care about land cells
	178	IF( ht_0(ji,jj) - ssh_ref > zdepwd ) CYCLE ! and cells which are unlikely to dry
	179	!
	180	zflxp(ji,jj) = MAX( zflxu(ji,jj) , 0._wp ) - MIN( zflxu(ji-1,jj ) , 0._wp ) &
	181	& + MAX( zflxv(ji,jj) , 0._wp ) - MIN( zflxv(ji, jj-1) , 0._wp )
	182	zflxn(ji,jj) = MIN( zflxu(ji,jj) , 0._wp ) - MAX( zflxu(ji-1,jj ) , 0._wp ) &
	183	& + MIN( zflxv(ji,jj) , 0._wp ) - MAX( zflxv(ji, jj-1) , 0._wp )
	184	!
	185	zdep2 = ht_0(ji,jj) + psshb1(ji,jj) - rn_wdmin1
	186	IF( zdep2 <= 0._wp ) THEN ! add more safty, but not necessary
	187	psshb1(ji,jj) = rn_wdmin1 - ht_0(ji,jj)
	188	IF(zflxu(ji, jj) > 0._wp) zwdlmtu(ji ,jj) = 0._wp
	189	IF(zflxu(ji-1,jj) < 0._wp) zwdlmtu(ji-1,jj) = 0._wp
	190	IF(zflxv(ji, jj) > 0._wp) zwdlmtv(ji ,jj) = 0._wp
	191	IF(zflxv(ji,jj-1) < 0._wp) zwdlmtv(ji,jj-1) = 0._wp
	192	wdmask(ji,jj) = 0._wp
	193	END IF
	194	END_2D
[9124]	195	!
	196	! ! HPG limiter from jholt
[12377]	197	wdramp(:,:) = min((ht_0(:,:) + psshb1(:,:) - rn_wdmin1)/(rn_wdmin0 - rn_wdmin1),1.0_wp)
[9124]	198	!jth assume don't need a lbc_lnk here
[13295]	199	DO_2D( 1, 0, 1, 0 )
[12377]	200	wdrampu(ji,jj) = MIN( wdramp(ji,jj) , wdramp(ji+1,jj) )
	201	wdrampv(ji,jj) = MIN( wdramp(ji,jj) , wdramp(ji,jj+1) )
	202	END_2D
[9124]	203	! ! end HPG limiter
	204	!
	205	!
	206	DO jk1 = 1, nn_wdit + 1 !== start limiter iterations ==!
	207	!
[9019]	208	zflxu1(:,:) = zflxu(:,:) * zwdlmtu(:,:)
	209	zflxv1(:,:) = zflxv(:,:) * zwdlmtv(:,:)
	210	jflag = 0 ! flag indicating if any further iterations are needed
[9124]	211	!
[13295]	212	DO_2D( 0, 1, 0, 1 )
[12377]	213	IF( tmask(ji, jj, 1) < 0.5_wp ) CYCLE
	214	IF( ht_0(ji,jj) > zdepwd ) CYCLE
	215	!
	216	ztmp = e1e2t(ji,jj)
	217	!
	218	zzflxp = MAX( zflxu1(ji,jj) , 0._wp ) - MIN( zflxu1(ji-1,jj ) , 0._wp) &
	219	& + MAX( zflxv1(ji,jj) , 0._wp ) - MIN( zflxv1(ji, jj-1) , 0._wp)
	220	zzflxn = MIN( zflxu1(ji,jj) , 0._wp ) - MAX( zflxu1(ji-1,jj ) , 0._wp) &
	221	& + MIN( zflxv1(ji,jj) , 0._wp ) - MAX( zflxv1(ji, jj-1) , 0._wp)
	222	!
	223	zdep1 = (zzflxp + zzflxn) * z2dt / ztmp
	224	zdep2 = ht_0(ji,jj) + psshb1(ji,jj) - rn_wdmin1 - z2dt * psshemp(ji,jj)
	225	!
	226	IF( zdep1 > zdep2 ) THEN
	227	wdmask(ji, jj) = 0._wp
	228	zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zflxp(ji,jj) * z2dt )
	229	!zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zzflxp * z2dt )
	230	! flag if the limiter has been used but stop flagging if the only
	231	! changes have zeroed the coefficient since further iterations will
	232	! not change anything
	233	IF( zcoef > 0._wp ) THEN ; jflag = 1
	234	ELSE ; zcoef = 0._wp
[9124]	235	ENDIF
[12377]	236	IF( jk1 > nn_wdit ) zcoef = 0._wp
	237	IF( zflxu1(ji ,jj ) > 0._wp ) zwdlmtu(ji ,jj ) = zcoef
	238	IF( zflxu1(ji-1,jj ) < 0._wp ) zwdlmtu(ji-1,jj ) = zcoef
	239	IF( zflxv1(ji ,jj ) > 0._wp ) zwdlmtv(ji ,jj ) = zcoef
	240	IF( zflxv1(ji ,jj-1) < 0._wp ) zwdlmtv(ji ,jj-1) = zcoef
	241	ENDIF
	242	END_2D
[14644]	243	CALL lbc_lnk( 'wet_dry', zwdlmtu, 'U', 1.0_wp, zwdlmtv, 'V', 1.0_wp )
[9124]	244	!
[10425]	245	CALL mpp_max('wet_dry', jflag) !max over the global domain
[9124]	246	!
	247	IF( jflag == 0 ) EXIT
	248	!
[9019]	249	END DO ! jk1 loop
[9124]	250	!
[9019]	251	DO jk = 1, jpkm1
[12377]	252	puu(:,:,jk,Kmm) = puu(:,:,jk,Kmm) * zwdlmtu(:,:)
	253	pvv(:,:,jk,Kmm) = pvv(:,:,jk,Kmm) * zwdlmtv(:,:)
[9019]	254	END DO
[12377]	255	uu_b(:,:,Kmm) = uu_b(:,:,Kmm) * zwdlmtu(:, :)
	256	vv_b(:,:,Kmm) = vv_b(:,:,Kmm) * zwdlmtv(:, :)
[9124]	257	!
	258	!!gm TO BE SUPPRESSED ? these lbc_lnk are useless since zwdlmtu and zwdlmtv are defined everywhere !
[14644]	259	CALL lbc_lnk( 'wet_dry', puu(:,:,:,Kmm) , 'U', -1.0_wp, pvv(:,:,:,Kmm) , 'V', -1.0_wp )
	260	CALL lbc_lnk( 'wet_dry', uu_b(:,:,Kmm), 'U', -1.0_wp, vv_b(:,:,Kmm), 'V', -1.0_wp )
[9124]	261	!!gm
[6152]	262	!
[9124]	263	IF(jflag == 1 .AND. lwp) WRITE(numout,*) 'Need more iterations in wad_lmt!!!'
[7646]	264	!
[12377]	265	!IF( ln_rnf ) CALL sbc_rnf_div( hdiv ) ! runoffs (update hdiv field)
[9023]	266	!
[9124]	267	IF( ln_timing ) CALL timing_stop('wad_lmt') !
[9023]	268	!
[6152]	269	END SUBROUTINE wad_lmt
	270
[7646]	271
[12489]	272	SUBROUTINE wad_lmt_bt( zflxu, zflxv, sshn_e, zssh_frc, rDt_e )
[6152]	273	!!----------------------------------------------------------------------
	274	!! * ROUTINE wad_lmt *
	275	!!
	276	!! ** Purpose : limiting flux in the barotropic stepping (dynspg_ts)
	277	!!
	278	!! ** Method : - Prevent negative depth occurring (Not ready for Agrif)
	279	!!
	280	!! ** Action : - calculate flux limiter and W/D flag
	281	!!----------------------------------------------------------------------
[12489]	282	REAL(wp) , INTENT(in ) :: rDt_e ! ocean time-step index
[14219]	283	REAL(dp), DIMENSION(:,:), INTENT(inout) :: sshn_e, zssh_frc
	284	REAL(wp), DIMENSION(:,:), INTENT(inout) :: zflxu, zflxv
[6152]	285	!
	286	INTEGER :: ji, jj, jk, jk1 ! dummy loop indices
[9124]	287	INTEGER :: jflag ! local integer
[6152]	288	REAL(wp) :: z2dt
	289	REAL(wp) :: zcoef, zdep1, zdep2 ! local scalars
	290	REAL(wp) :: zzflxp, zzflxn ! local scalars
	291	REAL(wp) :: zdepwd ! local scalar, always wet cell depth
	292	REAL(wp) :: ztmp ! local scalars
[9019]	293	REAL(wp), DIMENSION(jpi,jpj) :: zwdlmtu, zwdlmtv !: W/D flux limiters
	294	REAL(wp), DIMENSION(jpi,jpj) :: zflxp, zflxn ! local 2D workspace
	295	REAL(wp), DIMENSION(jpi,jpj) :: zflxu1, zflxv1 ! local 2D workspace
[6152]	296	!!----------------------------------------------------------------------
[9124]	297	IF( ln_timing ) CALL timing_start('wad_lmt_bt') !
[6152]	298	!
[9124]	299	jflag = 0
	300	zdepwd = 50._wp ! maximum depth that ocean cells can have W/D processes
[9023]	301	!
[12489]	302	z2dt = rDt_e
[9124]	303	!
[9019]	304	zflxp(:,:) = 0._wp
	305	zflxn(:,:) = 0._wp
	306	zwdlmtu(:,:) = 1._wp
	307	zwdlmtv(:,:) = 1._wp
[9124]	308	!
[13497]	309	DO_2D( 0, 1, 0, 1 ) ! Horizontal Flux in u and v direction
[12377]	310	!
	311	IF( tmask(ji, jj, 1 ) < 0.5_wp) CYCLE ! we don't care about land cells
	312	IF( ht_0(ji,jj) > zdepwd ) CYCLE ! and cells which are unlikely to dry
	313	!
	314	zflxp(ji,jj) = MAX( zflxu(ji,jj) , 0._wp ) - MIN( zflxu(ji-1,jj ) , 0._wp ) &
	315	& + MAX( zflxv(ji,jj) , 0._wp ) - MIN( zflxv(ji, jj-1) , 0._wp )
	316	zflxn(ji,jj) = MIN( zflxu(ji,jj) , 0._wp ) - MAX( zflxu(ji-1,jj ) , 0._wp ) &
	317	& + MIN( zflxv(ji,jj) , 0._wp ) - MAX( zflxv(ji, jj-1) , 0._wp )
	318	!
	319	zdep2 = ht_0(ji,jj) + sshn_e(ji,jj) - rn_wdmin1
	320	IF( zdep2 <= 0._wp ) THEN !add more safety, but not necessary
	321	sshn_e(ji,jj) = rn_wdmin1 - ht_0(ji,jj)
	322	IF( zflxu(ji ,jj ) > 0._wp) zwdlmtu(ji ,jj ) = 0._wp
	323	IF( zflxu(ji-1,jj ) < 0._wp) zwdlmtu(ji-1,jj ) = 0._wp
	324	IF( zflxv(ji ,jj ) > 0._wp) zwdlmtv(ji ,jj ) = 0._wp
	325	IF( zflxv(ji ,jj-1) < 0._wp) zwdlmtv(ji ,jj-1) = 0._wp
	326	ENDIF
	327	END_2D
[9124]	328	!
	329	DO jk1 = 1, nn_wdit + 1 !! start limiter iterations
	330	!
[9019]	331	zflxu1(:,:) = zflxu(:,:) * zwdlmtu(:,:)
	332	zflxv1(:,:) = zflxv(:,:) * zwdlmtv(:,:)
	333	jflag = 0 ! flag indicating if any further iterations are needed
[9124]	334	!
[13295]	335	DO_2D( 0, 1, 0, 1 )
[12377]	336	!
	337	IF( tmask(ji, jj, 1 ) < 0.5_wp ) CYCLE
	338	IF( ht_0(ji,jj) > zdepwd ) CYCLE
	339	!
	340	ztmp = e1e2t(ji,jj)
	341	!
	342	zzflxp = max(zflxu1(ji,jj), 0._wp) - min(zflxu1(ji-1,jj), 0._wp) &
	343	& + max(zflxv1(ji,jj), 0._wp) - min(zflxv1(ji, jj-1), 0._wp)
	344	zzflxn = min(zflxu1(ji,jj), 0._wp) - max(zflxu1(ji-1,jj), 0._wp) &
	345	& + min(zflxv1(ji,jj), 0._wp) - max(zflxv1(ji, jj-1), 0._wp)
	346
	347	zdep1 = (zzflxp + zzflxn) * z2dt / ztmp
	348	zdep2 = ht_0(ji,jj) + sshn_e(ji,jj) - rn_wdmin1 - z2dt * zssh_frc(ji,jj)
	349
	350	IF(zdep1 > zdep2) THEN
	351	zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zflxp(ji,jj) * z2dt )
	352	!zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zzflxp * z2dt )
	353	! flag if the limiter has been used but stop flagging if the only
	354	! changes have zeroed the coefficient since further iterations will
	355	! not change anything
	356	IF( zcoef > 0._wp ) THEN
	357	jflag = 1
	358	ELSE
	359	zcoef = 0._wp
	360	ENDIF
	361	IF(jk1 > nn_wdit) zcoef = 0._wp
	362	IF(zflxu1(ji, jj) > 0._wp) zwdlmtu(ji ,jj) = zcoef
	363	IF(zflxu1(ji-1,jj) < 0._wp) zwdlmtu(ji-1,jj) = zcoef
	364	IF(zflxv1(ji, jj) > 0._wp) zwdlmtv(ji ,jj) = zcoef
	365	IF(zflxv1(ji,jj-1) < 0._wp) zwdlmtv(ji,jj-1) = zcoef
	366	END IF
	367	END_2D
[9124]	368	!
[14644]	369	CALL lbc_lnk( 'wet_dry', zwdlmtu, 'U', 1.0_wp, zwdlmtv, 'V', 1.0_wp )
[9124]	370	!
[10425]	371	CALL mpp_max('wet_dry', jflag) !max over the global domain
[9124]	372	!
	373	IF(jflag == 0) EXIT
	374	!
[9019]	375	END DO ! jk1 loop
[9124]	376	!
[9019]	377	zflxu(:,:) = zflxu(:,:) * zwdlmtu(:, :)
	378	zflxv(:,:) = zflxv(:,:) * zwdlmtv(:, :)
[9124]	379	!
	380	!!gm THIS lbc_lnk is useless since it is already done at the end of the jk1-loop
[14644]	381	CALL lbc_lnk( 'wet_dry', zflxu, 'U', -1.0_wp, zflxv, 'V', -1.0_wp )
[9124]	382	!!gm end
	383	!
	384	IF( jflag == 1 .AND. lwp ) WRITE(numout,*) 'Need more iterations in wad_lmt_bt!!!'
	385	!
[12377]	386	!IF( ln_rnf ) CALL sbc_rnf_div( hdiv ) ! runoffs (update hdiv field)
[9019]	387	!
[9124]	388	IF( ln_timing ) CALL timing_stop('wad_lmt_bt') !
[9019]	389	!
[6152]	390	END SUBROUTINE wad_lmt_bt
[7646]	391
	392	!!==============================================================================
[6152]	393	END MODULE wet_dry

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

Download in other formats: