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 trunk/NEMO/OPA_SRC/OBC – NEMO

Context Navigation

source: trunk/NEMO/OPA_SRC/OBC/obcini.F90 @ 1152

Last change on this file since 1152 was 1152, checked in by rblod, 16 years ago
Convert cvs header to svn Id, step II
Property svn:eol-style set to `native` Property svn:keywords set to `Id`
File size: 30.9 KB

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

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

Download in other formats: