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.

obcini.F90 in branches/2011/DEV_r2739_STFC_dCSE/NEMOGCM/NEMO/OPA_SRC/OBC – NEMO

Context Navigation

source: branches/2011/DEV_r2739_STFC_dCSE/NEMOGCM/NEMO/OPA_SRC/OBC/obcini.F90 @ 4400

Last change on this file since 4400 was 3211, checked in by spickles2, 12 years ago

Stephen Pickles, 11 Dec 2011

Commit to bring the rest of the DCSE NEMO development branch
in line with the latest development version. This includes
array index re-ordering of all OPA_SRC/.

Property svn:keywords set to Id

File size: 30.2 KB

Line
1	MODULE obcini
2	!!======================================================================
3	!! * MODULE obcini *
4	!! OBC initial state : Open boundary initial state
5	!!======================================================================
6	!! History : 8.0 ! 97-07 (J.M. Molines, G. Madec) Original code
7	!! NEMO 1.0 ! 02-11 (C. Talandier, A-M. Treguier) Free surface, F90
8	!! 2.0 ! 05-11 (V. Garnier) Surface pressure gradient organization
9	!!----------------------------------------------------------------------
10	#if defined key_obc
11	!!----------------------------------------------------------------------
12	!! 'key_obc' Open Boundary Conditions
13	!!----------------------------------------------------------------------
14	!! obc_init : initialization for the open boundary condition
15	!!----------------------------------------------------------------------
16	USE oce ! ocean dynamics and tracers variables
17	USE dom_oce ! ocean space and time domain variables
18	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
19	USE phycst ! physical constants
20	USE obc_oce ! open boundary condition: ocean
21	USE obcdta ! open boundary condition: data
22	USE in_out_manager ! I/O units
23	USE lib_mpp ! MPP library
24	USE dynspg_oce ! flag lk_dynspg_flt
25
26	IMPLICIT NONE
27	PRIVATE
28
29	PUBLIC obc_init ! routine called by opa.F90
30
31	!! * Control permutation of array indices
32	# include "oce_ftrans.h90"
33	# include "dom_oce_ftrans.h90"
34	# include "obc_oce_ftrans.h90"
35
36	!! * Substitutions
37	# include "obc_vectopt_loop_substitute.h90"
38	!!----------------------------------------------------------------------
39	!! NEMO/OPA 3.3 , NEMO Consortium (2010)
40	!! $Id$
41	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
42	!!----------------------------------------------------------------------
43	CONTAINS
44
45	SUBROUTINE obc_init
46	!!----------------------------------------------------------------------
47	!! * ROUTINE obc_init *
48	!!
49	!! ** Purpose : Initialization of the dynamics and tracer fields at
50	!! the open boundaries.
51	!!
52	!! ** Method : initialization of open boundary variables
53	!! (u, v) over 3 time step and 3 rows
54	!! (t, s) over 2 time step and 2 rows
55	!! if ln_rstart = .FALSE. : no restart, fields set to zero
56	!! if ln_rstart = .TRUE. : restart, fields are read in a file
57	!! if rdpxxx = 0 then lfbc is set true for this boundary.
58	!!
59	!! ** Input : restart.obc file, restart file for open boundaries
60	!!----------------------------------------------------------------------
61	USE obcrst, ONLY : obc_rst_read ! Make obc_rst_read routine available
62	!!
63	INTEGER :: ji, jj, istop , inumfbc
64	INTEGER, DIMENSION(4) :: icorner
65	REAL(wp), DIMENSION(2) :: ztestmask
66	!!
67	NAMELIST/namobc/ rn_dpein, rn_dpwin, rn_dpnin, rn_dpsin, &
68	& rn_dpeob, rn_dpwob, rn_dpnob, rn_dpsob, &
69	& rn_volemp, nn_obcdta, cn_obcdta, &
70	& ln_obc_clim, ln_vol_cst, ln_obc_fla
71	!!----------------------------------------------------------------------
72
73	REWIND( numnam ) ! Namelist namobc : open boundaries
74	READ ( numnam, namobc )
75
76	! convert DOCTOR namelist name into the OLD names
77	nobc_dta = nn_obcdta
78	cffile = cn_obcdta
79	rdpein = rn_dpein
80	rdpwin = rn_dpwin
81	rdpsin = rn_dpsin
82	rdpnin = rn_dpnin
83	rdpeob = rn_dpeob
84	rdpwob = rn_dpwob
85	rdpsob = rn_dpsob
86	rdpnob = rn_dpnob
87	volemp = rn_volemp
88
89	! ! allocate obc arrays
90	IF( obc_oce_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'obc_init : unable to allocate obc_oce arrays' )
91	IF( obc_dta_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'obc_init : unable to allocate obc_dta arrays' )
92
93	! By security we set rdpxin and rdpxob respectively to 1. and 15. if the corresponding OBC is not activated
94	IF( .NOT.lp_obc_east ) THEN ; rdpein = 1. ; rdpeob = 15. ; END IF
95	IF( .NOT.lp_obc_west ) THEN ; rdpwin = 1. ; rdpwob = 15. ; END IF
96	IF( .NOT.lp_obc_north ) THEN ; rdpnin = 1. ; rdpnob = 15. ; END IF
97	IF( .NOT.lp_obc_south ) THEN ; rdpsin = 1. ; rdpsob = 15. ; END IF
98
99	! number of open boudaries and open boundary indicators
100	nbobc = 0
101	IF( lp_obc_east ) nbobc = nbobc + 1
102	IF( lp_obc_west ) nbobc = nbobc + 1
103	IF( lp_obc_north ) nbobc = nbobc + 1
104	IF( lp_obc_south ) nbobc = nbobc + 1
105
106	IF(lwp) WRITE(numout,*)
107	IF(lwp) WRITE(numout,*) 'obc_init : initialization of open boundaries'
108	IF(lwp) WRITE(numout,*) '~~~~~~~~'
109	IF(lwp) WRITE(numout,*) ' Number of open boundaries nbobc = ', nbobc
110	IF(lwp) WRITE(numout,*)
111
112	! control prints
113	IF(lwp) WRITE(numout,*) ' Namelist namobc'
114	IF(lwp) WRITE(numout,*) ' data in file (=1) or initial state used (=0) nn_obcdta = ', nn_obcdta
115	IF(lwp) WRITE(numout,*) ' climatology (true) or not ln_obc_clim = ', ln_obc_clim
116	IF(lwp) WRITE(numout,*) ' vol_cst (true) or not: ln_vol_cst = ', ln_vol_cst
117	IF(lwp) WRITE(numout,*) ' '
118	IF(lwp) WRITE(numout,*) ' WARNING '
119	IF(lwp) WRITE(numout,*) ' Flather"s algorithm is applied with explicit free surface scheme '
120	IF(lwp) WRITE(numout,*) ' or with free surface time-splitting scheme '
121	IF(lwp) WRITE(numout,*) ' Nor radiation neither relaxation is allowed with explicit free surface scheme: '
122	IF(lwp) WRITE(numout,*) ' Radiation and/or relaxation is allowed with free surface time-splitting scheme '
123	IF(lwp) WRITE(numout,*) ' depending of the choice of rdpXin = rdpXob = 0. for open boundaries '
124	IF(lwp) WRITE(numout,*)
125	IF(lwp) WRITE(numout,*) ' For the filtered free surface case, '
126	IF(lwp) WRITE(numout,*) ' radiation, relaxation or presciption of data can be applied '
127	IF(lwp) WRITE(numout,*)
128
129	IF( lwp.AND.lp_obc_east ) THEN
130	WRITE(numout,*) ' East open boundary :'
131	WRITE(numout,*) ' i index jpieob = ', jpieob
132	WRITE(numout,*) ' damping time scale (days) rn_dpeob = ', rn_dpeob
133	WRITE(numout,*) ' damping time scale (days) rn_dpein = ', rn_dpein
134	ENDIF
135
136	IF( lwp.AND.lp_obc_west ) THEN
137	WRITE(numout,*) ' West open boundary :'
138	WRITE(numout,*) ' i index jpiwob = ', jpiwob
139	WRITE(numout,*) ' damping time scale (days) rn_dpwob = ', rn_dpwob
140	WRITE(numout,*) ' damping time scale (days) rn_dpwin = ', rn_dpwin
141	ENDIF
142
143	IF( lwp.AND.lp_obc_north ) THEN
144	WRITE(numout,*) ' North open boundary :'
145	WRITE(numout,*) ' j index jpjnob = ', jpjnob
146	WRITE(numout,*) ' damping time scale (days) rn_dpnob = ', rn_dpnob
147	WRITE(numout,*) ' damping time scale (days) rn_dpnin = ', rn_dpnin
148	ENDIF
149
150	IF( lwp.AND.lp_obc_south ) THEN
151	WRITE(numout,*) ' South open boundary :'
152	WRITE(numout,*) ' j index jpjsob = ', jpjsob
153	WRITE(numout,*) ' damping time scale (days) rn_dpsob = ', rn_dpsob
154	WRITE(numout,*) ' damping time scale (days) rn_dpsin = ', rn_dpsin
155	WRITE(numout,*)
156	ENDIF
157
158	IF( nbobc >= 2 .AND. jperio /= 0 ) &
159	& CALL ctl_stop( ' Cyclic or symmetric, and open boundary condition are not compatible' )
160
161	! 1. Initialisation of constants
162	! ------------------------------
163	! ... convert rdp$ob in seconds
164	! Fixed Bdy flag inbound outbound
165	lfbceast = .FALSE. ; rdpein = rdpein * rday ; rdpeob = rdpeob * rday
166	lfbcwest = .FALSE. ; rdpwin = rdpwin * rday ; rdpwob = rdpwob * rday
167	lfbcnorth = .FALSE. ; rdpnin = rdpnin * rday ; rdpnob = rdpnob * rday
168	lfbcsouth = .FALSE. ; rdpsin = rdpsin * rday ; rdpsob = rdpsob * rday
169	inumfbc = 0
170	! ... look for Fixed Boundaries (rdp = 0 )
171	! ... When specified, lbcxxx flags are set to TRUE and rdpxxx are set to
172	! ... a small arbitrary value, (to avoid division by zero further on).
173	! ... rdpxxx is not used anymore.
174	IF( lp_obc_east ) THEN
175	IF( (rdpein+rdpeob) == 0 ) THEN
176	lfbceast = .TRUE. ; rdpein = 1e-3 ; rdpeob = 1e-3
177	inumfbc = inumfbc+1
178	ELSEIF ( (rdpein*rdpeob) == 0 ) THEN
179	CALL ctl_stop( 'obc_init : rn_dpein & rn_dpeob must be both zero or non zero' )
180	END IF
181	END IF
182
183	IF( lp_obc_west ) THEN
184	IF( (rdpwin + rdpwob) == 0 ) THEN
185	lfbcwest = .TRUE. ; rdpwin = 1e-3 ; rdpwob = 1e-3
186	inumfbc = inumfbc+1
187	ELSEIF ( (rdpwin*rdpwob) == 0 ) THEN
188	CALL ctl_stop( 'obc_init : rn_dpwin & rn_dpwob must be both zero or non zero' )
189	END IF
190	END IF
191	IF( lp_obc_north ) THEN
192	IF( (rdpnin + rdpnob) == 0 ) THEN
193	lfbcnorth = .TRUE. ; rdpnin = 1e-3 ; rdpnob = 1e-3
194	inumfbc = inumfbc+1
195	ELSEIF ( (rdpnin*rdpnob) == 0 ) THEN
196	CALL ctl_stop( 'obc_init : rn_dpnin & rn_dpnob must be both zero or non zero' )
197	END IF
198	END IF
199	IF( lp_obc_south ) THEN
200	IF( (rdpsin + rdpsob) == 0 ) THEN
201	lfbcsouth = .TRUE. ; rdpsin = 1e-3 ; rdpsob = 1e-3
202	inumfbc = inumfbc+1
203	ELSEIF ( (rdpsin*rdpsob) == 0 ) THEN
204	CALL ctl_stop( 'obc_init : rn_dpsin & rn_dpsob must be both zero or non zero' )
205	END IF
206	END IF
207
208	! 2. Clever mpp indices for loops on the open boundaries.
209	! The loops will be performed only on the processors
210	! that contain a given open boundary.
211	! --------------------------------------------------------
212
213	IF( lp_obc_east ) THEN
214	! ... mpp initialization
215	nie0 = max( 1, min(jpieob - nimpp+1, jpi ) )
216	nie1 = max( 0, min(jpieob - nimpp+1, jpi - 1 ) )
217	nie0p1 = max( 1, min(jpieob+1 - nimpp+1, jpi ) )
218	nie1p1 = max( 0, min(jpieob+1 - nimpp+1, jpi - 1 ) )
219	nie0m1 = max( 1, min(jpieob-1 - nimpp+1, jpi ) )
220	nie1m1 = max( 0, min(jpieob-1 - nimpp+1, jpi - 1 ) )
221	nje0 = max( 2, min(jpjed - njmpp+1, jpj ) )
222	nje1 = max( 0, min(jpjef - njmpp+1, jpj - 1 ) )
223	nje0p1 = max( 1, min(jpjedp1 - njmpp+1, jpj ) )
224	nje0m1 = max( 1, min(jpjed - njmpp+1, jpj ) )
225	nje1m1 = max( 0, min(jpjefm1 - njmpp+1, jpj - 1 ) )
226	nje1m2 = max( 0, min(jpjefm1-1- njmpp+1, jpj - 1 ) )
227	IF(lwp) THEN
228	IF( lfbceast ) THEN
229	WRITE(numout,*)' '
230	WRITE(numout,*)' Specified East Open Boundary'
231	ELSE
232	WRITE(numout,*)' '
233	WRITE(numout,*)' Radiative East Open Boundary'
234	END IF
235	END IF
236	END IF
237
238	IF( lp_obc_west ) THEN
239	! ... mpp initialization
240	niw0 = max( 1, min(jpiwob - nimpp+1, jpi ) )
241	niw1 = max( 0, min(jpiwob - nimpp+1, jpi - 1 ) )
242	niw0p1 = max( 1, min(jpiwob+1 - nimpp+1, jpi ) )
243	niw1p1 = max( 0, min(jpiwob+1 - nimpp+1, jpi - 1 ) )
244	njw0 = max( 2, min(jpjwd - njmpp+1, jpj ) )
245	njw1 = max( 0, min(jpjwf - njmpp+1, jpj - 1 ) )
246	njw0p1 = max( 1, min(jpjwdp1 - njmpp+1, jpj ) )
247	njw0m1 = max( 1, min(jpjwd - njmpp+1, jpj ) )
248	njw1m1 = max( 0, min(jpjwfm1 - njmpp+1, jpj - 1 ) )
249	njw1m2 = max( 0, min(jpjwfm1-1- njmpp+1, jpj - 1 ) )
250	IF(lwp) THEN
251	IF( lfbcwest ) THEN
252	WRITE(numout,*)' '
253	WRITE(numout,*)' Specified West Open Boundary'
254	ELSE
255	WRITE(numout,*)' '
256	WRITE(numout,*)' Radiative West Open Boundary'
257	END IF
258	END IF
259	END IF
260
261	IF( lp_obc_north ) THEN
262	! ... mpp initialization
263	nin0 = max( 2, min(jpind - nimpp+1, jpi ) )
264	nin1 = max( 0, min(jpinf - nimpp+1, jpi - 1 ) )
265	nin0p1 = max( 1, min(jpindp1 - nimpp+1, jpi ) )
266	nin0m1 = max( 1, min(jpind - nimpp+1, jpi ) )
267	nin1m1 = max( 0, min(jpinfm1 - nimpp+1, jpi - 1 ) )
268	nin1m2 = max( 0, min(jpinfm1-1- nimpp+1, jpi - 1 ) )
269	njn0 = max( 1, min(jpjnob - njmpp+1, jpj ) )
270	njn1 = max( 0, min(jpjnob - njmpp+1, jpj - 1 ) )
271	njn0p1 = max( 1, min(jpjnob+1 - njmpp+1, jpj ) )
272	njn1p1 = max( 0, min(jpjnob+1 - njmpp+1, jpj - 1 ) )
273	njn0m1 = max( 1, min(jpjnob-1 - njmpp+1, jpj ) )
274	njn1m1 = max( 0, min(jpjnob-1 - njmpp+1, jpj - 1 ) )
275	IF(lwp) THEN
276	IF( lfbcnorth ) THEN
277	WRITE(numout,*)' '
278	WRITE(numout,*)' Specified North Open Boundary'
279	ELSE
280	WRITE(numout,*)' '
281	WRITE(numout,*)' Radiative North Open Boundary'
282	END IF
283	END IF
284	END IF
285
286	IF( lp_obc_south ) THEN
287	! ... mpp initialization
288	nis0 = max( 2, min(jpisd - nimpp+1, jpi ) )
289	nis1 = max( 0, min(jpisf - nimpp+1, jpi - 1 ) )
290	nis0p1 = max( 1, min(jpisdp1 - nimpp+1, jpi ) )
291	nis0m1 = max( 1, min(jpisd - nimpp+1, jpi ) )
292	nis1m1 = max( 0, min(jpisfm1 - nimpp+1, jpi - 1 ) )
293	nis1m2 = max( 0, min(jpisfm1-1- nimpp+1, jpi - 1 ) )
294	njs0 = max( 1, min(jpjsob - njmpp+1, jpj ) )
295	njs1 = max( 0, min(jpjsob - njmpp+1, jpj - 1 ) )
296	njs0p1 = max( 1, min(jpjsob+1 - njmpp+1, jpj ) )
297	njs1p1 = max( 0, min(jpjsob+1 - njmpp+1, jpj - 1 ) )
298	IF(lwp) THEN
299	IF( lfbcsouth ) THEN
300	WRITE(numout,*)' '
301	WRITE(numout,*)' Specified South Open Boundary'
302	ELSE
303	WRITE(numout,*)' '
304	WRITE(numout,*)' Radiative South Open Boundary'
305	END IF
306	END IF
307	END IF
308
309	! 3. mask correction for OBCs
310	! ---------------------------
311
312	IF( lp_obc_east ) THEN
313	!... (jpjed,jpjefm1),jpieob
314	bmask(nie0p1:nie1p1,nje0:nje1m1) = 0.e0
315
316	! ... initilization to zero
317	uemsk(:,:) = 0.e0 ; vemsk(:,:) = 0.e0 ; temsk(:,:) = 0.e0
318
319	! ... set 2D mask on East OBC, Vopt
320	DO ji = fs_nie0, fs_nie1
321	DO jj = nje0, nje1
322	uemsk(jj,:) = umask(ji, jj,:) * tmask_i(ji,jj) * tmask_i(ji+1,jj)
323	vemsk(jj,:) = vmask(ji+1,jj,:) * tmask_i(ji+1,jj)
324	temsk(jj,:) = tmask(ji+1,jj,:) * tmask_i(ji+1,jj)
325	END DO
326	END DO
327
328	END IF
329
330	IF( lp_obc_west ) THEN
331	! ... (jpjwd,jpjwfm1),jpiwob
332	bmask(niw0:niw1,njw0:njw1m1) = 0.e0
333
334	! ... initilization to zero
335	uwmsk(:,:) = 0.e0 ; vwmsk(:,:) = 0.e0 ; twmsk(:,:) = 0.e0
336
337	! ... set 2D mask on West OBC, Vopt
338	DO ji = fs_niw0, fs_niw1
339	DO jj = njw0, njw1
340	uwmsk(jj,:) = umask(ji,jj,:) * tmask_i(ji,jj) * tmask_i(ji+1,jj)
341	vwmsk(jj,:) = vmask(ji,jj,:) * tmask_i(ji,jj)
342	twmsk(jj,:) = tmask(ji,jj,:) * tmask_i(ji,jj)
343	END DO
344	END DO
345
346	END IF
347
348	IF( lp_obc_north ) THEN
349	! ... jpjnob,(jpind,jpisfm1)
350	bmask(nin0:nin1m1,njn0p1:njn1p1) = 0.e0
351
352	! ... initilization to zero
353	unmsk(:,:) = 0.e0 ; vnmsk(:,:) = 0.e0 ; tnmsk(:,:) = 0.e0
354
355	! ... set 2D mask on North OBC, Vopt
356	DO jj = fs_njn0, fs_njn1
357	DO ji = nin0, nin1
358	unmsk(ji,:) = umask(ji,jj+1,:) * tmask_i(ji,jj+1)
359	vnmsk(ji,:) = vmask(ji,jj ,:) * tmask_i(ji,jj) * tmask_i(ji,jj+1)
360	tnmsk(ji,:) = tmask(ji,jj+1,:) * tmask_i(ji,jj+1)
361	END DO
362	END DO
363
364	END IF
365
366	IF( lp_obc_south ) THEN
367	! ... jpjsob,(jpisd,jpisfm1)
368	bmask(nis0:nis1m1,njs0:njs1) = 0.e0
369
370	! ... initilization to zero
371	usmsk(:,:) = 0.e0 ; vsmsk(:,:) = 0.e0 ; tsmsk(:,:) = 0.e0
372
373	! ... set 2D mask on South OBC, Vopt
374	DO jj = fs_njs0, fs_njs1
375	DO ji = nis0, nis1
376	usmsk(ji,:) = umask(ji,jj,:) * tmask_i(ji,jj)
377	vsmsk(ji,:) = vmask(ji,jj,:) * tmask_i(ji,jj) * tmask_i(ji,jj+1)
378	tsmsk(ji,:) = tmask(ji,jj,:) * tmask_i(ji,jj)
379	END DO
380	END DO
381
382	END IF
383
384	! ... Initialize obcumask and obcvmask for the Force filtering
385	! boundary condition in dynspg_flt
386	obcumask(:,:) = umask(:,:,1)
387	obcvmask(:,:) = vmask(:,:,1)
388
389	! ... Initialize obctmsk on overlap region and obcs. This mask
390	! is used in obcvol.F90 to calculate cumulate flux E-P.
391	! obc Tracer point are outside the domain ( U/V obc points) ==> masked by obctmsk
392	! - no flux E-P on obcs and overlap region (jpreci = jprecj = 1)
393	obctmsk(:,:) = tmask_i(:,:)
394
395	IF( lp_obc_east ) THEN
396	! ... East obc Force filtering mask for the grad D
397	obcumask(nie0 :nie1 ,nje0p1:nje1m1) = 0.e0
398	obcvmask(nie0p1:nie1p1,nje0p1:nje1m1) = 0.e0
399	! ... set to 0 on East OBC
400	obctmsk(nie0p1:nie1p1,nje0:nje1) = 0.e0
401	END IF
402
403	IF( lp_obc_west ) THEN
404	! ... West obc Force filtering mask for the grad D
405	obcumask(niw0:niw1,njw0:njw1) = 0.e0
406	obcvmask(niw0:niw1,njw0:njw1) = 0.e0
407	! ... set to 0 on West OBC
408	obctmsk(niw0:niw1,njw0:njw1) = 0.e0
409	END IF
410
411	IF( lp_obc_north ) THEN
412	! ... North obc Force filtering mask for the grad D
413	obcumask(nin0p1:nin1m1,njn0p1:njn1p1) = 0.e0
414	obcvmask(nin0p1:nin1m1,njn0 :njn1 ) = 0.e0
415	! ... set to 0 on North OBC
416	obctmsk(nin0:nin1,njn0p1:njn1p1) = 0.e0
417	END IF
418
419	IF( lp_obc_south ) THEN
420	! ... South obc Force filtering mask for the grad D
421	obcumask(nis0p1:nis1m1,njs0:njs1) = 0.e0
422	obcvmask(nis0p1:nis1m1,njs0:njs1) = 0.e0
423	! ... set to 0 on South OBC
424	obctmsk(nis0:nis1,njs0:njs1) = 0.e0
425	END IF
426
427	! 3.1 Total lateral surface
428	! -------------------------
429	obcsurftot = 0.e0
430
431	IF( lp_obc_east ) THEN ! ... East open boundary lateral surface
432	DO ji = nie0, nie1
433	DO jj = 1, jpj
434	obcsurftot = obcsurftot+hu(ji,jj)e2u(ji,jj)uemsk(jj,1) * MAX(obctmsk(ji,jj),obctmsk(ji+1,jj) )
435	END DO
436	END DO
437	END IF
438
439	IF( lp_obc_west ) THEN ! ... West open boundary lateral surface
440	DO ji = niw0, niw1
441	DO jj = 1, jpj
442	obcsurftot = obcsurftot+hu(ji,jj)e2u(ji,jj)uwmsk(jj,1) * MAX(obctmsk(ji,jj),obctmsk(ji+1,jj) )
443	END DO
444	END DO
445	END IF
446
447	IF( lp_obc_north ) THEN ! ... North open boundary lateral surface
448	DO jj = njn0, njn1
449	DO ji = 1, jpi
450	obcsurftot = obcsurftot+hv(ji,jj)e1v(ji,jj)vnmsk(ji,1) * MAX(obctmsk(ji,jj),obctmsk(ji,jj+1) )
451	END DO
452	END DO
453	END IF
454
455	IF( lp_obc_south ) THEN ! ... South open boundary lateral surface
456	DO jj = njs0, njs1
457	DO ji = 1, jpi
458	obcsurftot = obcsurftot+hv(ji,jj)e1v(ji,jj)vsmsk(ji,1) * MAX(obctmsk(ji,jj),obctmsk(ji,jj+1) )
459	END DO
460	END DO
461	END IF
462
463	IF( lk_mpp ) CALL mpp_sum( obcsurftot ) ! sum over the global domain
464
465	! 5. Control print on mask
466	! The extremities of the open boundaries must be in land
467	! or else correspond to an "ocean corner" between two open boundaries.
468	! corner 1 is southwest, 2 is south east, 3 is northeast, 4 is northwest.
469	! --------------------------------------------------------------------------
470
471	icorner(:)=0
472
473	! ... control of the west boundary
474	IF( lp_obc_west ) THEN
475	IF( jpiwob < 2 .OR. jpiwob >= jpiglo-2 ) THEN
476	WRITE(ctmp1,*) ' jpiwob exceed ', jpiglo-2, 'or less than 2'
477	CALL ctl_stop( ctmp1 )
478	END IF
479	ztestmask(:)=0.
480	DO ji=niw0,niw1
481	IF( (njw0 + njmpp - 1) == jpjwd ) ztestmask(1)=ztestmask(1)+ tmask(ji,njw0,1)
482	IF( (njw1 + njmpp - 1) == jpjwf ) ztestmask(2)=ztestmask(2)+ tmask(ji,njw1,1)
483	END DO
484	IF( lk_mpp ) CALL mpp_sum( ztestmask, 2 ) ! sum over the global domain
485
486	IF( ztestmask(1) /= 0. ) icorner(1)=icorner(1)+1
487	IF( ztestmask(2) /= 0. ) icorner(4)=icorner(4)+1
488	END IF
489
490	! ... control of the east boundary
491	IF( lp_obc_east ) THEN
492	IF( jpieob < 4 .OR. jpieob >= jpiglo ) THEN
493	WRITE(ctmp1,*) ' jpieob exceed ', jpiglo, ' or less than 4'
494	CALL ctl_stop( ctmp1 )
495	END IF
496	ztestmask(:)=0.
497	DO ji=nie0p1,nie1p1
498	IF( (nje0 + njmpp - 1) == jpjed ) ztestmask(1)=ztestmask(1)+ tmask(ji,nje0,1)
499	IF( (nje1 + njmpp - 1) == jpjef ) ztestmask(2)=ztestmask(2)+ tmask(ji,nje1,1)
500	END DO
501	IF( lk_mpp ) CALL mpp_sum( ztestmask, 2 ) ! sum over the global domain
502
503	IF( ztestmask(1) /= 0. ) icorner(2)=icorner(2)+1
504	IF( ztestmask(2) /= 0. ) icorner(3)=icorner(3)+1
505	END IF
506
507	! ... control of the north boundary
508	IF( lp_obc_north ) THEN
509	IF( jpjnob < 4 .OR. jpjnob >= jpjglo ) THEN
510	WRITE(ctmp1,*) 'jpjnob exceed ', jpjglo, ' or less than 4'
511	CALL ctl_stop( ctmp1 )
512	END IF
513	ztestmask(:)=0.
514	DO jj=njn0p1,njn1p1
515	IF( (nin0 + nimpp - 1) == jpind ) ztestmask(1)=ztestmask(1)+ tmask(nin0,jj,1)
516	IF( (nin1 + nimpp - 1) == jpinf ) ztestmask(2)=ztestmask(2)+ tmask(nin1,jj,1)
517	END DO
518	IF( lk_mpp ) CALL mpp_sum( ztestmask, 2 ) ! sum over the global domain
519
520	IF( ztestmask(1) /= 0. ) icorner(4)=icorner(4)+1
521	IF( ztestmask(2) /= 0. ) icorner(3)=icorner(3)+1
522	END IF
523
524	! ... control of the south boundary
525	IF( lp_obc_south ) THEN
526	IF( jpjsob < 2 .OR. jpjsob >= jpjglo-2 ) THEN
527	WRITE(ctmp1,*) ' jpjsob exceed ', jpjglo-2, ' or less than 2'
528	CALL ctl_stop( ctmp1 )
529	END IF
530	ztestmask(:)=0.
531	DO jj=njs0,njs1
532	IF( (nis0 + nimpp - 1) == jpisd ) ztestmask(1)=ztestmask(1)+ tmask(nis0,jj,1)
533	IF( (nis1 + nimpp - 1) == jpisf ) ztestmask(2)=ztestmask(2)+ tmask(nis1,jj,1)
534	END DO
535	IF( lk_mpp ) CALL mpp_sum( ztestmask, 2 ) ! sum over the global domain
536
537	IF( ztestmask(1) /= 0. ) icorner(1)=icorner(1)+1
538	IF( ztestmask(2) /= 0. ) icorner(2)=icorner(2)+1
539	END IF
540
541	IF( icorner(1) == 2 ) THEN
542	IF(lwp) WRITE(numout,*)
543	IF(lwp) WRITE(numout,*) ' South West ocean corner, two open boudaries'
544	IF(lwp) WRITE(numout,*) ' ========== '
545	IF(lwp) WRITE(numout,*)
546	IF( jpisd /= jpiwob.OR.jpjsob /= jpjwd ) &
547	& CALL ctl_stop( ' Open boundaries do not fit, we stop' )
548
549	ELSE IF( icorner(1) == 1 ) THEN
550	CALL ctl_stop( ' Open boundaries do not fit at SW corner, we stop' )
551	END IF
552
553	IF( icorner(2) == 2 ) THEN
554	IF(lwp) WRITE(numout,*)
555	IF(lwp) WRITE(numout,*) ' South East ocean corner, two open boudaries'
556	IF(lwp) WRITE(numout,*) ' ========== '
557	IF(lwp) WRITE(numout,*)
558	IF( jpisf /= jpieob+1.OR.jpjsob /= jpjed ) &
559	& CALL ctl_stop( ' Open boundaries do not fit, we stop' )
560	ELSE IF( icorner(2) == 1 ) THEN
561	CALL ctl_stop( ' Open boundaries do not fit at SE corner, we stop' )
562	END IF
563
564	IF( icorner(3) == 2 ) THEN
565	IF(lwp) WRITE(numout,*)
566	IF(lwp) WRITE(numout,*) ' North East ocean corner, two open boudaries'
567	IF(lwp) WRITE(numout,*) ' ========== '
568	IF(lwp) WRITE(numout,*)
569	IF( jpinf /= jpieob+1 .OR. jpjnob+1 /= jpjef ) &
570	& CALL ctl_stop( ' Open boundaries do not fit, we stop' )
571	ELSE IF( icorner(3) == 1 ) THEN
572	CALL ctl_stop( ' Open boundaries do not fit at NE corner, we stop' )
573	END IF
574
575	IF( icorner(4) == 2 ) THEN
576	IF(lwp) WRITE(numout,*)
577	IF(lwp) WRITE(numout,*) ' North West ocean corner, two open boudaries'
578	IF(lwp) WRITE(numout,*) ' ========== '
579	IF(lwp) WRITE(numout,*)
580	IF( jpind /= jpiwob.OR.jpjnob+1 /= jpjwf ) &
581	& CALL ctl_stop( ' Open boundaries do not fit, we stop' )
582	ELSE IF( icorner(4) == 1 ) THEN
583	CALL ctl_stop( ' Open boundaries do not fit at NW corner, we stop' )
584	END IF
585
586	! 6. Initialization of open boundary variables (u, v, t, s)
587	! --------------------------------------------------------------
588	! only if at least one boundary is radiative
589	IF ( inumfbc < nbobc .AND. ln_rstart ) THEN
590	! Restart from restart.obc
591	CALL obc_rst_read
592	ELSE
593
594	! ! ... Initialization to zero of radiation arrays.
595	! ! Those have dimensions of local subdomains
596
597	uebnd(:,:,:,:) = 0.e0 ; unbnd(:,:,:,:) = 0.e0
598	vebnd(:,:,:,:) = 0.e0 ; vnbnd(:,:,:,:) = 0.e0
599	tebnd(:,:,:,:) = 0.e0 ; tnbnd(:,:,:,:) = 0.e0
600	sebnd(:,:,:,:) = 0.e0 ; snbnd(:,:,:,:) = 0.e0
601
602	uwbnd(:,:,:,:) = 0.e0 ; usbnd(:,:,:,:) = 0.e0
603	vwbnd(:,:,:,:) = 0.e0 ; vsbnd(:,:,:,:) = 0.e0
604	twbnd(:,:,:,:) = 0.e0 ; tsbnd(:,:,:,:) = 0.e0
605	swbnd(:,:,:,:) = 0.e0 ; ssbnd(:,:,:,:) = 0.e0
606
607	END IF
608
609	! 7. Control print
610	! -----------------------------------------------------------------
611
612	! ... control of the east boundary
613	IF( lp_obc_east ) THEN
614	istop = 0
615	IF( jpieob < 4 .OR. jpieob >= jpiglo ) THEN
616	IF(lwp) WRITE(numout,cform_err)
617	IF(lwp) WRITE(numout,*) ' jpieob exceed ', jpim1, ' or less than 4'
618	istop = istop + 1
619	END IF
620
621	IF( lk_mpp ) THEN
622	! ...
623	IF( nimpp > jpieob-5) THEN
624	IF(lwp) WRITE(numout,cform_err)
625	IF(lwp) WRITE(numout,*) ' A sub-domain is too close to the East OBC'
626	IF(lwp) WRITE(numout,*) ' nimpp must be < jpieob-5'
627	istop = istop + 1
628	ENDIF
629	ELSE
630
631	! ... stop if e r r o r (s) detected
632	IF( istop /= 0 ) THEN
633	WRITE(ctmp1,*) istop,' obcini : E R R O R (S) detected : stop'
634	CALL ctl_stop( ctmp1 )
635	ENDIF
636	ENDIF
637	ENDIF
638
639	! ... control of the west boundary
640	IF( lp_obc_west ) THEN
641	istop = 0
642	IF( jpiwob < 2 .OR. jpiwob >= jpiglo ) THEN
643	IF(lwp) WRITE(numout,cform_err)
644	IF(lwp) WRITE(numout,*) ' jpiwob exceed ', jpim1, ' or less than 2'
645	istop = istop + 1
646	END IF
647
648	IF( lk_mpp ) THEN
649	IF( (nimpp < jpiwob+5) .AND. (nimpp > 1) ) THEN
650	IF(lwp) WRITE(numout,cform_err)
651	IF(lwp) WRITE(numout,*) ' A sub-domain is too close to the West OBC'
652	IF(lwp) WRITE(numout,*) ' nimpp must be > jpiwob-5 or =1'
653	istop = istop + 1
654	ENDIF
655	ELSE
656
657	! ... stop if e r r o r (s) detected
658	IF( istop /= 0 ) THEN
659	WRITE(ctmp1,*) istop,' obcini : E R R O R (S) detected : stop'
660	CALL ctl_stop( ctmp1 )
661	ENDIF
662	ENDIF
663	ENDIF
664
665	! control of the north boundary
666	IF( lp_obc_north ) THEN
667	istop = 0
668	IF( jpjnob < 4 .OR. jpjnob >= jpjglo ) THEN
669	IF(lwp) WRITE(numout,cform_err)
670	IF(lwp) WRITE(numout,*) ' jpjnob exceed ', jpjm1,' or less than 4'
671	istop = istop + 1
672	END IF
673
674	IF( lk_mpp ) THEN
675	IF( njmpp > jpjnob-5) THEN
676	IF(lwp) WRITE(numout,cform_err)
677	IF(lwp) WRITE(numout,*) ' A sub-domain is too close to the North OBC'
678	IF(lwp) WRITE(numout,*) ' njmpp must be < jpjnob-5'
679	istop = istop + 1
680	ENDIF
681	ELSE
682
683	! ... stop if e r r o r (s) detected
684	IF( istop /= 0 ) THEN
685	WRITE(ctmp1,*) istop,' obcini : E R R O R (S) detected : stop'
686	CALL ctl_stop( ctmp1 )
687	ENDIF
688	ENDIF
689	ENDIF
690
691	! control of the south boundary
692	IF( lp_obc_south ) THEN
693	istop = 0
694	IF( jpjsob < 2 .OR. jpjsob >= jpjglo ) THEN
695	IF(lwp) WRITE(numout,cform_err)
696	IF(lwp) WRITE(numout,*) ' jpjsob exceed ', jpjm1,' or less than 2'
697	istop = istop + 1
698	END IF
699
700	IF( lk_mpp ) THEN
701	IF( (njmpp < jpjsob+5) .AND. (njmpp > 1) ) THEN
702	IF(lwp) WRITE(numout,cform_err)
703	IF(lwp) WRITE(numout,*) ' A sub-domain is too close to the South OBC'
704	IF(lwp) WRITE(numout,*) ' njmpp must be > jpjsob+5 or =1'
705	istop = istop + 1
706	ENDIF
707	ELSE
708
709	! ... stop if e r r o r (s) detected
710	IF( istop /= 0 ) THEN
711	WRITE(ctmp1,*) istop,' obcini : E R R O R (S) detected : stop'
712	CALL ctl_stop( ctmp1 )
713	ENDIF
714	ENDIF
715	ENDIF
716
717	END SUBROUTINE obc_init
718
719	#else
720	!!---------------------------------------------------------------------------------
721	!! Dummy module NO open boundaries
722	!!---------------------------------------------------------------------------------
723	CONTAINS
724	SUBROUTINE obc_init ! Dummy routine
725	END SUBROUTINE obc_init
726	#endif
727
728	!!=================================================================================
729	END MODULE obcini

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

Download in other formats: