Context Navigation

bdytides.F90 @ 15814

Last change on this file since 15814 was 11101, checked in by frrh, 5 years ago

Merge changes from Met Office GMED ticket 450 to reduce unnecessary
text output from NEMO.
This output, which is typically not switchable, is rarely of interest
in normal (non-debugging) runs and simply redunantley consumes extra
file space.
Further, the presence of this text output has been shown to
significantly degrade performance of models which are run during
Met Office HPC RAID (disk) checks.
The new code introduces switches which are configurable via the
changes made in the associated Met Office MOCI ticket 399.

File size: 30.1 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
[4292]	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 ??
[5132]	34	USE fldread
[4292]	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
[4292]	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
[4354]	53	!$AGRIF_DO_NOT_TREAT
[3651]	54	TYPE(TIDES_DATA), PUBLIC, DIMENSION(jp_bdy), TARGET :: tides !: External tidal harmonics data
[4354]	55	!$AGRIF_END_DO_NOT_TREAT
[4292]	56	TYPE(OBC_DATA) , PRIVATE, DIMENSION(jp_bdy) :: dta_bdy_s !: bdy external data (slow component)
[911]	57
[1125]	58	!!----------------------------------------------------------------------
[2528]	59	!! NEMO/OPA 3.3 , NEMO Consortium (2010)
[6486]	60	!! $Id$
[2528]	61	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
[1125]	62	!!----------------------------------------------------------------------
[911]	63	CONTAINS
	64
[3651]	65	SUBROUTINE bdytide_init
[1125]	66	!!----------------------------------------------------------------------
[3651]	67	!! * SUBROUTINE bdytide_init *
[1125]	68	!!
[3294]	69	!! ** Purpose : - Read in namelist for tides and initialise external
	70	!! tidal harmonics data
[911]	71	!!
	72	!!----------------------------------------------------------------------
[3294]	73	!! namelist variables
	74	!!-------------------
[3651]	75	CHARACTER(len=80) :: filtide !: Filename root for tidal input files
	76	LOGICAL :: ln_bdytide_2ddta !: If true, read 2d harmonic data
	77	LOGICAL :: ln_bdytide_conj !: If true, assume complex conjugate tidal data
[1125]	78	!!
[3651]	79	INTEGER :: ib_bdy, itide, ib !: dummy loop indices
	80	INTEGER :: ii, ij !: dummy loop indices
[3294]	81	INTEGER :: inum, igrd
[3651]	82	INTEGER, DIMENSION(3) :: ilen0 !: length of boundary data (from OBC arrays)
	83	INTEGER, POINTER, DIMENSION(:) :: nblen, nblenrim ! short cuts
[4147]	84	INTEGER :: ios ! Local integer output status for namelist read
[3651]	85	CHARACTER(len=80) :: clfile !: full file name for tidal input file
	86	REAL(wp),ALLOCATABLE, DIMENSION(:,:,:) :: dta_read !: work space to read in tidal harmonics data
	87	REAL(wp), POINTER, DIMENSION(:,:) :: ztr, zti !: " " " " " " " "
[3294]	88	!!
[3651]	89	TYPE(TIDES_DATA), POINTER :: td !: local short cut
[5132]	90	TYPE(MAP_POINTER), DIMENSION(jpbgrd) :: ibmap_ptr !: array of pointers to nbmap
[3294]	91	!!
[3651]	92	NAMELIST/nambdy_tide/filtide, ln_bdytide_2ddta, ln_bdytide_conj
[1125]	93	!!----------------------------------------------------------------------
[911]	94
[3651]	95	IF( nn_timing == 1 ) CALL timing_start('bdytide_init')
[3294]	96
[3651]	97	IF (nb_bdy>0) THEN
	98	IF(lwp) WRITE(numout,*)
	99	IF(lwp) WRITE(numout,*) 'bdytide_init : initialization of tidal harmonic forcing at open boundaries'
	100	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~'
[11101]	101	IF(lflush) CALL flush(numout)
[3651]	102	ENDIF
[911]	103
[4147]	104	REWIND(numnam_cfg)
[3651]	105
[3294]	106	DO ib_bdy = 1, nb_bdy
	107	IF( nn_dyn2d_dta(ib_bdy) .ge. 2 ) THEN
[2528]	108
[3294]	109	td => tides(ib_bdy)
[3651]	110	nblen => idx_bdy(ib_bdy)%nblen
	111	nblenrim => idx_bdy(ib_bdy)%nblenrim
[911]	112
[3294]	113	! Namelist nambdy_tide : tidal harmonic forcing at open boundaries
	114	filtide(:) = ''
[911]	115
[4147]	116	! Don't REWIND here - may need to read more than one of these namelists.
	117	READ ( numnam_ref, nambdy_tide, IOSTAT = ios, ERR = 901)
	118	901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'nambdy_tide in reference namelist', lwp )
	119	READ ( numnam_cfg, nambdy_tide, IOSTAT = ios, ERR = 902 )
	120	902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'nambdy_tide in configuration namelist', lwp )
[11101]	121	IF(lwm .AND. nprint > 2) WRITE ( numond, nambdy_tide )
[3651]	122	! ! Parameter control and print
[11101]	123	IF(lwp) THEN
	124	WRITE(numout,*) ' '
	125	WRITE(numout,*) ' Namelist nambdy_tide : tidal harmonic forcing at open boundaries'
	126	WRITE(numout,*) ' read tidal data in 2d files: ', ln_bdytide_2ddta
	127	WRITE(numout,*) ' assume complex conjugate : ', ln_bdytide_conj
	128	WRITE(numout,*) ' Number of tidal components to read: ', nb_harmo
	129	WRITE(numout,*) ' Tidal components: '
[3651]	130	DO itide = 1, nb_harmo
[5084]	131	WRITE(numout,*) ' ', Wave(ntide(itide))%cname_tide
[3651]	132	END DO
[11101]	133	WRITE(numout,*) ' '
	134	IF(lflush) CALL flush(numout)
	135	ENDIF
[1125]	136
[3651]	137	! Allocate space for tidal harmonics data - get size from OBC data arrays
	138	! -----------------------------------------------------------------------
[911]	139
[3651]	140	! JC: If FRS scheme is used, we assume that tidal is needed over the whole
	141	! relaxation area
[4292]	142	IF( cn_dyn2d(ib_bdy) == 'frs' ) THEN
[3651]	143	ilen0(:)=nblen(:)
[3294]	144	ELSE
[3651]	145	ilen0(:)=nblenrim(:)
[3294]	146	ENDIF
[911]	147
[3651]	148	ALLOCATE( td%ssh0( ilen0(1), nb_harmo, 2 ) )
	149	ALLOCATE( td%ssh ( ilen0(1), nb_harmo, 2 ) )
[911]	150
[3651]	151	ALLOCATE( td%u0( ilen0(2), nb_harmo, 2 ) )
	152	ALLOCATE( td%u ( ilen0(2), nb_harmo, 2 ) )
[911]	153
[3651]	154	ALLOCATE( td%v0( ilen0(3), nb_harmo, 2 ) )
	155	ALLOCATE( td%v ( ilen0(3), nb_harmo, 2 ) )
[911]	156
[4292]	157	td%ssh0(:,:,:) = 0._wp
	158	td%ssh (:,:,:) = 0._wp
	159	td%u0 (:,:,:) = 0._wp
	160	td%u (:,:,:) = 0._wp
	161	td%v0 (:,:,:) = 0._wp
	162	td%v (:,:,:) = 0._wp
[3294]	163
[3651]	164	IF (ln_bdytide_2ddta) THEN
	165	! It is assumed that each data file contains all complex harmonic amplitudes
	166	! given on the data domain (ie global, jpidta x jpjdta)
	167	!
	168	CALL wrk_alloc( jpi, jpj, zti, ztr )
	169	!
	170	! SSH fields
	171	clfile = TRIM(filtide)//'_grid_T.nc'
	172	CALL iom_open (clfile , inum )
	173	igrd = 1 ! Everything is at T-points here
	174	DO itide = 1, nb_harmo
	175	CALL iom_get ( inum, jpdom_data, TRIM(Wave(ntide(itide))%cname_tide)//'_z1', ztr(:,:) )
	176	CALL iom_get ( inum, jpdom_data, TRIM(Wave(ntide(itide))%cname_tide)//'_z2', zti(:,:) )
	177	DO ib = 1, ilen0(igrd)
	178	ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
	179	ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
	180	td%ssh0(ib,itide,1) = ztr(ii,ij)
	181	td%ssh0(ib,itide,2) = zti(ii,ij)
	182	END DO
	183	END DO
[3294]	184	CALL iom_close( inum )
[3651]	185	!
	186	! U fields
	187	clfile = TRIM(filtide)//'_grid_U.nc'
	188	CALL iom_open (clfile , inum )
	189	igrd = 2 ! Everything is at U-points here
	190	DO itide = 1, nb_harmo
	191	CALL iom_get ( inum, jpdom_data, TRIM(Wave(ntide(itide))%cname_tide)//'_u1', ztr(:,:) )
	192	CALL iom_get ( inum, jpdom_data, TRIM(Wave(ntide(itide))%cname_tide)//'_u2', zti(:,:) )
	193	DO ib = 1, ilen0(igrd)
	194	ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
	195	ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
	196	td%u0(ib,itide,1) = ztr(ii,ij)
	197	td%u0(ib,itide,2) = zti(ii,ij)
	198	END DO
	199	END DO
[3294]	200	CALL iom_close( inum )
[3651]	201	!
	202	! V fields
	203	clfile = TRIM(filtide)//'_grid_V.nc'
	204	CALL iom_open (clfile , inum )
	205	igrd = 3 ! Everything is at V-points here
	206	DO itide = 1, nb_harmo
	207	CALL iom_get ( inum, jpdom_data, TRIM(Wave(ntide(itide))%cname_tide)//'_v1', ztr(:,:) )
	208	CALL iom_get ( inum, jpdom_data, TRIM(Wave(ntide(itide))%cname_tide)//'_v2', zti(:,:) )
	209	DO ib = 1, ilen0(igrd)
	210	ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
	211	ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
	212	td%v0(ib,itide,1) = ztr(ii,ij)
	213	td%v0(ib,itide,2) = zti(ii,ij)
	214	END DO
	215	END DO
[3294]	216	CALL iom_close( inum )
	217	!
[3651]	218	CALL wrk_dealloc( jpi, jpj, ztr, zti )
	219	!
	220	ELSE
	221	!
	222	! Read tidal data only on bdy segments
	223	!
	224	ALLOCATE( dta_read( MAXVAL(ilen0(1:3)), 1, 1 ) )
[5132]	225	!
	226	! Set map structure
	227	ibmap_ptr(1)%ptr => idx_bdy(ib_bdy)%nbmap(:,1)
	228	ibmap_ptr(1)%ll_unstruc = ln_coords_file(ib_bdy)
	229	ibmap_ptr(2)%ptr => idx_bdy(ib_bdy)%nbmap(:,2)
	230	ibmap_ptr(2)%ll_unstruc = ln_coords_file(ib_bdy)
	231	ibmap_ptr(3)%ptr => idx_bdy(ib_bdy)%nbmap(:,3)
	232	ibmap_ptr(3)%ll_unstruc = ln_coords_file(ib_bdy)
[3294]	233
[3651]	234	! Open files and read in tidal forcing data
	235	! -----------------------------------------
[3294]	236
[3651]	237	DO itide = 1, nb_harmo
	238	! ! SSH fields
	239	clfile = TRIM(filtide)//TRIM(Wave(ntide(itide))%cname_tide)//'_grid_T.nc'
	240	CALL iom_open( clfile, inum )
[5132]	241	CALL fld_map( inum, 'z1' , dta_read(1:ilen0(1),1:1,1:1) , 1, ibmap_ptr(1) )
[3651]	242	td%ssh0(:,itide,1) = dta_read(1:ilen0(1),1,1)
[5132]	243	CALL fld_map( inum, 'z2' , dta_read(1:ilen0(1),1:1,1:1) , 1, ibmap_ptr(1) )
[3651]	244	td%ssh0(:,itide,2) = dta_read(1:ilen0(1),1,1)
	245	CALL iom_close( inum )
	246	! ! U fields
	247	clfile = TRIM(filtide)//TRIM(Wave(ntide(itide))%cname_tide)//'_grid_U.nc'
	248	CALL iom_open( clfile, inum )
[5132]	249	CALL fld_map( inum, 'u1' , dta_read(1:ilen0(2),1:1,1:1) , 1, ibmap_ptr(2) )
[3651]	250	td%u0(:,itide,1) = dta_read(1:ilen0(2),1,1)
[5132]	251	CALL fld_map( inum, 'u2' , dta_read(1:ilen0(2),1:1,1:1) , 1, ibmap_ptr(2) )
[3651]	252	td%u0(:,itide,2) = dta_read(1:ilen0(2),1,1)
	253	CALL iom_close( inum )
	254	! ! V fields
	255	clfile = TRIM(filtide)//TRIM(Wave(ntide(itide))%cname_tide)//'_grid_V.nc'
	256	CALL iom_open( clfile, inum )
[5132]	257	CALL fld_map( inum, 'v1' , dta_read(1:ilen0(3),1:1,1:1) , 1, ibmap_ptr(3) )
[3651]	258	td%v0(:,itide,1) = dta_read(1:ilen0(3),1,1)
[5132]	259	CALL fld_map( inum, 'v2' , dta_read(1:ilen0(3),1:1,1:1) , 1, ibmap_ptr(3) )
[3651]	260	td%v0(:,itide,2) = dta_read(1:ilen0(3),1,1)
	261	CALL iom_close( inum )
[3294]	262	!
[3651]	263	END DO ! end loop on tidal components
[3294]	264	!
[3651]	265	DEALLOCATE( dta_read )
	266	ENDIF ! ln_bdytide_2ddta=.true.
[1125]	267	!
[3651]	268	IF ( ln_bdytide_conj ) THEN ! assume complex conjugate in data files
	269	td%ssh0(:,:,2) = - td%ssh0(:,:,2)
	270	td%u0 (:,:,2) = - td%u0 (:,:,2)
	271	td%v0 (:,:,2) = - td%v0 (:,:,2)
	272	ENDIF
	273	!
[4292]	274	IF ( lk_dynspg_ts ) THEN ! Allocate arrays to save slowly varying boundary data during
	275	! time splitting integration
	276	ALLOCATE( dta_bdy_s(ib_bdy)%ssh ( ilen0(1) ) )
	277	ALLOCATE( dta_bdy_s(ib_bdy)%u2d ( ilen0(2) ) )
	278	ALLOCATE( dta_bdy_s(ib_bdy)%v2d ( ilen0(3) ) )
	279	dta_bdy_s(ib_bdy)%ssh(:) = 0.e0
	280	dta_bdy_s(ib_bdy)%u2d(:) = 0.e0
	281	dta_bdy_s(ib_bdy)%v2d(:) = 0.e0
	282	ENDIF
	283	!
[3294]	284	ENDIF ! nn_dyn2d_dta(ib_bdy) .ge. 2
[1125]	285	!
[3294]	286	END DO ! loop on ib_bdy
[911]	287
[3651]	288	IF( nn_timing == 1 ) CALL timing_stop('bdytide_init')
[1125]	289
[3651]	290	END SUBROUTINE bdytide_init
[3294]	291
[3651]	292	SUBROUTINE bdytide_update ( kt, idx, dta, td, jit, time_offset )
[1125]	293	!!----------------------------------------------------------------------
[3651]	294	!! * SUBROUTINE bdytide_update *
[1125]	295	!!
[3294]	296	!! ** Purpose : - Add tidal forcing to ssh, u2d and v2d OBC data arrays.
[911]	297	!!
[1125]	298	!!----------------------------------------------------------------------
[3651]	299	INTEGER, INTENT( in ) :: kt ! Main timestep counter
	300	TYPE(OBC_INDEX), INTENT( in ) :: idx ! OBC indices
	301	TYPE(OBC_DATA), INTENT(inout) :: dta ! OBC external data
	302	TYPE(TIDES_DATA),INTENT( inout ) :: td ! tidal harmonics data
	303	INTEGER,INTENT(in),OPTIONAL :: jit ! Barotropic timestep counter (for timesplitting option)
	304	INTEGER,INTENT( in ), OPTIONAL :: time_offset ! time offset in units of timesteps. NB. if jit
	305	! is present then units = subcycle timesteps.
	306	! time_offset = 0 => get data at "now" time level
	307	! time_offset = -1 => get data at "before" time level
	308	! time_offset = +1 => get data at "after" time level
	309	! etc.
[1125]	310	!!
[3651]	311	INTEGER, DIMENSION(3) :: ilen0 !: length of boundary data (from OBC arrays)
	312	INTEGER :: itide, igrd, ib ! dummy loop indices
	313	INTEGER :: time_add ! time offset in units of timesteps
	314	REAL(wp) :: z_arg, z_sarg, zflag, zramp
	315	REAL(wp), DIMENSION(jpmax_harmo) :: z_sist, z_cost
[1125]	316	!!----------------------------------------------------------------------
[911]	317
[3651]	318	IF( nn_timing == 1 ) CALL timing_start('bdytide_update')
[3294]	319
[3651]	320	ilen0(1) = SIZE(td%ssh(:,1,1))
	321	ilen0(2) = SIZE(td%u(:,1,1))
	322	ilen0(3) = SIZE(td%v(:,1,1))
	323
	324	zflag=1
	325	IF ( PRESENT(jit) ) THEN
	326	IF ( jit /= 1 ) zflag=0
	327	ENDIF
	328
[4292]	329	IF ( nsec_day == NINT(0.5_wp * rdttra(1)) .AND. zflag==1 ) THEN
[3651]	330	!
	331	kt_tide = kt
	332	!
	333	IF(lwp) THEN
	334	WRITE(numout,*)
	335	WRITE(numout,*) 'bdytide_update : (re)Initialization of the tidal bdy forcing at kt=',kt
	336	WRITE(numout,*) '~~~~~~~~~~~~~~ '
[11101]	337	IF(lflush) CALL flush(numout)
[3651]	338	ENDIF
	339	!
	340	CALL tide_init_elevation ( idx, td )
	341	CALL tide_init_velocities( idx, td )
	342	!
	343	ENDIF
	344
[3294]	345	time_add = 0
	346	IF( PRESENT(time_offset) ) THEN
	347	time_add = time_offset
	348	ENDIF
	349
	350	IF( PRESENT(jit) ) THEN
[4292]	351	z_arg = ((kt-kt_tide) * rdt + (jit+0.5_wp(time_add-1)) rdt / REAL(nn_baro,wp) )
[3294]	352	ELSE
[3651]	353	z_arg = ((kt-kt_tide)+time_add) * rdt
[1125]	354	ENDIF
[911]	355
[3651]	356	! Linear ramp on tidal component at open boundaries
[4292]	357	zramp = 1._wp
	358	IF (ln_tide_ramp) zramp = MIN(MAX( (z_arg + (kt_tide-nit000)rdt)/(rdttideramprday),0._wp),1._wp)
[3651]	359
	360	DO itide = 1, nb_harmo
	361	z_sarg = z_arg * omega_tide(itide)
[1125]	362	z_cost(itide) = COS( z_sarg )
	363	z_sist(itide) = SIN( z_sarg )
	364	END DO
[911]	365
[3651]	366	DO itide = 1, nb_harmo
	367	igrd=1 ! SSH on tracer grid
	368	DO ib = 1, ilen0(igrd)
	369	dta%ssh(ib) = dta%ssh(ib) + zramp(td%ssh(ib,itide,1)z_cost(itide) + td%ssh(ib,itide,2)*z_sist(itide))
[1125]	370	END DO
[3294]	371	igrd=2 ! U grid
[3651]	372	DO ib = 1, ilen0(igrd)
	373	dta%u2d(ib) = dta%u2d(ib) + zramp(td%u (ib,itide,1)z_cost(itide) + td%u (ib,itide,2)*z_sist(itide))
[1125]	374	END DO
[3294]	375	igrd=3 ! V grid
[3651]	376	DO ib = 1, ilen0(igrd)
	377	dta%v2d(ib) = dta%v2d(ib) + zramp(td%v (ib,itide,1)z_cost(itide) + td%v (ib,itide,2)*z_sist(itide))
[1125]	378	END DO
	379	END DO
	380	!
[3651]	381	IF( nn_timing == 1 ) CALL timing_stop('bdytide_update')
[3294]	382	!
[3651]	383	END SUBROUTINE bdytide_update
[911]	384
[4292]	385	SUBROUTINE bdy_dta_tides( kt, kit, time_offset )
	386	!!----------------------------------------------------------------------
	387	!! * SUBROUTINE bdy_dta_tides *
	388	!!
	389	!! ** Purpose : - Add tidal forcing to ssh, u2d and v2d OBC data arrays.
	390	!!
	391	!!----------------------------------------------------------------------
	392	INTEGER, INTENT( in ) :: kt ! Main timestep counter
	393	INTEGER, INTENT( in ),OPTIONAL :: kit ! Barotropic timestep counter (for timesplitting option)
	394	INTEGER, INTENT( in ),OPTIONAL :: time_offset ! time offset in units of timesteps. NB. if kit
	395	! is present then units = subcycle timesteps.
	396	! time_offset = 0 => get data at "now" time level
	397	! time_offset = -1 => get data at "before" time level
	398	! time_offset = +1 => get data at "after" time level
	399	! etc.
	400	!!
	401	LOGICAL :: lk_first_btstp ! =.TRUE. if time splitting and first barotropic step
	402	INTEGER, DIMENSION(jpbgrd) :: ilen0
	403	INTEGER, DIMENSION(1:jpbgrd) :: nblen, nblenrim ! short cuts
	404	INTEGER :: itide, ib_bdy, ib, igrd ! loop indices
	405	INTEGER :: time_add ! time offset in units of timesteps
	406	REAL(wp) :: z_arg, z_sarg, zramp, zoff, z_cost, z_sist
	407	!!----------------------------------------------------------------------
	408
	409	IF( nn_timing == 1 ) CALL timing_start('bdy_dta_tides')
	410
	411	lk_first_btstp=.TRUE.
	412	IF ( PRESENT(kit).AND.( kit /= 1 ) ) THEN ; lk_first_btstp=.FALSE. ; ENDIF
	413
	414	time_add = 0
	415	IF( PRESENT(time_offset) ) THEN
	416	time_add = time_offset
	417	ENDIF
	418
	419	! Absolute time from model initialization:
	420	IF( PRESENT(kit) ) THEN
[6487]	421	z_arg = ( kt + (kit+time_add-1) / REAL(nn_baro,wp) ) * rdt
[4292]	422	ELSE
	423	z_arg = ( kt + time_add ) * rdt
	424	ENDIF
	425
	426	! Linear ramp on tidal component at open boundaries
	427	zramp = 1.
	428	IF (ln_tide_ramp) zramp = MIN(MAX( (z_arg - nit000rdt)/(rdttideramprday),0.),1.)
	429
	430	DO ib_bdy = 1,nb_bdy
	431
	432	IF ( nn_dyn2d_dta(ib_bdy) .ge. 2 ) THEN
	433
	434	nblen(1:jpbgrd) = idx_bdy(ib_bdy)%nblen(1:jpbgrd)
	435	nblenrim(1:jpbgrd) = idx_bdy(ib_bdy)%nblenrim(1:jpbgrd)
	436
	437	IF( cn_dyn2d(ib_bdy) == 'frs' ) THEN
	438	ilen0(:)=nblen(:)
	439	ELSE
	440	ilen0(:)=nblenrim(:)
	441	ENDIF
	442
	443	! We refresh nodal factors every day below
	444	! This should be done somewhere else
	445	IF ( nsec_day == NINT(0.5_wp * rdttra(1)) .AND. lk_first_btstp ) THEN
	446	!
	447	kt_tide = kt
	448	!
	449	IF(lwp) THEN
	450	WRITE(numout,*)
	451	WRITE(numout,*) 'bdy_tide_dta : Refresh nodal factors for tidal open bdy data at kt=',kt
	452	WRITE(numout,*) '~~~~~~~~~~~~~~ '
[11101]	453	IF(lflush) CALL flush(numout)
[4292]	454	ENDIF
	455	!
	456	CALL tide_init_elevation ( idx=idx_bdy(ib_bdy), td=tides(ib_bdy) )
	457	CALL tide_init_velocities( idx=idx_bdy(ib_bdy), td=tides(ib_bdy) )
	458	!
	459	ENDIF
	460	zoff = -kt_tide * rdt ! time offset relative to nodal factor computation time
	461	!
	462	! If time splitting, save data at first barotropic iteration
	463	IF ( PRESENT(kit) ) THEN
	464	IF ( lk_first_btstp ) THEN ! Save slow varying open boundary data:
[4758]	465	IF ( dta_bdy(ib_bdy)%ll_ssh ) dta_bdy_s(ib_bdy)%ssh(1:ilen0(1)) = dta_bdy(ib_bdy)%ssh(1:ilen0(1))
	466	IF ( dta_bdy(ib_bdy)%ll_u2d ) dta_bdy_s(ib_bdy)%u2d(1:ilen0(2)) = dta_bdy(ib_bdy)%u2d(1:ilen0(2))
	467	IF ( dta_bdy(ib_bdy)%ll_v2d ) dta_bdy_s(ib_bdy)%v2d(1:ilen0(3)) = dta_bdy(ib_bdy)%v2d(1:ilen0(3))
[4292]	468
	469	ELSE ! Initialize arrays from slow varying open boundary data:
[4758]	470	IF ( dta_bdy(ib_bdy)%ll_ssh ) dta_bdy(ib_bdy)%ssh(1:ilen0(1)) = dta_bdy_s(ib_bdy)%ssh(1:ilen0(1))
	471	IF ( dta_bdy(ib_bdy)%ll_u2d ) dta_bdy(ib_bdy)%u2d(1:ilen0(2)) = dta_bdy_s(ib_bdy)%u2d(1:ilen0(2))
	472	IF ( dta_bdy(ib_bdy)%ll_v2d ) dta_bdy(ib_bdy)%v2d(1:ilen0(3)) = dta_bdy_s(ib_bdy)%v2d(1:ilen0(3))
[4292]	473	ENDIF
	474	ENDIF
	475	!
	476	! Update open boundary data arrays:
	477	DO itide = 1, nb_harmo
	478	!
	479	z_sarg = (z_arg + zoff) * omega_tide(itide)
	480	z_cost = zramp * COS( z_sarg )
	481	z_sist = zramp * SIN( z_sarg )
	482	!
[4758]	483	IF ( dta_bdy(ib_bdy)%ll_ssh ) THEN
	484	igrd=1 ! SSH on tracer grid
	485	DO ib = 1, ilen0(igrd)
	486	dta_bdy(ib_bdy)%ssh(ib) = dta_bdy(ib_bdy)%ssh(ib) + &
	487	& ( tides(ib_bdy)%ssh(ib,itide,1)*z_cost + &
	488	& tides(ib_bdy)%ssh(ib,itide,2)*z_sist )
	489	END DO
	490	ENDIF
[4292]	491	!
[4758]	492	IF ( dta_bdy(ib_bdy)%ll_u2d ) THEN
	493	igrd=2 ! U grid
	494	DO ib = 1, ilen0(igrd)
	495	dta_bdy(ib_bdy)%u2d(ib) = dta_bdy(ib_bdy)%u2d(ib) + &
	496	& ( tides(ib_bdy)%u(ib,itide,1)*z_cost + &
	497	& tides(ib_bdy)%u(ib,itide,2)*z_sist )
	498	END DO
	499	ENDIF
[4292]	500	!
[4758]	501	IF ( dta_bdy(ib_bdy)%ll_v2d ) THEN
	502	igrd=3 ! V grid
	503	DO ib = 1, ilen0(igrd)
	504	dta_bdy(ib_bdy)%v2d(ib) = dta_bdy(ib_bdy)%v2d(ib) + &
	505	& ( tides(ib_bdy)%v(ib,itide,1)*z_cost + &
	506	& tides(ib_bdy)%v(ib,itide,2)*z_sist )
	507	END DO
	508	ENDIF
	509	END DO
[4292]	510	END IF
	511	END DO
	512	!
	513	IF( nn_timing == 1 ) CALL timing_stop('bdy_dta_tides')
	514	!
	515	END SUBROUTINE bdy_dta_tides
	516
[3651]	517	SUBROUTINE tide_init_elevation( idx, td )
[1125]	518	!!----------------------------------------------------------------------
[3651]	519	!! * ROUTINE tide_init_elevation *
[1125]	520	!!----------------------------------------------------------------------
[3651]	521	TYPE(OBC_INDEX), INTENT( in ) :: idx ! OBC indices
	522	TYPE(TIDES_DATA),INTENT( inout ) :: td ! tidal harmonics data
	523	!! * Local declarations
	524	INTEGER, DIMENSION(1) :: ilen0 !: length of boundary data (from OBC arrays)
	525	REAL(wp),ALLOCATABLE, DIMENSION(:) :: mod_tide, phi_tide
	526	INTEGER :: itide, igrd, ib ! dummy loop indices
[1125]	527
[3651]	528	igrd=1
	529	! SSH on tracer grid.
	530
	531	ilen0(1) = SIZE(td%ssh0(:,1,1))
	532
	533	ALLOCATE(mod_tide(ilen0(igrd)),phi_tide(ilen0(igrd)))
	534
	535	DO itide = 1, nb_harmo
	536	DO ib = 1, ilen0(igrd)
	537	mod_tide(ib)=SQRT(td%ssh0(ib,itide,1)2.+td%ssh0(ib,itide,2)2.)
	538	phi_tide(ib)=ATAN2(-td%ssh0(ib,itide,2),td%ssh0(ib,itide,1))
[1125]	539	END DO
[3651]	540	DO ib = 1 , ilen0(igrd)
	541	mod_tide(ib)=mod_tide(ib)*ftide(itide)
	542	phi_tide(ib)=phi_tide(ib)+v0tide(itide)+utide(itide)
	543	ENDDO
	544	DO ib = 1 , ilen0(igrd)
	545	td%ssh(ib,itide,1)= mod_tide(ib)*COS(phi_tide(ib))
	546	td%ssh(ib,itide,2)=-mod_tide(ib)*SIN(phi_tide(ib))
	547	ENDDO
	548	END DO
[911]	549
[3651]	550	DEALLOCATE(mod_tide,phi_tide)
[911]	551
[3651]	552	END SUBROUTINE tide_init_elevation
	553
	554	SUBROUTINE tide_init_velocities( idx, td )
[1125]	555	!!----------------------------------------------------------------------
[3651]	556	!! * ROUTINE tide_init_elevation *
[1125]	557	!!----------------------------------------------------------------------
[3651]	558	TYPE(OBC_INDEX), INTENT( in ) :: idx ! OBC indices
	559	TYPE(TIDES_DATA),INTENT( inout ) :: td ! tidal harmonics data
	560	!! * Local declarations
	561	INTEGER, DIMENSION(3) :: ilen0 !: length of boundary data (from OBC arrays)
	562	REAL(wp),ALLOCATABLE, DIMENSION(:) :: mod_tide, phi_tide
	563	INTEGER :: itide, igrd, ib ! dummy loop indices
[911]	564
[3651]	565	ilen0(2) = SIZE(td%u0(:,1,1))
	566	ilen0(3) = SIZE(td%v0(:,1,1))
[1125]	567
[3651]	568	igrd=2 ! U 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%u0(ib,itide,1)2.+td%u0(ib,itide,2)2.)
	575	phi_tide(ib)=ATAN2(-td%u0(ib,itide,2),td%u0(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%u(ib,itide,1)= mod_tide(ib)*COS(phi_tide(ib))
	583	td%u(ib,itide,2)=-mod_tide(ib)*SIN(phi_tide(ib))
	584	ENDDO
	585	END DO
[911]	586
[3651]	587	DEALLOCATE(mod_tide,phi_tide)
[911]	588
[3651]	589	igrd=3 ! V grid.
[911]	590
[3651]	591	ALLOCATE(mod_tide(ilen0(igrd)),phi_tide(ilen0(igrd)))
[911]	592
[3651]	593	DO itide = 1, nb_harmo
	594	DO ib = 1, ilen0(igrd)
	595	mod_tide(ib)=SQRT(td%v0(ib,itide,1)2.+td%v0(ib,itide,2)2.)
	596	phi_tide(ib)=ATAN2(-td%v0(ib,itide,2),td%v0(ib,itide,1))
	597	END DO
	598	DO ib = 1, ilen0(igrd)
	599	mod_tide(ib)=mod_tide(ib)*ftide(itide)
	600	phi_tide(ib)=phi_tide(ib)+v0tide(itide)+utide(itide)
	601	ENDDO
	602	DO ib = 1, ilen0(igrd)
	603	td%v(ib,itide,1)= mod_tide(ib)*COS(phi_tide(ib))
	604	td%v(ib,itide,2)=-mod_tide(ib)*SIN(phi_tide(ib))
	605	ENDDO
	606	END DO
[1125]	607
[3651]	608	DEALLOCATE(mod_tide,phi_tide)
[911]	609
[3651]	610	END SUBROUTINE tide_init_velocities
[911]	611	#else
[1125]	612	!!----------------------------------------------------------------------
	613	!! Dummy module NO Unstruct Open Boundary Conditions for tides
	614	!!----------------------------------------------------------------------
[911]	615	CONTAINS
[3651]	616	SUBROUTINE bdytide_init ! Empty routine
[9583]	617	IMPLICIT NONE
[3651]	618	WRITE(,) 'bdytide_init: You should not have seen this print! error?'
	619	END SUBROUTINE bdytide_init
	620	SUBROUTINE bdytide_update( kt, jit ) ! Empty routine
[9583]	621	IMPLICIT NONE
	622	INTEGER, INTENT( in ) :: kt ! Main timestep counter
	623	INTEGER,INTENT(in),OPTIONAL :: jit ! Barotropic timestep counter (for timesplitting option)
[3651]	624	WRITE(,) 'bdytide_update: You should not have seen this print! error?', kt, jit
	625	END SUBROUTINE bdytide_update
[4292]	626	SUBROUTINE bdy_dta_tides( kt, kit, time_offset ) ! Empty routine
[9583]	627	IMPLICIT NONE
[4292]	628	INTEGER, INTENT( in ) :: kt ! Dummy argument empty routine
	629	INTEGER, INTENT( in ),OPTIONAL :: kit ! Dummy argument empty routine
	630	INTEGER, INTENT( in ),OPTIONAL :: time_offset ! Dummy argument empty routine
[9583]	631	WRITE(,) 'bdy_dta_tides: You should not have seen this print! error?', kt, kit
[4292]	632	END SUBROUTINE bdy_dta_tides
[911]	633	#endif
	634
[1125]	635	!!======================================================================
[911]	636	END MODULE bdytides
[4292]	637

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

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

Context Navigation

source: branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/BDY/bdytides.F90 @ 15814

Download in other formats: