/[MITgcm]/MITgcm/pkg/ctrl/ctrl_getobcse.F
ViewVC logotype

Annotation of /MITgcm/pkg/ctrl/ctrl_getobcse.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.9 - (hide annotations) (download)
Mon Mar 7 09:24:10 2011 UTC (13 years, 3 months ago) by mlosch
Branch: MAIN
Changes since 1.8: +40 -39 lines
remove global declaration of tmpfldx/yz and tmpfldx/yz2 in order to
clean up and remove storage requirements a little. Fortunately,
we do no have any tests for the numerous cpp-flagged option of obcs
as control parameters so we will never know if this is an
improvement (but at least now things compile for the flags that I found)

1 mlosch 1.9 C $Header: /u/gcmpack/MITgcm/pkg/ctrl/ctrl_getobcse.F,v 1.8 2011/01/19 08:42:06 mlosch Exp $
2 jmc 1.7 C $Name: $
3 heimbach 1.2
4     #include "CTRL_CPPOPTIONS.h"
5     #ifdef ALLOW_OBCS
6     # include "OBCS_OPTIONS.h"
7     #endif
8    
9    
10     subroutine ctrl_getobcse(
11     I mytime,
12     I myiter,
13     I mythid
14     & )
15    
16     c ==================================================================
17     c SUBROUTINE ctrl_getobcse
18     c ==================================================================
19     c
20     c o Get eastern obc of the control vector and add it
21     c to dyn. fields
22     c
23     c started: heimbach@mit.edu, 29-Aug-2001
24     c
25     c ==================================================================
26     c SUBROUTINE ctrl_getobcse
27     c ==================================================================
28    
29     implicit none
30    
31     #ifdef ALLOW_OBCSE_CONTROL
32    
33     c == global variables ==
34    
35     #include "EEPARAMS.h"
36     #include "SIZE.h"
37     #include "PARAMS.h"
38     #include "GRID.h"
39     #include "OBCS.h"
40    
41     #include "ctrl.h"
42     #include "ctrl_dummy.h"
43     #include "optim.h"
44    
45     c == routine arguments ==
46    
47     _RL mytime
48     integer myiter
49     integer mythid
50    
51     c == local variables ==
52    
53     integer bi,bj
54     integer i,j,k
55     integer itlo,ithi
56     integer jtlo,jthi
57     integer jmin,jmax
58     integer imin,imax
59     integer ilobcse
60     integer iobcs
61    
62     _RL dummy
63     _RL obcsefac
64     logical obcsefirst
65     logical obcsechanged
66     integer obcsecount0
67     integer obcsecount1
68     integer ip1
69    
70     cgg _RL maskyz (1-oly:sny+oly,nr,nsx,nsy)
71 mlosch 1.9 _RL tmpfldyz (1-oly:sny+oly,nr,nsx,nsy)
72 heimbach 1.2
73     logical doglobalread
74     logical ladinit
75    
76     character*(80) fnameobcse
77    
78     cgg( Variables for splitting barotropic/baroclinic vels.
79     _RL ubaro
80     _RL utop
81     cgg)
82    
83     c == external functions ==
84    
85     integer ilnblnk
86     external ilnblnk
87    
88    
89     c == end of interface ==
90    
91     jtlo = mybylo(mythid)
92     jthi = mybyhi(mythid)
93     itlo = mybxlo(mythid)
94     ithi = mybxhi(mythid)
95     jmin = 1-oly
96     jmax = sny+oly
97     imin = 1-olx
98     imax = snx+olx
99     ip1 = 0
100    
101     cgg( Initialize variables for balancing volume flux.
102     ubaro = 0.d0
103     utop = 0.d0
104     cgg)
105    
106     c-- Now, read the control vector.
107     doglobalread = .false.
108     ladinit = .false.
109    
110     if (optimcycle .ge. 0) then
111     ilobcse=ilnblnk( xx_obcse_file )
112 jmc 1.7 write(fnameobcse(1:80),'(2a,i10.10)')
113 heimbach 1.2 & xx_obcse_file(1:ilobcse), '.', optimcycle
114     endif
115    
116     c-- Get the counters, flags, and the interpolation factor.
117     call ctrl_get_gen_rec(
118     I xx_obcsestartdate, xx_obcseperiod,
119     O obcsefac, obcsefirst, obcsechanged,
120     O obcsecount0,obcsecount1,
121     I mytime, myiter, mythid )
122    
123     do iobcs = 1,nobcs
124    
125     if ( obcsefirst ) then
126 jmc 1.7 call active_read_yz( fnameobcse, tmpfldyz,
127 heimbach 1.2 & (obcsecount0-1)*nobcs+iobcs,
128     & doglobalread, ladinit, optimcycle,
129     & mythid, xx_obcse_dummy )
130    
131 mlosch 1.5 #ifdef ALLOW_CTRL_OBCS_BALANCE
132    
133 jmc 1.7 if ( optimcycle .gt. 0) then
134 heimbach 1.2 if (iobcs .eq. 3) then
135     cgg Special attention is needed for the normal velocity.
136     cgg For the north, this is the v velocity, iobcs = 4.
137     cgg This is done on a columnwise basis here.
138     do bj = jtlo,jthi
139     do bi = itlo, ithi
140     do j = jmin,jmax
141     i = OB_Ie(J,bi,bj)
142    
143     cgg The barotropic velocity is stored in the level 1.
144     ubaro = tmpfldyz(j,1,bi,bj)
145     tmpfldyz(j,1,bi,bj) = 0.d0
146     utop = 0.d0
147    
148     do k = 1,Nr
149     cgg If cells are not full, this should be modified with hFac.
150 jmc 1.7 cgg
151 mlosch 1.9 cgg The xx field (tmpfldyz) does not contain the velocity at the
152 jmc 1.7 cgg surface level. This velocity is not independent; it must
153 heimbach 1.2 cgg exactly balance the volume flux, since we are dealing with
154 jmc 1.7 cgg the baroclinic velocity structure..
155 heimbach 1.2 utop = tmpfldyz(j,k,bi,bj)*
156 heimbach 1.4 & maskW(i+ip1,j,k,bi,bj) * delR(k) + utop
157 jmc 1.7 cgg Add the barotropic velocity component.
158 heimbach 1.2 if (maskW(i+ip1,j,k,bi,bj) .ne. 0.) then
159     tmpfldyz(j,k,bi,bj) = tmpfldyz(j,k,bi,bj)+ ubaro
160     endif
161     enddo
162     cgg Compute the baroclinic velocity at level 1. Should balance flux.
163 jmc 1.7 tmpfldyz(j,1,bi,bj) = tmpfldyz(j,1,bi,bj)
164 heimbach 1.4 & - utop / delR(1)
165 heimbach 1.2 enddo
166     enddo
167     enddo
168     endif
169     if (iobcs .eq. 4) then
170     cgg Special attention is needed for the normal velocity.
171     cgg For the north, this is the v velocity, iobcs = 4.
172     cgg This is done on a columnwise basis here.
173     do bj = jtlo,jthi
174     do bi = itlo, ithi
175     do j = jmin,jmax
176     i = OB_Ie(J,bi,bj)
177    
178     cgg The barotropic velocity is stored in the level 1.
179     ubaro = tmpfldyz(j,1,bi,bj)
180     tmpfldyz(j,1,bi,bj) = 0.d0
181     utop = 0.d0
182    
183     do k = 1,Nr
184     cgg If cells are not full, this should be modified with hFac.
185 jmc 1.7 cgg
186 mlosch 1.9 cgg The xx field (tmpfldyz) does not contain the velocity at the
187 jmc 1.7 cgg surface level. This velocity is not independent; it must
188 heimbach 1.2 cgg exactly balance the volume flux, since we are dealing with
189 jmc 1.7 cgg the baroclinic velocity structure..
190 heimbach 1.2 utop = tmpfldyz(j,k,bi,bj)*
191 heimbach 1.4 & maskS(i,j,k,bi,bj) * delR(k) + utop
192 jmc 1.7 cgg Add the barotropic velocity component.
193 heimbach 1.2 if (maskS(i,j,k,bi,bj) .ne. 0.) then
194     tmpfldyz(j,k,bi,bj) = tmpfldyz(j,k,bi,bj)+ ubaro
195     endif
196     enddo
197     cgg Compute the baroclinic velocity at level 1. Should balance flux.
198 jmc 1.7 tmpfldyz(j,1,bi,bj) = tmpfldyz(j,1,bi,bj)
199 heimbach 1.4 & - utop / delR(1)
200 heimbach 1.2 enddo
201     enddo
202     enddo
203     endif
204     endif
205    
206 mlosch 1.5 #endif /* ALLOW_CTRL_OBCS_BALANCE */
207    
208 heimbach 1.2 do bj = jtlo,jthi
209     do bi = itlo,ithi
210     do k = 1,nr
211     do j = jmin,jmax
212     xx_obcse1(j,k,bi,bj,iobcs) = tmpfldyz (j,k,bi,bj)
213     cgg & * maskyz (j,k,bi,bj)
214     enddo
215     enddo
216     enddo
217     enddo
218     endif
219    
220     if ( (obcsefirst) .or. (obcsechanged)) then
221 jmc 1.7
222 heimbach 1.2 do bj = jtlo,jthi
223 mlosch 1.8 do bi = itlo,ithi
224     do j = jmin,jmax
225     do k = 1,nr
226     xx_obcse0(j,k,bi,bj,iobcs) = xx_obcse1(j,k,bi,bj,iobcs)
227     tmpfldyz (j,k,bi,bj) = 0. _d 0
228     enddo
229 heimbach 1.2 enddo
230 mlosch 1.8 enddo
231 heimbach 1.2 enddo
232    
233 jmc 1.7 call active_read_yz( fnameobcse, tmpfldyz,
234 heimbach 1.2 & (obcsecount1-1)*nobcs+iobcs,
235     & doglobalread, ladinit, optimcycle,
236     & mythid, xx_obcse_dummy )
237    
238 mlosch 1.5 #ifdef ALLOW_CTRL_OBCS_BALANCE
239    
240 jmc 1.7 if ( optimcycle .gt. 0) then
241 heimbach 1.2 if (iobcs .eq. 3) then
242     cgg Special attention is needed for the normal velocity.
243     cgg For the north, this is the v velocity, iobcs = 4.
244     cgg This is done on a columnwise basis here.
245     do bj = jtlo,jthi
246     do bi = itlo, ithi
247     do j = jmin,jmax
248     i = OB_Ie(J,bi,bj)
249    
250     cgg The barotropic velocity is stored in the level 1.
251     ubaro = tmpfldyz(j,1,bi,bj)
252     tmpfldyz(j,1,bi,bj) = 0.d0
253     utop = 0.d0
254    
255     do k = 1,Nr
256     cgg If cells are not full, this should be modified with hFac.
257 jmc 1.7 cgg
258 mlosch 1.9 cgg The xx field (tmpfldyz) does not contain the velocity at the
259 jmc 1.7 cgg surface level. This velocity is not independent; it must
260 heimbach 1.2 cgg exactly balance the volume flux, since we are dealing with
261 jmc 1.7 cgg the baroclinic velocity structure..
262 heimbach 1.2 utop = tmpfldyz(j,k,bi,bj)*
263 heimbach 1.4 & maskW(i+ip1,j,k,bi,bj) * delR(k) + utop
264 jmc 1.7 cgg Add the barotropic velocity component.
265 heimbach 1.2 if (maskW(i+ip1,j,k,bi,bj) .ne. 0.) then
266     tmpfldyz(j,k,bi,bj) = tmpfldyz(j,k,bi,bj)+ ubaro
267     endif
268     enddo
269     cgg Compute the baroclinic velocity at level 1. Should balance flux.
270 jmc 1.7 tmpfldyz(j,1,bi,bj) = tmpfldyz(j,1,bi,bj)
271 heimbach 1.4 & - utop / delR(1)
272 heimbach 1.2 enddo
273     enddo
274     enddo
275     endif
276     if (iobcs .eq. 4) then
277     cgg Special attention is needed for the normal velocity.
278     cgg For the north, this is the v velocity, iobcs = 4.
279     cgg This is done on a columnwise basis here.
280     do bj = jtlo,jthi
281     do bi = itlo, ithi
282     do j = jmin,jmax
283     i = OB_Ie(J,bi,bj)
284    
285     cgg The barotropic velocity is stored in the level 1.
286     ubaro = tmpfldyz(j,1,bi,bj)
287     tmpfldyz(j,1,bi,bj) = 0.d0
288     utop = 0.d0
289    
290     do k = 1,Nr
291     cgg If cells are not full, this should be modified with hFac.
292 jmc 1.7 cgg
293 mlosch 1.9 cgg The xx field (tmpfldyz) does not contain the velocity at the
294 jmc 1.7 cgg surface level. This velocity is not independent; it must
295 heimbach 1.2 cgg exactly balance the volume flux, since we are dealing with
296 jmc 1.7 cgg the baroclinic velocity structure..
297 heimbach 1.2 utop = tmpfldyz(j,k,bi,bj)*
298 heimbach 1.4 & maskS(i,j,k,bi,bj) * delR(k) + utop
299 jmc 1.7 cgg Add the barotropic velocity component.
300 heimbach 1.2 if (maskS(i,j,k,bi,bj) .ne. 0.) then
301     tmpfldyz(j,k,bi,bj) = tmpfldyz(j,k,bi,bj)+ ubaro
302     endif
303     enddo
304     cgg Compute the baroclinic velocity at level 1. Should balance flux.
305 jmc 1.7 tmpfldyz(j,1,bi,bj) = tmpfldyz(j,1,bi,bj)
306 heimbach 1.4 & - utop / delR(1)
307 heimbach 1.2 enddo
308     enddo
309     enddo
310     endif
311     endif
312    
313 mlosch 1.5 #endif /* ALLOW_CTRL_OBCS_BALANCE */
314    
315 heimbach 1.2 do bj = jtlo,jthi
316     do bi = itlo,ithi
317     do k = 1,nr
318     do j = jmin,jmax
319     xx_obcse1 (j,k,bi,bj,iobcs) = tmpfldyz (j,k,bi,bj)
320     cgg & * maskyz (j,k,bi,bj)
321     enddo
322     enddo
323     enddo
324     enddo
325     endif
326    
327 mlosch 1.9 c-- Add control to model variable.
328 heimbach 1.2 do bj = jtlo,jthi
329 mlosch 1.9 do bi = itlo,ithi
330     c-- Calculate mask for tracer cells (0 => land, 1 => water).
331     do k = 1,nr
332     do j = 1,sny
333     i = OB_Ie(j,bi,bj)
334     if (iobcs .EQ. 1) then
335     OBEt(j,k,bi,bj) = OBEt(j,k,bi,bj)
336     & + obcsefac *xx_obcse0(j,k,bi,bj,iobcs)
337     & + (1. _d 0 - obcsefac)*xx_obcse1(j,k,bi,bj,iobcs)
338     OBEt(j,k,bi,bj) = OBEt(j,k,bi,bj)
339     & *maskW(i+ip1,j,k,bi,bj)
340     else if (iobcs .EQ. 2) then
341     OBEs(j,k,bi,bj) = OBEs(j,k,bi,bj)
342     & + obcsefac *xx_obcse0(j,k,bi,bj,iobcs)
343     & + (1. _d 0 - obcsefac)*xx_obcse1(j,k,bi,bj,iobcs)
344     OBEs(j,k,bi,bj) = OBEs(j,k,bi,bj)
345     & *maskW(i+ip1,j,k,bi,bj)
346     else if (iobcs .EQ. 3) then
347     OBEu(j,k,bi,bj) = OBEu(j,k,bi,bj)
348     & + obcsefac *xx_obcse0(j,k,bi,bj,iobcs)
349     & + (1. _d 0 - obcsefac)*xx_obcse1(j,k,bi,bj,iobcs)
350     OBEu(j,k,bi,bj) = OBEu(j,k,bi,bj)
351     & *maskW(i+ip1,j,k,bi,bj)
352     else if (iobcs .EQ. 4) then
353     OBEv(j,k,bi,bj) = OBEv(j,k,bi,bj)
354     & + obcsefac *xx_obcse0(j,k,bi,bj,iobcs)
355     & + (1. _d 0 - obcsefac)*xx_obcse1(j,k,bi,bj,iobcs)
356     OBEv(j,k,bi,bj) = OBEv(j,k,bi,bj)
357     & *maskS(i,j,k,bi,bj)
358     endif
359 heimbach 1.2 enddo
360 mlosch 1.9 enddo
361     enddo
362 heimbach 1.2 enddo
363 mlosch 1.9
364 heimbach 1.2 C-- End over iobcs loop
365     enddo
366    
367     #else /* ALLOW_OBCSE_CONTROL undefined */
368    
369     c == routine arguments ==
370    
371     _RL mytime
372     integer myiter
373     integer mythid
374    
375     c-- CPP flag ALLOW_OBCSE_CONTROL undefined.
376    
377     #endif /* ALLOW_OBCSE_CONTROL */
378    
379     end
380    

  ViewVC Help
Powered by ViewVC 1.1.22