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.

obcdta.F90 in tags/nemo_v2_2/NEMO/OPA_SRC/OBC – NEMO

Context Navigation

source: tags/nemo_v2_2/NEMO/OPA_SRC/OBC/obcdta.F90 @ 4310

Last change on this file since 4310 was 591, checked in by opalod, 17 years ago

nemo_v2_bugfix_003 : CT : - add declaration of idvar

add initialization of itobce, itobcw, itobcn, itobcs variables

Property svn:eol-style set to native
Property svn:keywords set to Author Date Id Revision

File size: 51.9 KB

Line
1	MODULE obcdta
2	!!==============================================================================
3	!! * MODULE obcdta *
4	!! Open boundary data : read the data for the open boundaries.
5	!!==============================================================================
6	!! History :
7	!! ! 98-05 (J.M. Molines) Original code
8	!! 8.5 ! 02-10 (C. Talandier, A-M. Treguier) Free surface, F90
9	!! 9.0 ! 04-06 (F. Durand, A-M. Treguier) Netcdf BC files on input
10	!!------------------------------------------------------------------------------
11	#if defined key_obc
12	!!------------------------------------------------------------------------------
13	!! 'key_obc' : Open Boundary Conditions
14	!!------------------------------------------------------------------------------
15	!! obc_dta : read u, v, t, s data along each open boundary
16	!! obc_dta_psi : read psi data along each open boundary (rigid lid only)
17	!!------------------------------------------------------------------------------
18	!! * Modules used
19	USE oce ! ocean dynamics and tracers
20	USE dom_oce ! ocean space and time domain
21	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
22	USE phycst ! physical constants
23	USE obc_oce ! ocean open boundary conditions
24	USE daymod ! calendar
25	USE in_out_manager ! I/O logical units
26	USE lib_mpp ! distributed memory computing
27	USE dynspg_oce ! choice/control of key cpp for surface pressure gradient
28	USE iom
29	# if defined key_dynspg_rl
30	USE obccli ! climatological obc, use only in rigid-lid case
31	# endif
32
33	IMPLICIT NONE
34	PRIVATE
35
36	!! * Accessibility
37	PUBLIC obc_dta ! routines called by step.F90
38	PUBLIC obc_dta_bt ! routines called by dynspg_ts.F90
39
40	!! * Shared module variables
41	INTEGER :: &
42	nlecto, & ! switch for the first read
43	ntobc1, & ! first record used
44	ntobc2, & ! second record used
45	ntobc ! number of time steps in OBC files
46
47	REAL(wp), DIMENSION(:), ALLOCATABLE :: tcobc ! time_counter variable of BCs
48
49	!! * Substitutions
50	# include "domzgr_substitute.h90"
51	# include "obc_vectopt_loop_substitute.h90"
52	!!---------------------------------------------------------------------------------
53	!! OPA 9.0 , LODYC-IPSL (2003)
54	!! $Header$
55	!! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt
56	!!---------------------------------------------------------------------------------
57
58	CONTAINS
59
60	SUBROUTINE obc_dta( kt )
61	!!--------------------------------------------------------------------
62	!! * SUBROUTINE obc_dta *
63	!!
64	!! ** Purpose : Find the climatological boundary arrays for the specified date,
65	!! The boundary arrays are netcdf files. Three possible cases:
66	!! - one time frame only in the file (time dimension = 1).
67	!! in that case the boundary data does not change in time.
68	!! - many time frames. In that case, if we have 12 frames
69	!! we assume monthly fields.
70	!! Else, we assume that time_counter is in seconds
71	!! since the beginning of either the current year or a reference
72	!! year given in the namelist.
73	!! (no check is done so far but one would have to check the "unit"
74	!! attribute of variable time_counter).
75	!!
76	!! History :
77	!! ! 98-05 (J.M. Molines) Original code
78	!! 8.5 ! 02-10 (C. Talandier, A-M. Treguier) Free surface, F90
79	!! 9.0 ! 04-06 (F. Durand, A-M. Treguier) Netcdf BC files on input
80	!!--------------------------------------------------------------------
81	!! * Arguments
82	INTEGER, INTENT( in ) :: kt ! ocean time-step index
83
84	!! * Local declarations
85	INTEGER :: ji, jj, jk, ii, ij ! dummy loop indices
86	INTEGER :: itimo, iman, imois
87	INTEGER :: i15
88	REAL(wp) :: zxy
89	!! * Ajouts FD
90	INTEGER :: isrel ! number of seconds since 1/1/1992
91	INTEGER, DIMENSION(1) :: itobce, itobcw, & ! number of time steps in OBC files
92	itobcs, itobcn ! " " " "
93	INTEGER :: istop
94	INTEGER :: iprint ! frequency for printouts.
95	INTEGER :: idvar, id_e, id_w, id_n, id_s ! file identifiers
96	LOGICAL :: llnot_done
97	CHARACTER(LEN=25) :: cl_vname
98	!!--------------------------------------------------------------------
99
100	IF( lk_dynspg_rl ) THEN
101	CALL obc_dta_psi (kt) ! update bsf data at open boundaries
102	IF ( nobc_dta == 1 .AND. kt == nit000 ) CALL ctl_stop( 'obcdta : time-variable psi boundary data not allowed yet' )
103	ENDIF
104
105	! 1. First call: check time frames available in files.
106	! -------------------------------------------------------
107
108	IF ( kt == nit000 ) THEN
109
110	nlecto = 0
111	itobce(1) = 0 ; itobcw(1) = 0
112	itobcn(1) = 0 ; itobcs(1) = 0
113
114	IF (lwp) WRITE(numout,*)
115	IF (lwp) WRITE(numout,*) 'obc_dta : find boundary data'
116	IF (lwp) WRITE(numout,*) '~~~~~~~'
117
118	IF ( nobc_dta == 0 ) THEN
119	IF(lwp) WRITE(numout,*) ' OBC data taken from initial conditions.'
120	ntobc1 = 1
121	ntobc2 = 1
122	ELSE
123	IF (lwp) WRITE(numout,*) ' OBC data taken from netcdf files.'
124	IF (lwp) WRITE(numout,*) ' climatology (T/F):',ln_obc_clim
125	! check the number of time steps in the files.
126	cl_vname = 'time_counter'
127	IF ( lp_obc_east ) THEN
128	CALL iom_open ( 'obceast_TS.nc' , id_e )
129	idvar = iom_varid( id_e, cl_vname, kdimsz = itobce )
130	ENDIF
131	IF ( lp_obc_west ) THEN
132	CALL iom_open ( 'obcwest_TS.nc' , id_w )
133	idvar = iom_varid( id_w, cl_vname, kdimsz = itobcw )
134	ENDIF
135	IF ( lp_obc_north ) THEN
136	CALL iom_open ( 'obcnorth_TS.nc', id_n )
137	idvar = iom_varid( id_n, cl_vname, kdimsz = itobcn )
138	ENDIF
139	IF ( lp_obc_south ) THEN
140	CALL iom_open ( 'obcsouth_TS.nc', id_s )
141	idvar = iom_varid( id_s, cl_vname, kdimsz = itobcs )
142	ENDIF
143
144	ntobc = MAX(itobce(1),itobcw(1),itobcn(1),itobcs(1))
145	istop = 0
146	IF ( lp_obc_east .AND. itobce(1) /= ntobc ) istop = 1
147	IF ( lp_obc_west .AND. itobcw(1) /= ntobc ) istop = 1
148	IF ( lp_obc_north .AND. itobcn(1) /= ntobc ) istop = 1
149	IF ( lp_obc_south .AND. itobcs(1) /= ntobc ) istop = 1
150	IF ( istop /= 0 ) THEN
151	WRITE(ctmp1,*) ' east, west, north, south: ', itobce(1), itobcw(1), itobcn(1), itobcs(1)
152	CALL ctl_stop( 'obcdta : all files must have the same number of time steps', ctmp1 )
153	ENDIF
154	IF ( ntobc == 1 ) THEN
155	IF ( lwp ) WRITE(numout,*) ' obcdta found one time step only in the OBC files'
156	ELSE
157	ALLOCATE (tcobc(ntobc))
158	llnot_done = .TRUE.
159	IF ( lp_obc_east ) THEN
160	IF ( llnot_done ) THEN
161	CALL iom_gettime ( id_e, TRIM(cl_vname), tcobc )
162	llnot_done = .FALSE.
163	ENDIF
164	CALL iom_close (id_e)
165	ENDIF
166	IF ( lp_obc_west ) THEN
167	IF ( llnot_done ) THEN
168	CALL iom_gettime ( id_w, TRIM(cl_vname), tcobc )
169	llnot_done = .FALSE.
170	ENDIF
171	CALL iom_close (id_w)
172	ENDIF
173	IF ( lp_obc_north ) THEN
174	IF ( llnot_done ) THEN
175	CALL iom_gettime ( id_n, TRIM(cl_vname), tcobc )
176	llnot_done = .FALSE.
177	ENDIF
178	CALL iom_close (id_n)
179	ENDIF
180	IF ( lp_obc_south ) THEN
181	IF ( llnot_done ) THEN
182	CALL iom_gettime ( id_s, TRIM(cl_vname), tcobc )
183	llnot_done = .FALSE.
184	ENDIF
185	CALL iom_close (id_s)
186	ENDIF
187	IF ( lwp ) WRITE(numout,*) ' obcdta found', ntobc,' time steps in the OBC files'
188	IF ( .NOT. ln_obc_clim .AND. ntobc == 12 ) THEN
189	IF ( lwp ) WRITE(numout,*) ' WARNING: With monthly data we assume climatology'
190	ln_obc_clim = .true.
191	ENDIF
192	ENDIF
193	ENDIF
194
195	! 1.1 Tangential velocities set to zero
196	! --------------------------------------
197	IF( lp_obc_east ) vfoe = 0.e0
198	IF( lp_obc_west ) vfow = 0.e0
199	IF( lp_obc_south ) ufos = 0.e0
200	IF( lp_obc_north ) ufon = 0.e0
201
202	! 1.2 Data temperature, salinity, normal velocities set to zero
203	! or initial conditions if nobc_dta == 0
204	! --------------------------------------------------------------
205
206	IF( lp_obc_east ) THEN
207	! initialisation to zero
208	sedta(:,:,:) = 0.e0
209	tedta(:,:,:) = 0.e0
210	uedta(:,:,:) = 0.e0
211	! ! ================== !
212	IF( nobc_dta == 0 ) THEN ! initial state used
213	! ! ================== !
214	! Fills sedta, tedta, uedta (global arrays)
215	! Remark: this works for njzoom = 1.
216	! Should the definition of ij include njzoom?
217	DO ji = nie0, nie1
218	DO jk = 1, jpkm1
219	DO jj = nje0p1, nje1m1
220	ij = jj -1 + njmpp
221	sedta(ij,jk,1) = sn(ji,jj,jk)*tmask(ji,jj,jk)
222	tedta(ij,jk,1) = tn(ji,jj,jk)*tmask(ji,jj,jk)
223	uedta(ij,jk,1) = un(ji,jj,jk)*umask(ji,jj,jk)
224	END DO
225	END DO
226	END DO
227	ENDIF
228	ENDIF
229
230	IF( lp_obc_west ) THEN
231	! initialisation to zero
232	swdta(:,:,:) = 0.e0
233	twdta(:,:,:) = 0.e0
234	uwdta(:,:,:) = 0.e0
235	! ! ================== !
236	IF( nobc_dta == 0 ) THEN ! initial state used !
237	! ! ================== !
238	! Fills swdta, twdta, uwdta (global arrays)
239	! Remark: this works for njzoom = 1.
240	! Should the definition of ij include njzoom?
241	DO ji = niw0, niw1
242	DO jk = 1, jpkm1
243	DO jj = njw0p1, njw1m1
244	ij = jj -1 + njmpp
245	swdta(ij,jk,1) = sn(ji,jj,jk)*tmask(ji,jj,jk)
246	twdta(ij,jk,1) = tn(ji,jj,jk)*tmask(ji,jj,jk)
247	uwdta(ij,jk,1) = un(ji,jj,jk)*umask(ji,jj,jk)
248	END DO
249	END DO
250	END DO
251	ENDIF
252	ENDIF
253
254	IF( lp_obc_north) THEN
255	! initialisation to zero
256	sndta(:,:,:) = 0.e0
257	tndta(:,:,:) = 0.e0
258	vndta(:,:,:) = 0.e0
259	! ! ================== !
260	IF( nobc_dta == 0 ) THEN ! initial state used
261	! ! ================== !
262	! Fills sndta, tndta, vndta (global arrays)
263	! Remark: this works for njzoom = 1.
264	! Should the definition of ij include njzoom?
265	DO jj = njn0, njn1
266	DO jk = 1, jpkm1
267	DO ji = nin0p1, nin1m1
268	ii = ji -1 + nimpp
269	sndta(ii,jk,1) = sn(ji,jj,jk)*tmask(ji,jj,jk)
270	tndta(ii,jk,1) = tn(ji,jj,jk)*tmask(ji,jj,jk)
271	vndta(ii,jk,1) = vn(ji,jj,jk)*vmask(ji,jj,jk)
272	END DO
273	END DO
274	END DO
275	ENDIF
276	ENDIF
277
278	IF( lp_obc_south ) THEN
279	! initialisation to zero
280	ssdta(:,:,:) = 0.e0
281	tsdta(:,:,:) = 0.e0
282	vsdta(:,:,:) = 0.e0
283	! ! ================== !
284	IF( nobc_dta == 0 ) THEN ! initial state used
285	! ! ================== !
286	! Fills ssdta, tsdta, vsdta (global arrays)
287	! Remark: this works for njzoom = 1.
288	! Should the definition of ij include njzoom?
289	DO jj = njs0, njs1
290	DO jk = 1, jpkm1
291	DO ji = nis0p1, nis1m1
292	ii = ji -1 + nimpp
293	ssdta(ii,jk,1) = sn(ji,jj,jk)*tmask(ji,jj,jk)
294	tsdta(ii,jk,1) = tn(ji,jj,jk)*tmask(ji,jj,jk)
295	vsdta(ii,jk,1) = vn(ji,jj,jk)*vmask(ji,jj,jk)
296	END DO
297	END DO
298	END DO
299	ENDIF
300	ENDIF
301
302	ENDIF ! end if kt == nit000
303
304	! 2. Initialize the time we are at.
305	! Does this every time the routine is called,
306	! excepted when nobc_dta = 0
307	!---------------------------------------------------------------------
308	IF( nobc_dta == 0 ) THEN
309	itimo = 1
310	zxy = 0
311	ELSE
312	IF( ntobc == 1 ) THEN
313	itimo = 1
314	ELSE IF( ntobc == 12 ) THEN ! BC are monthly
315	! we assume we have climatology in that case
316	iman = 12
317	i15 = nday / 16
318	imois = nmonth + i15 - 1
319	IF( imois == 0 ) imois = iman
320	itimo = imois
321	ELSE
322	IF(lwp) WRITE(numout,*) 'data other than constant or monthly', kt
323	iman = ntobc
324	itimo = FLOOR( kt*rdt / (tcobc(2)-tcobc(1)) )
325	isrel = kt*rdt
326	ENDIF
327	ENDIF
328
329	! 2.1 Read two records in the file if necessary
330	! ---------------------------------------------
331	IF( ( nobc_dta == 1 ) .AND. ( ( kt == nit000 .AND. nlecto == 0 ) .OR. itimo /= ntobc1 ) ) THEN
332	nlecto = 1
333
334	! Calendar computation
335	IF( ntobc == 1 ) THEN ! BC are constant in time
336	ntobc1 = 1
337	ntobc2 = 1
338	ELSE IF( ntobc == 12 ) THEN ! BC are monthly
339	ntobc1 = itimo ! first file record used
340	ntobc2 = ntobc1 + 1 ! last file record used
341	ntobc1 = MOD( ntobc1, iman )
342	IF( ntobc1 == 0 ) ntobc1 = iman
343	ntobc2 = MOD( ntobc2, iman )
344	IF( ntobc2 == 0 ) ntobc2 = iman
345	IF( lwp ) THEN
346	WRITE(numout,*) ' read monthly obc first record file used ntobc1 ', ntobc1
347	WRITE(numout,*) ' read monthly obc last record file used ntobc2 ', ntobc2
348	ENDIF
349	ELSE
350	isrel=kt*rdt
351	ntobc1 = itimo ! first file record used
352	ntobc2 = ntobc1 + 1 ! last file record used
353	ntobc1 = MOD( ntobc1, iman )
354	IF( ntobc1 == 0 ) ntobc1 = iman
355	ntobc2 = MOD( ntobc2, iman )
356	IF( ntobc2 == 0 ) ntobc2 = iman
357	IF(lwp) WRITE(numout,*) ' read obc first record file used ntobc1 ', ntobc1
358	IF(lwp) WRITE(numout,*) ' read obc last record file used ntobc2 ', ntobc2
359	ENDIF
360	! ======================= !
361	! BCs read !
362	! ! ======================= !
363	IF( lp_obc_east ) THEN
364	! ... Read datafile and set temperature, salinity and normal velocity
365	! ... initialise the sedta, tedta, uedta arrays
366	CALL iom_open ( 'obceast_TS.nc' , id_e )
367	CALL iom_get ( id_e, jpdom_unknown, 'votemper', tedta(:,:,1), ktime=ntobc1 )
368	CALL iom_get ( id_e, jpdom_unknown, 'votemper', tedta(:,:,2), ktime=ntobc2 )
369	CALL iom_get ( id_e, jpdom_unknown, 'vosaline', sedta(:,:,1), ktime=ntobc1 )
370	CALL iom_get ( id_e, jpdom_unknown, 'vosaline', sedta(:,:,2), ktime=ntobc2 )
371	CALL iom_close (id_e)
372	!
373	CALL iom_open ( 'obceast_U.nc' , id_e )
374	CALL iom_get ( id_e, jpdom_unknown, 'vozocrtx', uedta(:,:,1), ktime=ntobc1 )
375	CALL iom_get ( id_e, jpdom_unknown, 'vozocrtx', uedta(:,:,2), ktime=ntobc2 )
376	CALL iom_close ( id_e )
377	! Usually printout is done only once at kt = nit000,
378	! unless nprint (namelist) > 1
379	IF( lwp .AND. ( kt == nit000 .OR. nprint /= 0 ) ) THEN
380	WRITE(numout,*)
381	WRITE(numout,*) ' Read East OBC data records ', ntobc1, ntobc2
382	iprint = (jpjef-jpjed+1)/20 +1
383	WRITE(numout,*) ' Temperature record 1 - printout every 3 level'
384	CALL prihre( tedta(:,:,1),jpjef-jpjed+1,jpk,1,jpjef-jpjed+1,iprint, &
385	& jpk, 1, -3, 1., numout )
386	WRITE(numout,*)
387	WRITE(numout,*) ' Salinity record 1 - printout every 3 level'
388	CALL prihre( sedta(:,:,1), jpjef-jpjed+1, jpk, 1, jpjef-jpjed+1, iprint, &
389	& jpk, 1, -3, 1., numout )
390	WRITE(numout,*)
391	WRITE(numout,*) ' Normal velocity U record 1 - printout every 3 level'
392	CALL prihre( uedta(:,:,1), jpjef-jpjed+1, jpk, 1, jpjef-jpjed+1, iprint, &
393	& jpk, 1, -3, 1., numout )
394	ENDIF
395	ENDIF
396
397	IF( lp_obc_west ) THEN
398	! ... Read datafile and set temperature, salinity and normal velocity
399	! ... initialise the swdta, twdta, uwdta arrays
400	CALL iom_open ( 'obcwest_TS.nc' , id_w )
401	CALL iom_get ( id_w, jpdom_unknown, 'votemper', twdta(:,:,1), ktime=ntobc1 )
402	CALL iom_get ( id_w, jpdom_unknown, 'votemper', twdta(:,:,2), ktime=ntobc2 )
403	CALL iom_get ( id_w, jpdom_unknown, 'vosaline', swdta(:,:,1), ktime=ntobc1 )
404	CALL iom_get ( id_w, jpdom_unknown, 'vosaline', swdta(:,:,2), ktime=ntobc2 )
405	CALL iom_close (id_w)
406	!
407	CALL iom_open ( 'obcwest_U.nc' , id_w )
408	CALL iom_get ( id_w, jpdom_unknown, 'vozocrtx', uwdta(:,:,1), ktime=ntobc1 )
409	CALL iom_get ( id_w, jpdom_unknown, 'vozocrtx', uwdta(:,:,2), ktime=ntobc2 )
410	CALL iom_close ( id_w )
411	!
412	IF( lwp .AND. ( kt == nit000 .OR. nprint /= 0 ) ) THEN
413	WRITE(numout,*)
414	WRITE(numout,*) ' Read West OBC data records ', ntobc1, ntobc2
415	iprint = (jpjwf-jpjwd+1)/20 +1
416	WRITE(numout,*) ' Temperature record 1 - printout every 3 level'
417	CALL prihre( twdta(:,:,1), jpjwf-jpjwd+1, jpk, 1, jpjwf-jpjwd+1, iprint, &
418	& jpk, 1, -3, 1., numout )
419	WRITE(numout,*)
420	WRITE(numout,*) ' Salinity record 1 - printout every 3 level'
421	CALL prihre( swdta(:,:,1), jpjwf-jpjwd+1, jpk, 1, jpjwf-jpjwd+1, iprint, &
422	& jpk, 1, -3, 1., numout )
423	WRITE(numout,*)
424	WRITE(numout,*) ' Normal velocity U record 1 - printout every 3 level'
425	CALL prihre( uwdta(:,:,1), jpjwf-jpjwd+1, jpk, 1, jpjwf-jpjwd+1, iprint, &
426	& jpk, 1, -3, 1., numout )
427	ENDIF
428	ENDIF
429
430	IF( lp_obc_north ) THEN
431	CALL iom_open ( 'obcnorth_TS.nc', id_n )
432	CALL iom_get ( id_n, jpdom_unknown, 'votemper', tndta(:,:,1), ktime=ntobc1 )
433	CALL iom_get ( id_n, jpdom_unknown, 'votemper', tndta(:,:,2), ktime=ntobc2 )
434	CALL iom_get ( id_n, jpdom_unknown, 'vosaline', sndta(:,:,1), ktime=ntobc1 )
435	CALL iom_get ( id_n, jpdom_unknown, 'vosaline', sndta(:,:,2), ktime=ntobc2 )
436	CALL iom_close ( id_n )
437	!
438	CALL iom_open ( 'obcnorth_V.nc', id_n )
439	CALL iom_get ( id_n, jpdom_unknown, 'vomecrty', vndta(:,:,1), ktime=ntobc1 )
440	CALL iom_get ( id_n, jpdom_unknown ,'vomecrty', vndta(:,:,2), ktime=ntobc2 )
441	CALL iom_close ( id_n )
442	!
443	IF( lwp .AND. ( kt == nit000 .OR. nprint /= 0 ) ) THEN
444	WRITE(numout,*)
445	WRITE(numout,*) ' Read North OBC data records ', ntobc1, ntobc2
446	iprint = (jpinf-jpind+1)/20 +1
447	WRITE(numout,*) ' Temperature record 1 - printout every 3 level'
448	CALL prihre( tndta(:,:,1), jpinf-jpind+1, jpk, 1, jpinf-jpind+1, iprint, &
449	& jpk, 1, -3, 1., numout )
450	WRITE(numout,*)
451	WRITE(numout,*) ' Salinity record 1 - printout every 3 level'
452	CALL prihre( sndta(:,:,1), jpinf-jpind+1, jpk, 1, jpinf-jpind+1, iprint, &
453	& jpk, 1, -3, 1., numout )
454	WRITE(numout,*)
455	WRITE(numout,*) ' Normal velocity V record 1 - printout every 3 level'
456	CALL prihre( vndta(:,:,1), jpinf-jpind+1, jpk, 1, jpinf-jpind+1, iprint, &
457	& jpk, 1, -3, 1., numout )
458	ENDIF
459	ENDIF
460
461	IF( lp_obc_south ) THEN
462	CALL iom_open ( 'obcsouth_TS.nc', id_s )
463	CALL iom_get ( id_s, jpdom_unknown, 'votemper', tsdta(:,:,1), ktime=ntobc1 )
464	CALL iom_get ( id_s, jpdom_unknown, 'votemper', tsdta(:,:,2), ktime=ntobc2 )
465	CALL iom_get ( id_s, jpdom_unknown, 'vosaline', ssdta(:,:,1), ktime=ntobc1 )
466	CALL iom_get ( id_s, jpdom_unknown, 'vosaline', ssdta(:,:,2), ktime=ntobc2 )
467	CALL iom_close ( id_s )
468	!
469	CALL iom_open ( 'obcsouth_V.nc', id_s )
470	CALL iom_get ( id_s, jpdom_unknown, 'vomecrty', vsdta(:,:,1), ktime=ntobc1 )
471	CALL iom_get ( id_s, jpdom_unknown ,'vomecrty', vsdta(:,:,2), ktime=ntobc2 )
472	CALL iom_close ( id_s )
473	!
474	IF( lwp .AND. ( kt == nit000 .OR. nprint /= 0 ) ) THEN
475	WRITE(numout,*)
476	WRITE(numout,*) ' Read South OBC data records ', ntobc1, ntobc2
477	iprint = (jpisf-jpisd+1)/20 +1
478	WRITE(numout,*) ' Temperature record 1 - printout every 3 level'
479	CALL prihre( tsdta(:,:,1), jpisf-jpisd+1, jpk, 1, jpisf-jpisd+1, iprint, &
480	& jpk, 1, -3, 1., numout )
481	WRITE(numout,*)
482	WRITE(numout,*) ' Salinity record 1 - printout every 3 level'
483	CALL prihre( ssdta(:,:,1), jpisf-jpisd+1, jpk, 1, jpisf-jpisd+1, iprint, &
484	& jpk, 1, -3, 1., numout )
485	WRITE(numout,*)
486	WRITE(numout,*) ' Normal velocity V record 1 - printout every 3 level'
487	CALL prihre( vsdta(:,:,1), jpisf-jpisd+1, jpk, 1, jpisf-jpisd+1, iprint, &
488	& jpk, 1, -3, 1., numout )
489	ENDIF
490	ENDIF
491
492	ELSE
493
494	nlecto = 0 ! no reading of OBC barotropic data
495
496	ENDIF ! end of the test on the condition to read or not the files
497
498	! 3. Call at every time step :
499	! Linear interpolation of BCs to current time step
500	! ----------------------------------------------------
501
502	IF( ntobc == 1 .OR. nobc_dta == 0 ) THEN
503	zxy = 0.
504	ELSE IF( ntobc == 12 ) THEN
505	zxy = FLOAT( nday + 15 - 30 * i15 ) / 30.
506	ELSE
507	zxy = (tcobc(ntobc1)-FLOAT(isrel))/(tcobc(ntobc1)-tcobc(ntobc2))
508	ENDIF
509
510	IF( lp_obc_east ) THEN
511	! fills sfoe, tfoe, ufoe (local to each processor)
512	DO jk = 1, jpkm1
513	DO jj = nje0p1, nje1m1
514	ij = jj -1 + njmpp
515	sfoe(jj,jk) = ( zxy * sedta(ij,jk,2) + &
516	& (1.-zxy) * sedta(ij,jk,1) ) * temsk(jj,jk)
517	tfoe(jj,jk) = ( zxy * tedta(ij,jk,2) + &
518	& (1.-zxy) * tedta(ij,jk,1) ) * temsk(jj,jk)
519	ufoe(jj,jk) = ( zxy * uedta(ij,jk,2) + &
520	& (1.-zxy) * uedta(ij,jk,1) ) * uemsk(jj,jk)
521	END DO
522	END DO
523	ENDIF
524
525	IF( lp_obc_west ) THEN
526	! fills sfow, tfow, ufow (local to each processor)
527	DO jk = 1, jpkm1
528	DO jj = njw0p1, njw1m1
529	ij = jj -1 + njmpp
530	sfow(jj,jk) = ( zxy * swdta(ij,jk,2) + &
531	& (1.-zxy) * swdta(ij,jk,1) ) * twmsk(jj,jk)
532	tfow(jj,jk) = ( zxy * twdta(ij,jk,2) + &
533	& (1.-zxy) * twdta(ij,jk,1) ) * twmsk(jj,jk)
534	ufow(jj,jk) = ( zxy * uwdta(ij,jk,2) + &
535	& (1.-zxy) * uwdta(ij,jk,1) ) * uwmsk(jj,jk)
536	END DO
537	END DO
538	ENDIF
539
540	IF( lp_obc_north ) THEN
541	! fills sfon, tfon, vfon (local to each processor)
542	DO jk = 1, jpkm1
543	DO ji = nin0p1, nin1m1
544	ii = ji -1 + nimpp
545	sfon(ji,jk) = ( zxy * sndta(ii,jk,2) + &
546	& (1.-zxy) * sndta(ii,jk,1) ) * tnmsk(ji,jk)
547	tfon(ji,jk) = ( zxy * tndta(ii,jk,2) + &
548	& (1.-zxy) * tndta(ii,jk,1) ) * tnmsk(ji,jk)
549	vfon(ji,jk) = ( zxy * vndta(ii,jk,2) + &
550	& (1.-zxy) * vndta(ii,jk,1) ) * vnmsk(ji,jk)
551	END DO
552	END DO
553	ENDIF
554
555	IF( lp_obc_south ) THEN
556	! fills sfos, tfos, vfos (local to each processor)
557	DO jk = 1, jpkm1
558	DO ji = nis0p1, nis1m1
559	ii = ji -1 + nimpp
560	sfos(ji,jk) = ( zxy * ssdta(ii,jk,2) + &
561	& (1.-zxy) * ssdta(ii,jk,1) ) * tsmsk(ji,jk)
562	tfos(ji,jk) = ( zxy * tsdta(ii,jk,2) + &
563	& (1.-zxy) * tsdta(ii,jk,1) ) * tsmsk(ji,jk)
564	vfos(ji,jk) = ( zxy * vsdta(ii,jk,2) + &
565	& (1.-zxy) * vsdta(ii,jk,1) ) * vsmsk(ji,jk)
566	END DO
567	END DO
568	ENDIF
569
570	END SUBROUTINE obc_dta
571
572	# if defined key_dynspg_rl
573	!!-----------------------------------------------------------------------------
574	!! Rigid-lid
575	!!-----------------------------------------------------------------------------
576
577	SUBROUTINE obc_dta_psi ( kt )
578	!!-----------------------------------------------------------------------------
579	!! * SUBROUTINE obc_dta_psi *
580	!!
581	!! ** Purpose :
582	!! Update the climatological streamfunction OBC at each time step.
583	!! Depends on the user's configuration. Here data are read only once
584	!! at the beginning of the run.
585	!!
586	!! ** Method :
587	!! 1. initialization
588	!! kbsfstart: number of time steps over which increase bsf
589	!! during initialization. This is provided for a smooth start
590	!! in cases where the transport is large (like on the Antarctic
591	!! continent). also note that when kbfstart=1, the transport
592	!! increases a lot in one time step and the precision usually
593	!! required for the solver may not be enough.
594	!! 2. set the time evolution of the climatological barotropic streamfunction
595	!! along the isolated coastlines ( gcbic(jnic) ).
596	!! 3. set the climatological barotropic streamfunction at the boundary.
597	!!
598	!! The last two steps are done only at first step (nit000) or if kt <= kbfstart
599	!!
600	!! History :
601	!! ! 97-08 (G. Madec, J.M. Molines)
602	!! 8.5 ! 02-10 (C. Talandier, A-M. Treguier) Free surface, F90
603	!! 9.0 ! 05-11 (V. Garnier) Surface pressure gradient organization
604	!!----------------------------------------------------------------------------
605	!! * Arguments
606	INTEGER, INTENT( in ) :: kt ! ocean time-step index
607
608	!! * Local declarations
609	INTEGER :: ji, jj, jnic, jip ! dummy loop indices
610	INTEGER :: inum = 11 ! temporary logical unit
611	INTEGER :: ip, ii, ij, iii, ijj
612	INTEGER :: kbsfstart
613	REAL(wp) :: zsver1, zsver2, zsver3, z2dtr, zcoef
614	!!----------------------------------------------------------------------------
615
616	! 1. initialisation
617	! -----------------
618
619	kbsfstart = 1
620	zsver1 = bsfic0(1)
621	zsver2 = zsver1
622	IF( kt <= kbsfstart ) THEN
623	zcoef = float(kt)/float(kbsfstart)
624	ELSE
625	zcoef = 1.
626	END IF
627	bsfic(1) = zsver1*zcoef
628	IF( lwp .AND. ( kt <= kbsfstart ) ) THEN
629	IF(lwp) WRITE(numout,*)' '
630	IF(lwp) WRITE(numout,*)'obcdta: spinup phase in obc_dta_psi routine'
631	IF(lwp) WRITE(numout,*)'~~~~~~ it=',kt,' OBC: spinup coef: ', &
632	zcoef, ' and transport: ',bsfic(1)
633	END IF
634
635	zsver2 = bsfic(1)-bsfic(2)
636	zsver3 = bsfic(2)
637
638	! 2. Right hand side of the barotropic elliptic equation (isolated coastlines)
639	! ----------------------------------------------------------------------------
640
641	IF( ( neuler == 0 ) .AND. ( kt == nit000 ) ) THEN
642	z2dtr = 1./rdt
643	ELSE
644	z2dtr = 1./2./rdt
645	END IF
646	! ... bsfb(ii,ij) should be constant but due to the Asselin filter it
647	! ... converges asymptotically towards bsfic(jnic)
648	! ... However, bsfb(ii,ij) is constant along the same coastlines
649	! ... ---> can be improved using an extra array for storing bsficb (before)
650	IF( nbobc > 1 ) THEN
651	DO jnic = 1,nbobc - 1
652	gcbic(jnic) = 0.e0
653	ip=mnic(0,jnic)
654	DO jip = 1,ip
655	ii = miic(jip,0,jnic)
656	ij = mjic(jip,0,jnic)
657	IF( ii >= nldi+ nimpp - 1 .AND. ii <= nlei+ nimpp - 1 .AND. &
658	ij >= nldj+ njmpp - 1 .AND. ij <= nlej+ njmpp - 1 ) THEN
659	iii=ii-nimpp+1
660	ijj=ij-njmpp+1
661	gcbic(jnic) = ( bsfic(jnic) - bsfb(iii,ijj) ) * z2dtr
662	END IF
663	END DO
664	END DO
665	END IF
666
667	IF( lk_mpp ) CALL mpp_isl( gcbic, 3 )
668
669	! 3. Update the climatological barotropic function at the boundary
670	! ----------------------------------------------------------------
671
672	IF( lpeastobc ) THEN
673
674	IF( kt == nit000 .OR. kt <= kbsfstart ) THEN
675	OPEN(inum,file='obceastbsf.dta')
676	READ(inum,*)
677	READ(inum,*)
678	READ(inum,*)
679	READ(inum,*)
680	READ(inum,*)
681	READ(inum,*) (bfoe(jj),jj=jpjed, jpjef)
682	CLOSE(inum)
683	END IF
684	DO jj=jpjed, jpjefm1
685	bfoe(jj)=bfoe(jj)*zcoef
686	END DO
687
688	END IF
689
690	IF( lpwestobc) THEN
691
692	IF( kt == nit000 .OR. kt <= kbsfstart ) then
693	OPEN(inum,file='obcwestbsf.dta')
694	READ(inum,*)
695	READ(inum,*)
696	READ(inum,*)
697	READ(inum,*)
698	READ(inum,*)
699	READ(inum,*) (bfow(jj),jj=jpjwd, jpjwf)
700	CLOSE(inum)
701	END IF
702	DO jj=jpjwd, jpjwfm1
703	bfow(jj)=bfow(jj)*zcoef
704	END DO
705
706	END IF
707
708	IF( lpsouthobc) THEN
709
710	IF( kt == nit000 .OR. kt <= kbsfstart ) THEN
711	OPEN(inum,file='obcsouthbsf.dta')
712	READ(inum,*)
713	READ(inum,*)
714	READ(inum,*)
715	READ(inum,*)
716	READ(inum,*)
717	READ(inum,*) (bfos(jj),jj=jpisd, jpisf)
718	CLOSE(inum)
719	END IF
720	DO ji=jpisd, jpisfm1
721	bfos(ji)=bfos(ji)*zcoef
722	END DO
723
724	END IF
725
726	IF( lpnorthobc) THEN
727	IF( kt == nit000 .OR. kt <= kbsfstart ) THEN
728	OPEN(inum,file='obcnorthbsf.dta')
729	READ(inum,*)
730	READ(inum,*)
731	READ(inum,*)
732	READ(inum,*)
733	READ(inum,*)
734	READ(inum,*) (bfon(jj),jj=jpind, jpinf)
735	CLOSE(inum)
736	END IF
737	DO ji=jpind, jpinfm1
738	bfon(ji)=bfon(ji)*zcoef
739	END DO
740
741	END IF
742
743	END SUBROUTINE obc_dta_psi
744	#else
745	!!-----------------------------------------------------------------------------
746	!! Default option
747	!!-----------------------------------------------------------------------------
748	SUBROUTINE obc_dta_psi ( kt ) ! Empty routine
749	!! * Arguments
750	INTEGER,INTENT(in) :: kt
751	WRITE(,) 'obc_dta_psi: You should not have seen this print! error?', kt
752	END SUBROUTINE obc_dta_psi
753	# endif
754
755
756	#if defined key_dynspg_ts \|\| defined key_dynspg_exp
757	SUBROUTINE obc_dta_bt( kt, kbt )
758	!!---------------------------------------------------------------------------
759	!! * SUBROUTINE obc_dta *
760	!!
761	!! ** Purpose : time interpolation of barotropic data for time-splitting scheme
762	!! Data at the boundary must be in m2/s
763	!!
764	!! History :
765	!! 9.0 ! 05-11 (V. garnier) Original code
766	!!---------------------------------------------------------------------------
767	!! * Arguments
768	INTEGER, INTENT( in ) :: kt ! ocean time-step index
769	INTEGER, INTENT( in ) :: kbt ! barotropic ocean time-step index
770
771	!! * Local declarations
772	INTEGER :: ji, jj, jk, ii, ij ! dummy loop indices
773	INTEGER :: id_e, id_w, id_n, id_s, fid ! file identifiers
774	INTEGER :: itimo, iman, imois, i15
775	INTEGER :: itobcm, itobcp, itimom, itimop
776	REAL(wp) :: zxy
777	INTEGER :: isrel, ikt ! number of seconds since 1/1/1992
778	INTEGER :: iprint ! frequency for printouts.
779
780	!!---------------------------------------------------------------------------
781
782	! 1. First call: check time frames available in files.
783	! -------------------------------------------------------
784
785	IF( kt == nit000 ) THEN
786
787	! 1.1 Barotropic tangential velocities set to zero
788	! -------------------------------------------------
789	IF( lp_obc_east ) vbtfoe(:) = 0.e0
790	IF( lp_obc_west ) vbtfow(:) = 0.e0
791	IF( lp_obc_south ) ubtfos(:) = 0.e0
792	IF( lp_obc_north ) ubtfon(:) = 0.e0
793
794	! 1.2 Sea surface height and normal barotropic velocities set to zero
795	! or initial conditions if nobc_dta == 0
796	! --------------------------------------------------------------------
797
798	IF( lp_obc_east ) THEN
799	! initialisation to zero
800	sshedta(:,:) = 0.e0
801	ubtedta(:,:) = 0.e0
802	! ! ================== !
803	IF( nobc_dta == 0 ) THEN ! initial state used !
804	! ! ================== !
805	! Fills sedta, tedta, uedta (global arrays)
806	! Remark: this works for njzoom = 1. Should the definition of ij include njzoom?
807	DO ji = nie0, nie1
808	DO jj = nje0p1, nje1m1
809	ij = jj -1 + njmpp
810	sshedta(ij,1) = sshn(ji+1,jj) * tmask(ji+1,jj,1)
811	END DO
812	END DO
813	ENDIF
814	ENDIF
815
816	IF( lp_obc_west) THEN
817	! initialisation to zero
818	sshwdta(:,:) = 0.e0
819	ubtwdta(:,:) = 0.e0
820	! ! ================== !
821	IF( nobc_dta == 0 ) THEN ! initial state used !
822	! ! ================== !
823	! Fills swdta, twdta, uwdta (global arrays)
824	! Remark: this works for njzoom = 1. Should the definition of ij include njzoom?
825	DO ji = niw0, niw1
826	DO jj = njw0p1, njw1m1
827	ij = jj -1 + njmpp
828	sshwdta(ij,1) = sshn(ji,jj) * tmask(ji,jj,1)
829	END DO
830	END DO
831	ENDIF
832	ENDIF
833
834	IF( lp_obc_north) THEN
835	! initialisation to zero
836	sshndta(:,:) = 0.e0
837	vbtndta(:,:) = 0.e0
838	! ! ================== !
839	IF( nobc_dta == 0 ) THEN ! initial state used !
840	! ! ================== !
841	! Fills sndta, tndta, vndta (global arrays)
842	! Remark: this works for njzoom = 1. Should the definition of ij include njzoom?
843	DO jj = njn0, njn1
844	DO ji = nin0p1, nin1m1
845	DO jk = 1, jpkm1
846	ii = ji -1 + nimpp
847	vbtndta(ii,1) = vbtndta(ii,1) + vndta(ii,jk,1)*fse3v(ji,jj,jk)
848	END DO
849	sshndta(ii,1) = sshn(ii,jj+1) * tmask(ji,jj+1,1)
850	END DO
851	END DO
852	ENDIF
853	ENDIF
854
855	IF( lp_obc_south) THEN
856	! initialisation to zero
857	ssdta(:,:,:) = 0.e0
858	tsdta(:,:,:) = 0.e0
859	vsdta(:,:,:) = 0.e0
860	sshsdta(:,:) = 0.e0
861	vbtsdta(:,:) = 0.e0
862	! ! ================== !
863	IF( nobc_dta == 0 ) THEN ! initial state used !
864	! ! ================== !
865	! Fills ssdta, tsdta, vsdta (global arrays)
866	! Remark: this works for njzoom = 1. Should the definition of ij include njzoom?
867	DO jj = njs0, njs1
868	DO ji = nis0p1, nis1m1
869	DO jk = 1, jpkm1
870	ii = ji -1 + nimpp
871	vbtsdta(ii,1) = vbtsdta(ii,1) + vsdta(ii,jk,1)*fse3v(ji,jj,jk)
872	END DO
873	sshsdta(ii,1) = sshn(ji,jj) * tmask(ii,jj,1)
874	END DO
875	END DO
876	ENDIF
877	ENDIF
878
879	ENDIF ! END IF kt == nit000
880
881	!!------------------------------------------------------------------------------------
882	! 2. Initialize the time we are at. Does this every time the routine is called,
883	! excepted when nobc_dta = 0
884	!
885	IF( nobc_dta == 0) THEN
886	itimo = 1
887	zxy = 0
888	ELSE
889	IF(ntobc == 1) THEN
890	itimo = 1
891	ELSE IF (ntobc == 12) THEN ! BC are monthly
892	! we assume we have climatology in that case
893	iman = 12
894	i15 = nday / 16
895	imois = nmonth + i15 - 1
896	IF( imois == 0 ) imois = iman
897	itimo = imois
898	ELSE
899	IF(lwp) WRITE(numout,*) 'data other than constant or monthly',kt
900	iman = ntobc
901	itimo = FLOOR( kt*rdt / tcobc(1))
902	isrel=kt*rdt
903	ENDIF
904	ENDIF
905
906	! 2. Read two records in the file if necessary
907	! ---------------------------------------------
908
909	IF( nobc_dta == 1 .AND. nlecto == 1 ) THEN
910
911	IF( lp_obc_east ) THEN
912	! ... Read datafile and set sea surface height and barotropic velocity
913	! ... initialise the sshedta, ubtedta arrays
914	sshedta(:,0) = sshedta(:,1)
915	ubtedta(:,0) = ubtedta(:,1)
916	CALL iom_open ( 'obceast_TS.nc', id_e )
917	CALL iom_get ( id_e, jpdom_unknown, 'vossurfh', sshedta(:,1), ktime=ntobc1 )
918	CALL iom_get ( id_e, jpdom_unknown, 'vossurfh', sshedta(:,2), ktime=ntobc2 )
919	IF( lk_dynspg_ts ) THEN
920	CALL iom_get (id_e, jpdom_unknown, 'vossurfh', sshedta(:,3), ktime=ntobc2+1 )
921	ENDIF
922	CALL iom_close ( id_e )
923	!
924	CALL iom_open ( 'obceast_U.nc', id_e )
925	CALL iom_get ( id_e, jpdom_unknown, 'vozoubt', ubtedta(:,1), ktime=ntobc1 )
926	CALL iom_get ( id_e, jpdom_unknown, 'vozoubt', ubtedta(:,2), ktime=ntobc2 )
927	IF( lk_dynspg_ts ) THEN
928	CALL iom_get ( id_e, jpdom_unknown, 'vozoubt', ubtedta(:,3), ktime=ntobc2+1 )
929	ENDIF
930	CALL iom_close ( id_e )
931	! ... Usually printout is done only once at kt = nit000, unless nprint (namelist) > 1
932	IF( lwp .AND. ( kt == nit000 .OR. nprint /= 0 ) ) THEN
933	WRITE(numout,*)
934	WRITE(numout,*) ' Read East OBC barotropic data records ', ntobc1, ntobc2
935	iprint = (jpjef-jpjed+1)/20 +1
936	WRITE(numout,*)
937	WRITE(numout,*) ' Sea surface height record 1'
938	CALL prihre( sshedta(:,1), jpjef-jpjed+1, 1, 1, jpjef-jpjed+1, iprint, 1, 1, -3, 1., numout )
939	WRITE(numout,*)
940	WRITE(numout,*) ' Normal transport (m2/s) record 1',iprint
941	CALL prihre( ubtedta(:,1), jpjef-jpjed+1, 1, 1, jpjef-jpjed+1, iprint, 1, 1, -3, 1., numout )
942	ENDIF
943	ENDIF
944
945	IF( lp_obc_west ) THEN
946	! ... Read datafile and set temperature, salinity and normal velocity
947	! ... initialise the swdta, twdta, uwdta arrays
948	sshwdta(:,0) = sshwdta(:,1)
949	ubtwdta(:,0) = ubtwdta(:,1)
950	CALL iom_open ( 'obcwest_TS.nc', id_w )
951	CALL iom_get ( id_w, jpdom_unknown, 'vossurfh', sshwdta(:,1), ktime=ntobc1 )
952	CALL iom_get ( id_w, jpdom_unknown, 'vossurfh', sshwdta(:,2), ktime=ntobc2 )
953	IF( lk_dynspg_ts ) THEN
954	CALL ( id_w, jpdom_unknown, 'vossurfh', sshwdta(:,3), ktime=ntobc2+1 )
955	ENDIF
956	CALL iom_close ( id_w )
957	!
958	CALL iom_open ( 'obcwest_U.nc', id_w )
959	CALL iom_get ( id_w, jpdom_unknown, 'vozoubt', ubtwdta(:,1), ktime=ntobc1 )
960	CALL iom_get ( id_w, jpdom_unknown, 'vozoubt', ubtwdta(:,2), ktime=ntobc2 )
961	IF( lk_dynspg_ts ) THEN
962	CALL iom_get ( id_w, jpdom_unknown, 'vozoubt', ubtwdta(:,3), ktime=ntobc2+1 )
963	ENDIF
964	CALL iom_close ( id_w )
965	! ... Usually printout is done only once at kt = nit000, unless nprint (namelist) > 1
966	IF( lwp .AND. ( kt == nit000 .OR. nprint /= 0 ) ) THEN
967	WRITE(numout,*)
968	WRITE(numout,*) ' Read West OBC barotropic data records ', ntobc1, ntobc2
969	iprint = (jpjwf-jpjwd+1)/20 +1
970	WRITE(numout,*)
971	WRITE(numout,*) ' Sea surface height record 1 - printout surface level'
972	CALL prihre( sshwdta(:,1), jpjwf-jpjwd+1, 1, 1, jpjwf-jpjwd+1, iprint, 1, 1, -3, 1., numout )
973	WRITE(numout,*)
974	WRITE(numout,*) ' Normal transport (m2/s) record 1'
975	CALL prihre( ubtwdta(:,1), jpjwf-jpjwd+1, 1, 1, jpjwf-jpjwd+1, iprint, 1, 1, -3, 1., numout )
976	ENDIF
977	ENDIF
978
979	IF( lp_obc_north) THEN
980	! ... Read datafile and set sea surface height and barotropic velocity
981	! ... initialise the sshndta, ubtndta arrays
982	sshndta(:,0) = sshndta(:,1)
983	vbtndta(:,0) = vbtndta(:,1)
984	CALL iom_open ( 'obcnorth_TS.nc', id_n )
985	CALL iom_get (id_n, jpdom_unknown, 'vossurfh', sshndta(:,1), ktime=ntobc1 )
986	CALL iom_get (id_n, jpdom_unknown, 'vossurfh', sshndta(:,2), ktime=ntobc2 )
987	IF( lk_dynspg_ts ) THEN
988	CALL iom_get (id_n, jpdom_unknown, 'vossurfh', sshndta(:,3), ktime=ntobc2+1 )
989	ENDIF
990	CALL iom_close ( id_n )
991
992	CALL iom_open ( 'obcnorth_V.nc', id_n )
993	CALL iom_get (id_n, jpdom_unknown, 'vomevbt', vbtndta(:,1), ktime=ntobc1 )
994	CALL iom_get (id_n, jpdom_unknown, 'vomevbt', vbtndta(:,2), ktime=ntobc2 )
995	IF( lk_dynspg_ts ) THEN
996	CALL iom_get (id_n, jpdom_unknown, 'vomevbt', vbtndta(:,3), ktime=ntobc2+1 )
997	ENDIF
998	CALL iom_close ( id_n )
999
1000	! ... Usually printout is done only once at kt = nit000, unless nprint (namelist) > 1
1001	IF( lwp .AND. ( kt == nit000 .OR. nprint /= 0 ) ) THEN
1002	WRITE(numout,*)
1003	WRITE(numout,*) ' Read North OBC barotropic data records ', ntobc1, ntobc2
1004	iprint = (jpinf-jpind+1)/20 +1
1005	WRITE(numout,*)
1006	WRITE(numout,*) ' Sea surface height record 1 - printout surface level'
1007	CALL prihre( sshndta(:,1), jpinf-jpind+1, 1, 1, jpinf-jpind+1, iprint, 1, 1, -3, 1., numout )
1008	WRITE(numout,*)
1009	WRITE(numout,*) ' Normal transport (m2/s) record 1'
1010	CALL prihre( vbtndta(:,1), jpinf-jpind+1, 1, 1, jpinf-jpind+1, iprint, 1, 1, -3, 1., numout )
1011	ENDIF
1012	ENDIF
1013
1014	IF( lp_obc_south) THEN
1015	! ... Read datafile and set sea surface height and barotropic velocity
1016	! ... initialise the sshsdta, ubtsdta arrays
1017	sshsdta(:,0) = sshsdta(:,1)
1018	vbtsdta(:,0) = vbtsdta(:,1)
1019	CALL iom_open ( 'obcsouth_TS.nc', id_s )
1020	CALL iom_get ( id_s, jpdom_unknown, 'vossurfh', sshsdta(:,1), ktime=ntobc1 )
1021	CALL iom_get ( id_s, jpdom_unknown, 'vossurfh', sshsdta(:,2), ktime=ntobc2 )
1022	IF( lk_dynspg_ts ) THEN
1023	CALL iom_get ( id_s, jpdom_unknown, 'vossurfh', sshsdta(:,3), ktime=ntobc2+1 )
1024	ENDIF
1025	CALL iom_close ( id_s )
1026
1027	CALL iom_open ( 'obcsouth_V.nc', id_s )
1028	CALL iom_get ( id_s, jpdom_unknown, 'vomevbt', vbtsdta(:,1), ktime=ntobc1 )
1029	CALL iom_get ( id_s, jpdom_unknown, 'vomevbt', vbtsdta(:,2), ktime=ntobc2 )
1030	IF( lk_dynspg_ts ) THEN
1031	CALL iom_get ( id_s, jpdom_unknown, 'vomevbt', vbtsdta(:,3), ktime=ntobc2+1 )
1032	ENDIF
1033	CALL iom_close ( id_s )
1034
1035	! ... Usually printout is done only once at kt = nit000, unless nprint (namelist) > 1
1036	IF( lwp .AND. ( kt == nit000 .OR. nprint /= 0 ) ) THEN
1037	WRITE(numout,*)
1038	WRITE(numout,*) ' Read South OBC barotropic data records ', ntobc1, ntobc2
1039	iprint = (jpisf-jpisd+1)/20 +1
1040	WRITE(numout,*)
1041	WRITE(numout,*) ' Sea surface height record 1 - printout surface level'
1042	CALL prihre( sshsdta(:,1), jpisf-jpisd+1, 1, 1, jpisf-jpisd+1, iprint, 1, 1, -3, 1., numout )
1043	WRITE(numout,*)
1044	WRITE(numout,*) ' Normal transport (m2/s) record 1'
1045	CALL prihre( vbtsdta(:,1), jpisf-jpisd+1, 1, 1, jpisf-jpisd+1, iprint, 1, 1, -3, 1., numout )
1046	ENDIF
1047	ENDIF
1048
1049	ENDIF ! end of the test on the condition to read or not the files
1050
1051	! 3. Call at every time step : Linear interpolation of BCs to current time step
1052	! ----------------------------------------------------------------------
1053
1054	IF( lk_dynspg_ts ) THEN
1055	isrel = (kt-1)rdt + kbtrdtbt
1056
1057	IF( nobc_dta == 1 ) THEN
1058	isrel = (kt-1)rdt + kbtrdtbt
1059	itimo = FLOOR( kt*rdt / (tcobc(2)-tcobc(1)) )
1060	itimom = FLOOR( (kt-1)*rdt / (tcobc(2)-tcobc(1)) )
1061	itimop = FLOOR( (kt+1)*rdt / (tcobc(2)-tcobc(1)) )
1062	IF( itimom == itimo .AND. itimop == itimo ) THEN
1063	itobcm = ntobc1
1064	itobcp = ntobc2
1065
1066	ELSEIF ( itimom <= itimo .AND. itimop == itimo ) THEN
1067	IF( FLOOR( isrel / (tcobc(2)-tcobc(1)) ) < itimo ) THEN
1068	itobcm = ntobc1-1
1069	itobcp = ntobc2-1
1070	ELSE
1071	itobcm = ntobc1
1072	itobcp = ntobc2
1073	ENDIF
1074
1075	ELSEIF ( itimom == itimo .AND. itimop >= itimo ) THEN
1076	IF( FLOOR( isrel / (tcobc(2)-tcobc(1)) ) < itimop ) THEN
1077	itobcm = ntobc1
1078	itobcp = ntobc2
1079	ELSE
1080	itobcm = ntobc1+1
1081	itobcp = ntobc2+1
1082	ENDIF
1083
1084	ELSEIF ( itimom == itimo-1 .AND. itimop == itimo+1 ) THEN
1085	IF( FLOOR( isrel / (tcobc(2)-tcobc(1)) ) < itimo ) THEN
1086	itobcm = ntobc1-1
1087	itobcp = ntobc2-1
1088	ELSEIF ( FLOOR( isrel / (tcobc(2)-tcobc(1)) ) < itimop ) THEN
1089	itobcm = ntobc1
1090	itobcp = ntobc2
1091	ELSEIF ( FLOOR( isrel / (tcobc(2)-tcobc(1)) ) == itimop ) THEN
1092	itobcm = ntobc1+1
1093	itobcp = ntobc2+2
1094	ELSE
1095	IF(lwp) WRITE(numout, *) 'obc_dta_bt: You should not have seen this print! error 1?'
1096	ENDIF
1097	ELSE
1098	IF(lwp) WRITE(numout, *) 'obc_dta_bt: You should not have seen this print! error 2?'
1099	ENDIF
1100
1101	ENDIF
1102
1103	ELSE IF( lk_dynspg_exp ) THEN
1104	isrel=kt*rdt
1105	itobcm = ntobc1
1106	itobcp = ntobc2
1107	ENDIF
1108
1109	IF( ntobc == 1 .OR. nobc_dta == 0 ) THEN
1110	zxy = 0.e0
1111	ELSE IF( ntobc == 12 ) THEN
1112	zxy = FLOAT( nday + 15 - 30 * i15 ) / 30.
1113	ELSE
1114	zxy = (tcobc(itobcm)-FLOAT(isrel)) / (tcobc(itobcm)-tcobc(itobcp))
1115	ENDIF
1116
1117	IF( lp_obc_east ) THEN ! fills sshfoe, ubtfoe (local to each processor)
1118	DO jj = nje0p1, nje1m1
1119	ij = jj -1 + njmpp
1120	sshfoe(jj) = ( zxy * sshedta(ij,2) + (1.-zxy) * sshedta(ij,1) ) * temsk(jj,1)
1121	ubtfoe(jj) = ( zxy * ubtedta(ij,2) + (1.-zxy) * ubtedta(ij,1) ) * uemsk(jj,1)
1122	END DO
1123	ENDIF
1124
1125	IF( lp_obc_west) THEN ! fills sshfow, ubtfow (local to each processor)
1126	DO jj = njw0p1, njw1m1
1127	ij = jj -1 + njmpp
1128	sshfow(jj) = ( zxy * sshwdta(ij,2) + (1.-zxy) * sshwdta(ij,1) ) * twmsk(jj,1)
1129	ubtfow(jj) = ( zxy * ubtwdta(ij,2) + (1.-zxy) * ubtwdta(ij,1) ) * uwmsk(jj,1)
1130	END DO
1131	ENDIF
1132
1133	IF( lp_obc_north) THEN ! fills sshfon, vbtfon (local to each processor)
1134	DO ji = nin0p1, nin1m1
1135	ii = ji -1 + nimpp
1136	sshfon(ji) = ( zxy * sshndta(ii,2) + (1.-zxy) * sshndta(ii,1) ) * tnmsk(ji,1)
1137	vbtfon(ji) = ( zxy * vbtndta(ii,2) + (1.-zxy) * vbtndta(ii,1) ) * vnmsk(ji,1)
1138	END DO
1139	ENDIF
1140
1141	IF( lp_obc_south) THEN ! fills sshfos, vbtfos (local to each processor)
1142	DO ji = nis0p1, nis1m1
1143	ii = ji -1 + nimpp
1144	sshfos(ji) = ( zxy * sshsdta(ii,2) + (1.-zxy) * sshsdta(ii,1) ) * tsmsk(ji,1)
1145	vbtfos(ji) = ( zxy * vbtsdta(ii,2) + (1.-zxy) * vbtsdta(ii,1) ) * vsmsk(ji,1)
1146	END DO
1147	ENDIF
1148
1149	END SUBROUTINE obc_dta_bt
1150
1151	#else
1152	!!-----------------------------------------------------------------------------
1153	!! Default option
1154	!!-----------------------------------------------------------------------------
1155	SUBROUTINE obc_dta_bt ( kt, kbt ) ! Empty routine
1156	!! * Arguments
1157	INTEGER,INTENT(in) :: kt
1158	INTEGER, INTENT( in ) :: kbt ! barotropic ocean time-step index
1159	WRITE(,) 'obc_dta_bt: You should not have seen this print! error?', kt
1160	WRITE(,) 'obc_dta_bt: You should not have seen this print! error?', kbt
1161	END SUBROUTINE obc_dta_bt
1162	#endif
1163
1164	#else
1165	!!--------------------------------------------------------------------
1166	!! default option : Dummy module NO Open Boundary Conditions
1167	!!--------------------------------------------------------------------
1168	CONTAINS
1169	SUBROUTINE obc_dta( kt ) ! Dummy routine
1170	INTEGER, INTENT (in) :: kt
1171	WRITE(,) 'obc_dta: You should not have seen this print! error?', kt
1172	END SUBROUTINE obc_dta
1173	SUBROUTINE obc_dta_bt( kt, jn) ! Dummy routine
1174	INTEGER, INTENT (in) :: kt, jn
1175	WRITE(,) 'obc_dta_bt: You should not have seen this print! error?', kt
1176	WRITE(,) 'obc_dta_bt: You should not have seen this print! error?', jn
1177	END SUBROUTINE obc_dta_bt
1178	#endif
1179
1180	!!=====================================================================
1181	END MODULE obcdta

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

Download in other formats: