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

Annotation of /MITgcm_contrib/SOSE/BoxAdj/code_ad/ctrl_getobcsw.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_getobcsw(
9     I mytime,
10     I myiter,
11     I mythid
12     & )
13    
14     c ==================================================================
15     c SUBROUTINE ctrl_getobcsw
16     c ==================================================================
17     c
18     c o Get western 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 modified: gebbie@mit.edu, 18-Mar-2003
24     c ==================================================================
25     c SUBROUTINE ctrl_getobcsw
26     c ==================================================================
27    
28     implicit none
29    
30     #ifdef ALLOW_OBCSW_CONTROL
31    
32     c == global variables ==
33    
34     #include "EEPARAMS.h"
35     #include "SIZE.h"
36     #include "PARAMS.h"
37     #include "GRID.h"
38     #include "OBCS.h"
39    
40     #include "ctrl.h"
41     #include "ctrl_dummy.h"
42     #include "optim.h"
43    
44     c == routine arguments ==
45    
46     _RL mytime
47     integer myiter
48     integer mythid
49    
50     c == local variables ==
51    
52     integer bi,bj
53     integer i,j,k
54     integer itlo,ithi
55     integer jtlo,jthi
56     integer jmin,jmax
57     integer imin,imax
58     integer ilobcsw
59     integer iobcs
60     integer ip1
61    
62     _RL dummy
63     _RL obcswfac
64     logical obcswfirst
65     logical obcswchanged
66     integer obcswcount0
67     integer obcswcount1
68     cih
69     Integer nk,nz
70     _RL tmpz (nr,nsx,nsy)
71     _RL stmp
72     character*(80) fnamein
73     cih
74     logical doglobalread
75     logical ladinit
76    
77     character*(80) fnameobcsw
78    
79    
80     c == external functions ==
81    
82     integer ilnblnk
83     external ilnblnk
84    
85    
86     c == end of interface ==
87    
88     jtlo = mybylo(mythid)
89     jthi = mybyhi(mythid)
90     itlo = mybxlo(mythid)
91     ithi = mybxhi(mythid)
92     jmin = 1-oly
93     jmax = sny+oly
94     imin = 1-olx
95     imax = snx+olx
96     ip1 = 1
97    
98    
99     c-- Now, read the control vector.
100     doglobalread = .false.
101     ladinit = .false.
102    
103     if (optimcycle .ge. 0) then
104     ilobcsw=ilnblnk( xx_obcsw_file )
105     write(fnameobcsw(1:80),'(2a,i10.10)')
106     & xx_obcsw_file(1:ilobcsw), '.', optimcycle
107     endif
108    
109     c-- Get the counters, flags, and the interpolation factor.
110     call ctrl_get_gen_rec(
111     I xx_obcswstartdate, xx_obcswperiod,
112     O obcswfac, obcswfirst, obcswchanged,
113     O obcswcount0,obcswcount1,
114     I mytime, myiter, mythid )
115     cih
116     do iobcs = 1,nobcs
117     if ( obcswfirst ) then
118     call active_read_yz( fnameobcsw, tmpfldyz,
119     & (obcswcount0-1)*nobcs+iobcs,
120     & doglobalread, ladinit, optimcycle,
121     & mythid, xx_obcsw_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_Iw(j,bi,bj) .ne. 0. ) then
132     i = OB_Iw(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_Iw(j,bi,bj) .ne. 0. ) then
171     i = OB_Iw(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_obcsw1(j,k,bi,bj,iobcs) = tmpfldyz (j,k,bi,bj)
209     enddo
210     enddo
211     enddo
212     enddo
213     cih End if obcswfirst
214     endif
215     cih
216     if ( (obcswfirst) .or. (obcswchanged) ) then
217     cih
218     do bj = jtlo,jthi
219     do bi = itlo,ithi
220     do k = 1,nr
221     do j = jmin,jmax
222     tmpfldyz(j,k,bi,bj) = xx_obcsw1(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 k = 1,nr
233     do j = jmin,jmax
234     xx_obcsw0(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( fnameobcsw, tmpfldyz,
241     & (obcswcount1-1)*nobcs+iobcs,
242     & doglobalread, ladinit, optimcycle,
243     & mythid, xx_obcsw_dummy )
244     cih
245     if ( optimcycle .ge. 0) then
246    
247     cih If normal velocity
248     if (iobcs .eq. 3) then
249     cih Begin loop over y-points.
250     do bj = jtlo,jthi
251     do bi = itlo, ithi
252     do j = jmin,jmax
253     cih If open boundary.
254     if ( OB_Iw(j,bi,bj) .ne. 0. ) then
255     i = OB_Iw(j,bi,bj)
256     cih Determine number of open vertical layers.
257     nz = 0
258     do k = 1,Nr
259     nz = nz + maskW(i+ip1,j,k,bi,bj)
260     end do
261     cih Compute absolute velocities from the barotropic-baroclinic modes.
262     do k = 1,Nr
263     if (k.le.nz) then
264     stmp = 0.
265     do nk = 1,nz
266     stmp = stmp +
267     & modesv(k,nk,nz)*tmpfldyz(j,nk,bi,bj)
268     end do
269     tmpz(k,bi,bj) = stmp
270     else
271     tmpz(k,bi,bj) = 0.
272     end if
273     end do
274     do k = 1,Nr
275     tmpfldyz(j,k,bi,bj) = tmpz(k,bi,bj)
276     end do
277     cih End open boundary.
278     end if
279     cih End loop over x-points
280     end do
281     end do
282     end do
283     cih End if iobcs = 3.
284     end if
285     cih
286     cih If tangential velocity
287     if (iobcs .eq. 4) then
288     cih Begin loop over y-points.
289     do bj = jtlo,jthi
290     do bi = itlo, ithi
291     do j = jmin,jmax
292     cih If open boundary.
293     if ( OB_Iw(j,bi,bj) .ne. 0. ) then
294     i = OB_Iw(j,bi,bj)
295     cih Determine number of open vertical layers.
296     nz = 0
297     do k = 1,Nr
298     nz = nz + maskS(i,j,k,bi,bj)
299     end do
300     cih Compute absolute velocities from the barotropic-baroclinic modes.
301     do k = 1,Nr
302     if (k.le.nz) then
303     stmp = 0.
304     do nk = 1,nz
305     stmp = stmp +
306     & modesv(k,nk,nz)*tmpfldyz(j,nk,bi,bj)
307     end do
308     tmpz(k,bi,bj) = stmp
309     else
310     tmpz(k,bi,bj) = 0.
311     end if
312     end do
313     do k = 1,Nr
314     tmpfldyz(j,k,bi,bj) = tmpz(k,bi,bj)
315     end do
316     cih End open boundary.
317     end if
318     cih End loop over x-points
319     end do
320     end do
321     end do
322     cih End if iobcs = 4.
323     end if
324     cih End if optimcycle > 0 .
325     end if
326     cih
327     do bj = jtlo,jthi
328     do bi = itlo,ithi
329     do k = 1,nr
330     do j = jmin,jmax
331     xx_obcsw1 (j,k,bi,bj,iobcs) = tmpfldyz (j,k,bi,bj)
332     enddo
333     enddo
334     enddo
335     enddo
336     cih End if obcswfirst .or. obcswchanged
337     endif
338    
339     c-- Add control to model variable.
340     do bj = jtlo, jthi
341     do bi = itlo, ithi
342     c-- Calculate mask for tracer cells (0 => land, 1 => water).
343     do k = 1,nr
344     do j = 1,sny
345     i = OB_Iw(j,bi,bj)
346     if (iobcs .EQ. 1) then
347     OBWt(j,k,bi,bj) = OBWt (j,k,bi,bj)
348     & + obcswfac *xx_obcsw0(j,k,bi,bj,iobcs)
349     & + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)
350     OBWt(j,k,bi,bj) = OBWt(j,k,bi,bj)
351     & *maskW(i+ip1,j,k,bi,bj)
352     else if (iobcs .EQ. 2) then
353     OBWs(j,k,bi,bj) = OBWs (j,k,bi,bj)
354     & + obcswfac *xx_obcsw0(j,k,bi,bj,iobcs)
355     & + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)
356     OBWs(j,k,bi,bj) = OBWs(j,k,bi,bj)
357     & *maskW(i+ip1,j,k,bi,bj)
358     else if (iobcs .EQ. 3) then
359     OBWu(j,k,bi,bj) = OBWu (j,k,bi,bj)
360     & + obcswfac *xx_obcsw0(j,k,bi,bj,iobcs)
361     & + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)
362     OBWu(j,k,bi,bj) = OBWu(j,k,bi,bj)
363     & *maskW(i+ip1,j,k,bi,bj)
364     else if (iobcs .EQ. 4) then
365     OBWv(j,k,bi,bj) = OBWv (j,k,bi,bj)
366     & + obcswfac *xx_obcsw0(j,k,bi,bj,iobcs)
367     & + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)
368     OBWv(j,k,bi,bj) = OBWv(j,k,bi,bj)
369     & *maskS(i,j,k,bi,bj)
370     endif
371     enddo
372     enddo
373     enddo
374     enddo
375    
376     C-- End over iobcs loop
377     enddo
378    
379     #else /* ALLOW_OBCSW_CONTROL undefined */
380    
381     c == routine arguments ==
382    
383     _RL mytime
384     integer myiter
385     integer mythid
386    
387     c-- CPP flag ALLOW_OBCSW_CONTROL undefined.
388    
389     #endif /* ALLOW_OBCSW_CONTROL */
390    
391     end
392    

  ViewVC Help
Powered by ViewVC 1.1.22