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.

bdytides.F90 in branches/2013/dev_r3858_NOC_ZTC/NEMOGCM/NEMO/OPA_SRC/BDY – NEMO

Context Navigation

source: branches/2013/dev_r3858_NOC_ZTC/NEMOGCM/NEMO/OPA_SRC/BDY/bdytides.F90 @ 4267

Last change on this file since 4267 was 4263, checked in by acc, 11 years ago
Branch 2013/dev_r3858_NOC_ZTC, #863. Fixes to allow successful compilation of ORCA2_LIM with key_dynspg_ts
Property svn:keywords set to `Id`
File size: 28.7 KB

Rev	Line
[911]	1	MODULE bdytides
[1125]	2	!!======================================================================
[911]	3	!! * MODULE bdytides *
	4	!! Ocean dynamics: Tidal forcing at open boundaries
[1125]	5	!!======================================================================
	6	!! History : 2.0 ! 2007-01 (D.Storkey) Original code
	7	!! 2.3 ! 2008-01 (J.Holt) Add date correction. Origins POLCOMS v6.3 2007
	8	!! 3.0 ! 2008-04 (NEMO team) add in the reference version
[2528]	9	!! 3.3 ! 2010-09 (D.Storkey and E.O'Dea) bug fixes
[3651]	10	!! 3.4 ! 2012-09 (G. Reffray and J. Chanut) New inputs + mods
[4105]	11	!! 3.5 ! 2013-07 (J. Chanut) Compliant with time splitting changes
[1125]	12	!!----------------------------------------------------------------------
	13	#if defined key_bdy
	14	!!----------------------------------------------------------------------
[3294]	15	!! 'key_bdy' Open Boundary Condition
[1125]	16	!!----------------------------------------------------------------------
[911]	17	!! PUBLIC
[3651]	18	!! bdytide_init : read of namelist and initialisation of tidal harmonics data
[1125]	19	!! tide_update : calculation of tidal forcing at each timestep
	20	!!----------------------------------------------------------------------
[3294]	21	USE timing ! Timing
[911]	22	USE oce ! ocean dynamics and tracers
	23	USE dom_oce ! ocean space and time domain
	24	USE iom
	25	USE in_out_manager ! I/O units
	26	USE phycst ! physical constants
	27	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
	28	USE bdy_par ! Unstructured boundary parameters
	29	USE bdy_oce ! ocean open boundary conditions
[2528]	30	USE daymod ! calendar
[3651]	31	USE wrk_nemo ! Memory allocation
	32	USE tideini
	33	! USE tide_mod ! Useless ??
[3294]	34	USE fldread, ONLY: fld_map
[4105]	35	USE dynspg_oce, ONLY: lk_dynspg_ts
[1125]	36
[911]	37	IMPLICIT NONE
	38	PRIVATE
	39
[3651]	40	PUBLIC bdytide_init ! routine called in bdy_init
	41	PUBLIC bdytide_update ! routine called in bdy_dta
[4105]	42	PUBLIC bdy_dta_tides ! routine called in dyn_spg_ts
[911]	43
[3294]	44	TYPE, PUBLIC :: TIDES_DATA !: Storage for external tidal harmonics data
[3651]	45	REAL(wp), POINTER, DIMENSION(:,:,:) :: ssh0 !: Tidal constituents : SSH0 (read in file)
	46	REAL(wp), POINTER, DIMENSION(:,:,:) :: u0 !: Tidal constituents : U0 (read in file)
	47	REAL(wp), POINTER, DIMENSION(:,:,:) :: v0 !: Tidal constituents : V0 (read in file)
	48	REAL(wp), POINTER, DIMENSION(:,:,:) :: ssh !: Tidal constituents : SSH (after nodal cor.)
	49	REAL(wp), POINTER, DIMENSION(:,:,:) :: u !: Tidal constituents : U (after nodal cor.)
	50	REAL(wp), POINTER, DIMENSION(:,:,:) :: v !: Tidal constituents : V (after nodal cor.)
[3294]	51	END TYPE TIDES_DATA
[911]	52
[3651]	53	TYPE(TIDES_DATA), PUBLIC, DIMENSION(jp_bdy), TARGET :: tides !: External tidal harmonics data
[4105]	54	TYPE(OBC_DATA) , PRIVATE, DIMENSION(jp_bdy) :: dta_bdy_s !: bdy external data (slow component)
[911]	55
[1125]	56	!!----------------------------------------------------------------------
[2528]	57	!! NEMO/OPA 3.3 , NEMO Consortium (2010)
[1146]	58	!! $Id$
[2528]	59	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
[1125]	60	!!----------------------------------------------------------------------
[911]	61	CONTAINS
	62
[3651]	63	SUBROUTINE bdytide_init
[1125]	64	!!----------------------------------------------------------------------
[3651]	65	!! * SUBROUTINE bdytide_init *
[1125]	66	!!
[3294]	67	!! ** Purpose : - Read in namelist for tides and initialise external
	68	!! tidal harmonics data
[911]	69	!!
	70	!!----------------------------------------------------------------------
[3294]	71	!! namelist variables
	72	!!-------------------
[3651]	73	CHARACTER(len=80) :: filtide !: Filename root for tidal input files
	74	LOGICAL :: ln_bdytide_2ddta !: If true, read 2d harmonic data
	75	LOGICAL :: ln_bdytide_conj !: If true, assume complex conjugate tidal data
[1125]	76	!!
[3651]	77	INTEGER :: ib_bdy, itide, ib !: dummy loop indices
	78	INTEGER :: ii, ij !: dummy loop indices
[3294]	79	INTEGER :: inum, igrd
[3651]	80	INTEGER, DIMENSION(3) :: ilen0 !: length of boundary data (from OBC arrays)
	81	INTEGER, POINTER, DIMENSION(:) :: nblen, nblenrim ! short cuts
	82	CHARACTER(len=80) :: clfile !: full file name for tidal input file
	83	REAL(wp),ALLOCATABLE, DIMENSION(:,:,:) :: dta_read !: work space to read in tidal harmonics data
	84	REAL(wp), POINTER, DIMENSION(:,:) :: ztr, zti !: " " " " " " " "
[3294]	85	!!
[3651]	86	TYPE(TIDES_DATA), POINTER :: td !: local short cut
[3294]	87	!!
[3651]	88	NAMELIST/nambdy_tide/filtide, ln_bdytide_2ddta, ln_bdytide_conj
[1125]	89	!!----------------------------------------------------------------------
[911]	90
[3651]	91	IF( nn_timing == 1 ) CALL timing_start('bdytide_init')
[3294]	92
[3651]	93	IF (nb_bdy>0) THEN
	94	IF(lwp) WRITE(numout,*)
	95	IF(lwp) WRITE(numout,*) 'bdytide_init : initialization of tidal harmonic forcing at open boundaries'
	96	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~'
	97	ENDIF
[911]	98
[3651]	99	ln_bdytide_2ddta = .FALSE.
	100	ln_bdytide_conj = .FALSE.
	101
[3294]	102	REWIND(numnam)
	103	DO ib_bdy = 1, nb_bdy
	104	IF( nn_dyn2d_dta(ib_bdy) .ge. 2 ) THEN
[2528]	105
[3294]	106	td => tides(ib_bdy)
[3651]	107	nblen => idx_bdy(ib_bdy)%nblen
	108	nblenrim => idx_bdy(ib_bdy)%nblenrim
[911]	109
[3294]	110	! Namelist nambdy_tide : tidal harmonic forcing at open boundaries
	111	filtide(:) = ''
[911]	112
[3294]	113	! Don't REWIND here - may need to read more than one of these namelists.
	114	READ ( numnam, nambdy_tide )
[3651]	115	! ! Parameter control and print
	116	IF(lwp) WRITE(numout,*) ' '
	117	IF(lwp) WRITE(numout,*) ' Namelist nambdy_tide : tidal harmonic forcing at open boundaries'
	118	IF(lwp) WRITE(numout,*) ' read tidal data in 2d files: ', ln_bdytide_2ddta
	119	IF(lwp) WRITE(numout,*) ' assume complex conjugate : ', ln_bdytide_conj
	120	IF(lwp) WRITE(numout,*) ' Number of tidal components to read: ', nb_harmo
	121	IF(lwp) THEN
	122	WRITE(numout,*) ' Tidal cpt name - Phase speed (deg/hr)'
	123	DO itide = 1, nb_harmo
	124	WRITE(numout,*) ' ', Wave(ntide(itide))%cname_tide, omega_tide(itide)
	125	END DO
	126	ENDIF
	127	IF(lwp) WRITE(numout,*) ' '
[1125]	128
[3651]	129	! Allocate space for tidal harmonics data - get size from OBC data arrays
	130	! -----------------------------------------------------------------------
[911]	131
[3651]	132	! JC: If FRS scheme is used, we assume that tidal is needed over the whole
	133	! relaxation area
[4254]	134	IF( cn_dyn2d(ib_bdy) == 'frs' ) THEN
[3651]	135	ilen0(:)=nblen(:)
[3294]	136	ELSE
[3651]	137	ilen0(:)=nblenrim(:)
[3294]	138	ENDIF
[911]	139
[3651]	140	ALLOCATE( td%ssh0( ilen0(1), nb_harmo, 2 ) )
	141	ALLOCATE( td%ssh ( ilen0(1), nb_harmo, 2 ) )
[911]	142
[3651]	143	ALLOCATE( td%u0( ilen0(2), nb_harmo, 2 ) )
	144	ALLOCATE( td%u ( ilen0(2), nb_harmo, 2 ) )
[911]	145
[3651]	146	ALLOCATE( td%v0( ilen0(3), nb_harmo, 2 ) )
	147	ALLOCATE( td%v ( ilen0(3), nb_harmo, 2 ) )
[911]	148
[4105]	149	td%ssh0(:,:,:) = 0._wp
	150	td%ssh (:,:,:) = 0._wp
	151	td%u0 (:,:,:) = 0._wp
	152	td%u (:,:,:) = 0._wp
	153	td%v0 (:,:,:) = 0._wp
	154	td%v (:,:,:) = 0._wp
[3294]	155
[3651]	156	IF (ln_bdytide_2ddta) THEN
	157	! It is assumed that each data file contains all complex harmonic amplitudes
	158	! given on the data domain (ie global, jpidta x jpjdta)
	159	!
	160	CALL wrk_alloc( jpi, jpj, zti, ztr )
	161	!
	162	! SSH fields
	163	clfile = TRIM(filtide)//'_grid_T.nc'
	164	CALL iom_open (clfile , inum )
	165	igrd = 1 ! Everything is at T-points here
	166	DO itide = 1, nb_harmo
	167	CALL iom_get ( inum, jpdom_data, TRIM(Wave(ntide(itide))%cname_tide)//'_z1', ztr(:,:) )
	168	CALL iom_get ( inum, jpdom_data, TRIM(Wave(ntide(itide))%cname_tide)//'_z2', zti(:,:) )
	169	DO ib = 1, ilen0(igrd)
	170	ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
	171	ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
	172	td%ssh0(ib,itide,1) = ztr(ii,ij)
	173	td%ssh0(ib,itide,2) = zti(ii,ij)
	174	END DO
	175	END DO
[3294]	176	CALL iom_close( inum )
[3651]	177	!
	178	! U fields
	179	clfile = TRIM(filtide)//'_grid_U.nc'
	180	CALL iom_open (clfile , inum )
	181	igrd = 2 ! Everything is at U-points here
	182	DO itide = 1, nb_harmo
	183	CALL iom_get ( inum, jpdom_data, TRIM(Wave(ntide(itide))%cname_tide)//'_u1', ztr(:,:) )
	184	CALL iom_get ( inum, jpdom_data, TRIM(Wave(ntide(itide))%cname_tide)//'_u2', zti(:,:) )
	185	DO ib = 1, ilen0(igrd)
	186	ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
	187	ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
	188	td%u0(ib,itide,1) = ztr(ii,ij)
	189	td%u0(ib,itide,2) = zti(ii,ij)
	190	END DO
	191	END DO
[3294]	192	CALL iom_close( inum )
[3651]	193	!
	194	! V fields
	195	clfile = TRIM(filtide)//'_grid_V.nc'
	196	CALL iom_open (clfile , inum )
	197	igrd = 3 ! Everything is at V-points here
	198	DO itide = 1, nb_harmo
	199	CALL iom_get ( inum, jpdom_data, TRIM(Wave(ntide(itide))%cname_tide)//'_v1', ztr(:,:) )
	200	CALL iom_get ( inum, jpdom_data, TRIM(Wave(ntide(itide))%cname_tide)//'_v2', zti(:,:) )
	201	DO ib = 1, ilen0(igrd)
	202	ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
	203	ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
	204	td%v0(ib,itide,1) = ztr(ii,ij)
	205	td%v0(ib,itide,2) = zti(ii,ij)
	206	END DO
	207	END DO
[3294]	208	CALL iom_close( inum )
	209	!
[3651]	210	CALL wrk_dealloc( jpi, jpj, ztr, zti )
	211	!
	212	ELSE
	213	!
	214	! Read tidal data only on bdy segments
	215	!
	216	ALLOCATE( dta_read( MAXVAL(ilen0(1:3)), 1, 1 ) )
[3294]	217
[3651]	218	! Open files and read in tidal forcing data
	219	! -----------------------------------------
[3294]	220
[3651]	221	DO itide = 1, nb_harmo
	222	! ! SSH fields
	223	clfile = TRIM(filtide)//TRIM(Wave(ntide(itide))%cname_tide)//'_grid_T.nc'
	224	CALL iom_open( clfile, inum )
	225	CALL fld_map( inum, 'z1' , dta_read(1:ilen0(1),1:1,1:1) , 1, idx_bdy(ib_bdy)%nbmap(:,1) )
	226	td%ssh0(:,itide,1) = dta_read(1:ilen0(1),1,1)
	227	CALL fld_map( inum, 'z2' , dta_read(1:ilen0(1),1:1,1:1) , 1, idx_bdy(ib_bdy)%nbmap(:,1) )
	228	td%ssh0(:,itide,2) = dta_read(1:ilen0(1),1,1)
	229	CALL iom_close( inum )
	230	! ! U fields
	231	clfile = TRIM(filtide)//TRIM(Wave(ntide(itide))%cname_tide)//'_grid_U.nc'
	232	CALL iom_open( clfile, inum )
	233	CALL fld_map( inum, 'u1' , dta_read(1:ilen0(2),1:1,1:1) , 1, idx_bdy(ib_bdy)%nbmap(:,2) )
	234	td%u0(:,itide,1) = dta_read(1:ilen0(2),1,1)
	235	CALL fld_map( inum, 'u2' , dta_read(1:ilen0(2),1:1,1:1) , 1, idx_bdy(ib_bdy)%nbmap(:,2) )
	236	td%u0(:,itide,2) = dta_read(1:ilen0(2),1,1)
	237	CALL iom_close( inum )
	238	! ! V fields
	239	clfile = TRIM(filtide)//TRIM(Wave(ntide(itide))%cname_tide)//'_grid_V.nc'
	240	CALL iom_open( clfile, inum )
	241	CALL fld_map( inum, 'v1' , dta_read(1:ilen0(3),1:1,1:1) , 1, idx_bdy(ib_bdy)%nbmap(:,3) )
	242	td%v0(:,itide,1) = dta_read(1:ilen0(3),1,1)
	243	CALL fld_map( inum, 'v2' , dta_read(1:ilen0(3),1:1,1:1) , 1, idx_bdy(ib_bdy)%nbmap(:,3) )
	244	td%v0(:,itide,2) = dta_read(1:ilen0(3),1,1)
	245	CALL iom_close( inum )
[3294]	246	!
[3651]	247	END DO ! end loop on tidal components
[3294]	248	!
[3651]	249	DEALLOCATE( dta_read )
	250	ENDIF ! ln_bdytide_2ddta=.true.
[1125]	251	!
[3651]	252	IF ( ln_bdytide_conj ) THEN ! assume complex conjugate in data files
	253	td%ssh0(:,:,2) = - td%ssh0(:,:,2)
	254	td%u0 (:,:,2) = - td%u0 (:,:,2)
	255	td%v0 (:,:,2) = - td%v0 (:,:,2)
	256	ENDIF
	257	!
[4105]	258	IF ( lk_dynspg_ts ) THEN ! Allocate arrays to save slowly varying boundary data during
	259	! time splitting integration
	260	ALLOCATE( dta_bdy_s(ib_bdy)%ssh ( ilen0(1) ) )
	261	ALLOCATE( dta_bdy_s(ib_bdy)%u2d ( ilen0(2) ) )
	262	ALLOCATE( dta_bdy_s(ib_bdy)%v2d ( ilen0(3) ) )
	263	dta_bdy_s(ib_bdy)%ssh(:) = 0.e0
	264	dta_bdy_s(ib_bdy)%u2d(:) = 0.e0
	265	dta_bdy_s(ib_bdy)%v2d(:) = 0.e0
	266	ENDIF
	267	!
[3294]	268	ENDIF ! nn_dyn2d_dta(ib_bdy) .ge. 2
[1125]	269	!
[3294]	270	END DO ! loop on ib_bdy
[911]	271
[3651]	272	IF( nn_timing == 1 ) CALL timing_stop('bdytide_init')
[1125]	273
[3651]	274	END SUBROUTINE bdytide_init
[3294]	275
[3651]	276	SUBROUTINE bdytide_update ( kt, idx, dta, td, jit, time_offset )
[1125]	277	!!----------------------------------------------------------------------
[3651]	278	!! * SUBROUTINE bdytide_update *
[1125]	279	!!
[3294]	280	!! ** Purpose : - Add tidal forcing to ssh, u2d and v2d OBC data arrays.
[911]	281	!!
[1125]	282	!!----------------------------------------------------------------------
[3651]	283	INTEGER, INTENT( in ) :: kt ! Main timestep counter
	284	TYPE(OBC_INDEX), INTENT( in ) :: idx ! OBC indices
	285	TYPE(OBC_DATA), INTENT(inout) :: dta ! OBC external data
	286	TYPE(TIDES_DATA),INTENT( inout ) :: td ! tidal harmonics data
	287	INTEGER,INTENT(in),OPTIONAL :: jit ! Barotropic timestep counter (for timesplitting option)
	288	INTEGER,INTENT( in ), OPTIONAL :: time_offset ! time offset in units of timesteps. NB. if jit
	289	! is present then units = subcycle timesteps.
	290	! time_offset = 0 => get data at "now" time level
	291	! time_offset = -1 => get data at "before" time level
	292	! time_offset = +1 => get data at "after" time level
	293	! etc.
[1125]	294	!!
[3651]	295	INTEGER, DIMENSION(3) :: ilen0 !: length of boundary data (from OBC arrays)
	296	INTEGER :: itide, igrd, ib ! dummy loop indices
	297	INTEGER :: time_add ! time offset in units of timesteps
	298	REAL(wp) :: z_arg, z_sarg, zflag, zramp
	299	REAL(wp), DIMENSION(jpmax_harmo) :: z_sist, z_cost
[1125]	300	!!----------------------------------------------------------------------
[911]	301
[3651]	302	IF( nn_timing == 1 ) CALL timing_start('bdytide_update')
[3294]	303
[3651]	304	ilen0(1) = SIZE(td%ssh(:,1,1))
	305	ilen0(2) = SIZE(td%u(:,1,1))
	306	ilen0(3) = SIZE(td%v(:,1,1))
	307
	308	zflag=1
	309	IF ( PRESENT(jit) ) THEN
	310	IF ( jit /= 1 ) zflag=0
	311	ENDIF
	312
[4105]	313	IF ( nsec_day == NINT(0.5_wp * rdttra(1)) .AND. zflag==1 ) THEN
[3651]	314	!
	315	kt_tide = kt
	316	!
	317	IF(lwp) THEN
	318	WRITE(numout,*)
	319	WRITE(numout,*) 'bdytide_update : (re)Initialization of the tidal bdy forcing at kt=',kt
	320	WRITE(numout,*) '~~~~~~~~~~~~~~ '
	321	ENDIF
	322	!
	323	CALL tide_init_elevation ( idx, td )
	324	CALL tide_init_velocities( idx, td )
	325	!
	326	ENDIF
	327
[3294]	328	time_add = 0
	329	IF( PRESENT(time_offset) ) THEN
	330	time_add = time_offset
	331	ENDIF
	332
	333	IF( PRESENT(jit) ) THEN
[4105]	334	z_arg = ((kt-kt_tide) * rdt + (jit+0.5_wp(time_add-1)) rdt / REAL(nn_baro,wp) )
[3294]	335	ELSE
[3651]	336	z_arg = ((kt-kt_tide)+time_add) * rdt
[1125]	337	ENDIF
[911]	338
[3651]	339	! Linear ramp on tidal component at open boundaries
[4105]	340	zramp = 1._wp
	341	IF (ln_tide_ramp) zramp = MIN(MAX( (z_arg + (kt_tide-nit000)rdt)/(rdttideramprday),0._wp),1._wp)
[3651]	342
	343	DO itide = 1, nb_harmo
	344	z_sarg = z_arg * omega_tide(itide)
[1125]	345	z_cost(itide) = COS( z_sarg )
	346	z_sist(itide) = SIN( z_sarg )
	347	END DO
[911]	348
[3651]	349	DO itide = 1, nb_harmo
	350	igrd=1 ! SSH on tracer grid
	351	DO ib = 1, ilen0(igrd)
	352	dta%ssh(ib) = dta%ssh(ib) + zramp(td%ssh(ib,itide,1)z_cost(itide) + td%ssh(ib,itide,2)*z_sist(itide))
[1125]	353	END DO
[3294]	354	igrd=2 ! U grid
[3651]	355	DO ib = 1, ilen0(igrd)
	356	dta%u2d(ib) = dta%u2d(ib) + zramp(td%u (ib,itide,1)z_cost(itide) + td%u (ib,itide,2)*z_sist(itide))
[1125]	357	END DO
[3294]	358	igrd=3 ! V grid
[3651]	359	DO ib = 1, ilen0(igrd)
	360	dta%v2d(ib) = dta%v2d(ib) + zramp(td%v (ib,itide,1)z_cost(itide) + td%v (ib,itide,2)*z_sist(itide))
[1125]	361	END DO
	362	END DO
	363	!
[3651]	364	IF( nn_timing == 1 ) CALL timing_stop('bdytide_update')
[3294]	365	!
[3651]	366	END SUBROUTINE bdytide_update
[911]	367
[4105]	368	SUBROUTINE bdy_dta_tides( kt, kit, time_offset )
	369	!!----------------------------------------------------------------------
	370	!! * SUBROUTINE bdy_dta_tides *
	371	!!
	372	!! ** Purpose : - Add tidal forcing to ssh, u2d and v2d OBC data arrays.
	373	!!
	374	!!----------------------------------------------------------------------
	375	INTEGER, INTENT( in ) :: kt ! Main timestep counter
	376	INTEGER, INTENT( in ),OPTIONAL :: kit ! Barotropic timestep counter (for timesplitting option)
	377	INTEGER, INTENT( in ),OPTIONAL :: time_offset ! time offset in units of timesteps. NB. if kit
	378	! is present then units = subcycle timesteps.
	379	! time_offset = 0 => get data at "now" time level
	380	! time_offset = -1 => get data at "before" time level
	381	! time_offset = +1 => get data at "after" time level
	382	! etc.
	383	!!
	384	LOGICAL :: lk_first_btstp ! =.TRUE. if time splitting and first barotropic step
	385	INTEGER, DIMENSION(jpbgrd) :: ilen0
	386	INTEGER, DIMENSION(1:jpbgrd) :: nblen, nblenrim ! short cuts
	387	INTEGER :: itide, ib_bdy, ib, igrd ! loop indices
	388	INTEGER :: time_add ! time offset in units of timesteps
	389	REAL(wp) :: z_arg, z_sarg, zramp, zoff, z_cost, z_sist
	390	!!----------------------------------------------------------------------
	391
	392	IF( nn_timing == 1 ) CALL timing_start('bdy_dta_tides')
	393
	394	lk_first_btstp=.TRUE.
	395	IF ( PRESENT(kit).AND.( kit /= 1 ) ) THEN ; lk_first_btstp=.FALSE. ; ENDIF
	396
	397	time_add = 0
	398	IF( PRESENT(time_offset) ) THEN
	399	time_add = time_offset
	400	ENDIF
	401
	402	! Absolute time from model initialization:
	403	IF( PRESENT(kit) ) THEN
	404	z_arg = ( kt + (kit+0.5_wp(time_add-1)) / REAL(nn_baro,wp) ) rdt
	405	ELSE
	406	z_arg = ( kt + time_add ) * rdt
	407	ENDIF
	408
	409	! Linear ramp on tidal component at open boundaries
	410	zramp = 1.
	411	IF (ln_tide_ramp) zramp = MIN(MAX( (z_arg - nit000rdt)/(rdttideramprday),0.),1.)
	412
	413	DO ib_bdy = 1,nb_bdy
	414
	415	! line below should be simplified (runoff case)
[4254]	416	!! CHANUT: TO BE SORTED OUT
	417	!! IF (( nn_dyn2d_dta(ib_bdy) .ge. 2 ).AND.(nn_tra(ib_bdy).NE.4)) THEN
	418	IF ( nn_dyn2d_dta(ib_bdy) .ge. 2 ) THEN
[4105]	419
	420	nblen(1:jpbgrd) = idx_bdy(ib_bdy)%nblen(1:jpbgrd)
	421	nblenrim(1:jpbgrd) = idx_bdy(ib_bdy)%nblenrim(1:jpbgrd)
	422
[4254]	423	IF( cn_dyn2d(ib_bdy) == 'frs' ) THEN
[4105]	424	ilen0(:)=nblen(:)
	425	ELSE
	426	ilen0(:)=nblenrim(:)
	427	ENDIF
	428
	429	! We refresh nodal factors every day below
	430	! This should be done somewhere else
	431	IF ( nsec_day == NINT(0.5_wp * rdttra(1)) .AND. lk_first_btstp ) THEN
	432	!
	433	kt_tide = kt
	434	!
	435	IF(lwp) THEN
	436	WRITE(numout,*)
	437	WRITE(numout,*) 'bdy_tide_dta : Refresh nodal factors for tidal open bdy data at kt=',kt
	438	WRITE(numout,*) '~~~~~~~~~~~~~~ '
	439	ENDIF
	440	!
	441	CALL tide_init_elevation ( idx=idx_bdy(ib_bdy), td=tides(ib_bdy) )
	442	CALL tide_init_velocities( idx=idx_bdy(ib_bdy), td=tides(ib_bdy) )
	443	!
	444	ENDIF
	445	zoff = -kt_tide * rdt ! time offset relative to nodal factor computation time
	446	!
	447	! If time splitting, save data at first barotropic iteration
	448	IF ( PRESENT(kit) ) THEN
	449	IF ( lk_first_btstp ) THEN ! Save slow varying open boundary data:
	450	dta_bdy_s(ib_bdy)%ssh(1:ilen0(1)) = dta_bdy(ib_bdy)%ssh(1:ilen0(1))
	451	dta_bdy_s(ib_bdy)%u2d(1:ilen0(2)) = dta_bdy(ib_bdy)%u2d(1:ilen0(2))
	452	dta_bdy_s(ib_bdy)%v2d(1:ilen0(3)) = dta_bdy(ib_bdy)%v2d(1:ilen0(3))
	453
	454	ELSE ! Initialize arrays from slow varying open boundary data:
	455	dta_bdy(ib_bdy)%ssh(1:ilen0(1)) = dta_bdy_s(ib_bdy)%ssh(1:ilen0(1))
	456	dta_bdy(ib_bdy)%u2d(1:ilen0(2)) = dta_bdy_s(ib_bdy)%u2d(1:ilen0(2))
	457	dta_bdy(ib_bdy)%v2d(1:ilen0(3)) = dta_bdy_s(ib_bdy)%v2d(1:ilen0(3))
	458	ENDIF
	459	ENDIF
	460	!
	461	! Update open boundary data arrays:
	462	DO itide = 1, nb_harmo
	463	!
	464	z_sarg = (z_arg + zoff) * omega_tide(itide)
	465	z_cost = zramp * COS( z_sarg )
	466	z_sist = zramp * SIN( z_sarg )
	467	!
	468	igrd=1 ! SSH on tracer grid
	469	DO ib = 1, ilen0(igrd)
	470	dta_bdy(ib_bdy)%ssh(ib) = dta_bdy(ib_bdy)%ssh(ib) + &
	471	& ( tides(ib_bdy)%ssh(ib,itide,1)*z_cost + &
	472	& tides(ib_bdy)%ssh(ib,itide,2)*z_sist )
	473	END DO
	474	!
	475	igrd=2 ! U grid
	476	DO ib = 1, ilen0(igrd)
	477	dta_bdy(ib_bdy)%u2d(ib) = dta_bdy(ib_bdy)%u2d(ib) + &
	478	& ( tides(ib_bdy)%u(ib,itide,1)*z_cost + &
	479	& tides(ib_bdy)%u(ib,itide,2)*z_sist )
	480	END DO
	481	!
	482	igrd=3 ! V grid
	483	DO ib = 1, ilen0(igrd)
	484	dta_bdy(ib_bdy)%v2d(ib) = dta_bdy(ib_bdy)%v2d(ib) + &
	485	& ( tides(ib_bdy)%v(ib,itide,1)*z_cost + &
	486	& tides(ib_bdy)%v(ib,itide,2)*z_sist )
	487	END DO
	488	END DO
	489	END IF
	490	END DO
	491	!
	492	IF( nn_timing == 1 ) CALL timing_stop('bdy_dta_tides')
	493	!
	494	END SUBROUTINE bdy_dta_tides
	495
[3651]	496	SUBROUTINE tide_init_elevation( idx, td )
[1125]	497	!!----------------------------------------------------------------------
[3651]	498	!! * ROUTINE tide_init_elevation *
[1125]	499	!!----------------------------------------------------------------------
[3651]	500	TYPE(OBC_INDEX), INTENT( in ) :: idx ! OBC indices
	501	TYPE(TIDES_DATA),INTENT( inout ) :: td ! tidal harmonics data
	502	!! * Local declarations
	503	INTEGER, DIMENSION(1) :: ilen0 !: length of boundary data (from OBC arrays)
	504	REAL(wp),ALLOCATABLE, DIMENSION(:) :: mod_tide, phi_tide
	505	INTEGER :: itide, igrd, ib ! dummy loop indices
[1125]	506
[3651]	507	igrd=1
	508	! SSH on tracer grid.
	509
	510	ilen0(1) = SIZE(td%ssh0(:,1,1))
	511
	512	ALLOCATE(mod_tide(ilen0(igrd)),phi_tide(ilen0(igrd)))
	513
	514	DO itide = 1, nb_harmo
	515	DO ib = 1, ilen0(igrd)
	516	mod_tide(ib)=SQRT(td%ssh0(ib,itide,1)2.+td%ssh0(ib,itide,2)2.)
	517	phi_tide(ib)=ATAN2(-td%ssh0(ib,itide,2),td%ssh0(ib,itide,1))
[1125]	518	END DO
[3651]	519	DO ib = 1 , ilen0(igrd)
	520	mod_tide(ib)=mod_tide(ib)*ftide(itide)
	521	phi_tide(ib)=phi_tide(ib)+v0tide(itide)+utide(itide)
	522	ENDDO
	523	DO ib = 1 , ilen0(igrd)
	524	td%ssh(ib,itide,1)= mod_tide(ib)*COS(phi_tide(ib))
	525	td%ssh(ib,itide,2)=-mod_tide(ib)*SIN(phi_tide(ib))
	526	ENDDO
	527	END DO
[911]	528
[3651]	529	DEALLOCATE(mod_tide,phi_tide)
[911]	530
[3651]	531	END SUBROUTINE tide_init_elevation
	532
	533	SUBROUTINE tide_init_velocities( idx, td )
[1125]	534	!!----------------------------------------------------------------------
[3651]	535	!! * ROUTINE tide_init_elevation *
[1125]	536	!!----------------------------------------------------------------------
[3651]	537	TYPE(OBC_INDEX), INTENT( in ) :: idx ! OBC indices
	538	TYPE(TIDES_DATA),INTENT( inout ) :: td ! tidal harmonics data
	539	!! * Local declarations
	540	INTEGER, DIMENSION(3) :: ilen0 !: length of boundary data (from OBC arrays)
	541	REAL(wp),ALLOCATABLE, DIMENSION(:) :: mod_tide, phi_tide
	542	INTEGER :: itide, igrd, ib ! dummy loop indices
[911]	543
[3651]	544	ilen0(2) = SIZE(td%u0(:,1,1))
	545	ilen0(3) = SIZE(td%v0(:,1,1))
[1125]	546
[3651]	547	igrd=2 ! U grid.
[911]	548
[3651]	549	ALLOCATE(mod_tide(ilen0(igrd)),phi_tide(ilen0(igrd)))
[911]	550
[3651]	551	DO itide = 1, nb_harmo
	552	DO ib = 1, ilen0(igrd)
	553	mod_tide(ib)=SQRT(td%u0(ib,itide,1)2.+td%u0(ib,itide,2)2.)
	554	phi_tide(ib)=ATAN2(-td%u0(ib,itide,2),td%u0(ib,itide,1))
	555	END DO
	556	DO ib = 1, ilen0(igrd)
	557	mod_tide(ib)=mod_tide(ib)*ftide(itide)
	558	phi_tide(ib)=phi_tide(ib)+v0tide(itide)+utide(itide)
	559	ENDDO
	560	DO ib = 1, ilen0(igrd)
	561	td%u(ib,itide,1)= mod_tide(ib)*COS(phi_tide(ib))
	562	td%u(ib,itide,2)=-mod_tide(ib)*SIN(phi_tide(ib))
	563	ENDDO
	564	END DO
[911]	565
[3651]	566	DEALLOCATE(mod_tide,phi_tide)
[911]	567
[3651]	568	igrd=3 ! V grid.
[911]	569
[3651]	570	ALLOCATE(mod_tide(ilen0(igrd)),phi_tide(ilen0(igrd)))
[911]	571
[3651]	572	DO itide = 1, nb_harmo
	573	DO ib = 1, ilen0(igrd)
	574	mod_tide(ib)=SQRT(td%v0(ib,itide,1)2.+td%v0(ib,itide,2)2.)
	575	phi_tide(ib)=ATAN2(-td%v0(ib,itide,2),td%v0(ib,itide,1))
	576	END DO
	577	DO ib = 1, ilen0(igrd)
	578	mod_tide(ib)=mod_tide(ib)*ftide(itide)
	579	phi_tide(ib)=phi_tide(ib)+v0tide(itide)+utide(itide)
	580	ENDDO
	581	DO ib = 1, ilen0(igrd)
	582	td%v(ib,itide,1)= mod_tide(ib)*COS(phi_tide(ib))
	583	td%v(ib,itide,2)=-mod_tide(ib)*SIN(phi_tide(ib))
	584	ENDDO
	585	END DO
[1125]	586
[3651]	587	DEALLOCATE(mod_tide,phi_tide)
[911]	588
[3651]	589	END SUBROUTINE tide_init_velocities
[911]	590	#else
[1125]	591	!!----------------------------------------------------------------------
	592	!! Dummy module NO Unstruct Open Boundary Conditions for tides
	593	!!----------------------------------------------------------------------
[911]	594	CONTAINS
[3651]	595	SUBROUTINE bdytide_init ! Empty routine
	596	WRITE(,) 'bdytide_init: You should not have seen this print! error?'
	597	END SUBROUTINE bdytide_init
	598	SUBROUTINE bdytide_update( kt, jit ) ! Empty routine
	599	WRITE(,) 'bdytide_update: You should not have seen this print! error?', kt, jit
	600	END SUBROUTINE bdytide_update
[4263]	601	SUBROUTINE bdy_dta_tides( kt, kit, time_offset ) ! Empty routine
	602	INTEGER, INTENT( in ) :: kt ! Dummy argument empty routine
	603	INTEGER, INTENT( in ),OPTIONAL :: kit ! Dummy argument empty routine
	604	INTEGER, INTENT( in ),OPTIONAL :: time_offset ! Dummy argument empty routine
[4105]	605	WRITE(,) 'bdy_dta_tides: You should not have seen this print! error?', kt, jit
	606	END SUBROUTINE bdy_dta_tides
[911]	607	#endif
	608
[1125]	609	!!======================================================================
[911]	610	END MODULE bdytides
[4254]	611

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

Download in other formats: