/[MITgcm]/MITgcm_contrib/SOSE/BoxAdj/code_ad/ctrl_getobcse.F
ViewVC logotype

Annotation of /MITgcm_contrib/SOSE/BoxAdj/code_ad/ctrl_getobcse.F

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


Revision 1.1 - (hide annotations) (download)
Tue Jan 18 19:33:08 2011 UTC (14 years, 11 months ago) by mmazloff
Branch: MAIN
Adding in simple box model with code for controlling obcs with modal decomposition.

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

  ViewVC Help
Powered by ViewVC 1.1.22