/[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.2 - (hide annotations) (download)
Wed Jan 19 12:39:13 2011 UTC (14 years, 11 months ago) by mmazloff
Branch: MAIN
Changes since 1.1: +17 -6 lines
account for partial cells

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 mmazloff 1.2 cnma Here we go into the modes' space
123     #ifdef CTRL_OBCS_MODES
124 mmazloff 1.1 cih
125     if ( optimcycle .ge. 0) then
126     cih If normal velocity
127     if (iobcs .eq. 3) then
128     cih Begin loop over y-points.
129     do bj = jtlo,jthi
130     do bi = itlo, ithi
131     do j = jmin,jmax
132     cih If open boundary.
133     if ( OB_Iw(j,bi,bj) .ne. 0. ) then
134     i = OB_Iw(j,bi,bj)
135     cih Determine number of open vertical layers.
136     nz = 0
137     do k = 1,Nr
138     nz = nz + maskW(i+ip1,j,k,bi,bj)
139     end do
140     cih Compute absolute velocities from the barotropic-baroclinic modes.
141     do k = 1,Nr
142     if (k.le.nz) then
143     stmp = 0.
144     do nk = 1,nz
145     stmp = stmp +
146     & modesv(k,nk,nz)*tmpfldyz(j,nk,bi,bj)
147     end do
148     tmpz(k,bi,bj) = stmp
149     else
150     tmpz(k,bi,bj) = 0.
151     end if
152     end do
153     do k = 1,Nr
154 mmazloff 1.2 tmpfldyz(j,k,bi,bj) = tmpz(k,bi,bj)
155     & *recip_hFacW(i+ip1,j,k,bi,bj)
156 mmazloff 1.1 end do
157     cih End open boundary.
158     end if
159     cih End loop over x-points
160     end do
161     end do
162     end do
163     cih End if iobcs = 3.
164     end if
165     cih
166     cih If tangential velocity
167     if (iobcs .eq. 4) then
168     cih Begin loop over y-points.
169     do bj = jtlo,jthi
170     do bi = itlo, ithi
171     do j = jmin,jmax
172     cih If open boundary.
173     if ( OB_Iw(j,bi,bj) .ne. 0. ) then
174     i = OB_Iw(j,bi,bj)
175     cih Determine number of open vertical layers.
176     nz = 0
177     do k = 1,Nr
178     nz = nz + maskS(i,j,k,bi,bj)
179     end do
180     cih Compute absolute velocities from the barotropic-baroclinic modes.
181     do k = 1,Nr
182     if (k.le.nz) then
183     stmp = 0.
184     do nk = 1,nz
185     stmp = stmp +
186     & modesv(k,nk,nz)*tmpfldyz(j,nk,bi,bj)
187     end do
188     tmpz(k,bi,bj) = stmp
189     else
190     tmpz(k,bi,bj) = 0.
191     end if
192     end do
193     do k = 1,Nr
194 mmazloff 1.2 tmpfldyz(j,k,bi,bj) = tmpz(k,bi,bj)
195     & *recip_hFacS(i,j,k,bi,bj)
196 mmazloff 1.1 end do
197     cih End open boundary.
198     end if
199     cih End loop over x-points
200     end do
201     end do
202     end do
203     cih End if iobcs = 4.
204     end if
205     cih End if optimcycle > 0 .
206 mmazloff 1.2 end if
207     cnma
208     #endif
209 mmazloff 1.1 cih
210     do bj = jtlo,jthi
211     do bi = itlo,ithi
212     do k = 1,nr
213     do j = jmin,jmax
214     xx_obcsw1(j,k,bi,bj,iobcs) = tmpfldyz (j,k,bi,bj)
215     enddo
216     enddo
217     enddo
218     enddo
219     cih End if obcswfirst
220     endif
221     cih
222     if ( (obcswfirst) .or. (obcswchanged) ) then
223     cih
224     do bj = jtlo,jthi
225     do bi = itlo,ithi
226     do k = 1,nr
227     do j = jmin,jmax
228     tmpfldyz(j,k,bi,bj) = xx_obcsw1(j,k,bi,bj,iobcs)
229     enddo
230     enddo
231     enddo
232     enddo
233     cih
234     call exf_swapffields_yz( tmpfldyz2, tmpfldyz, mythid)
235     cih
236     do bj = jtlo,jthi
237     do bi = itlo,ithi
238     do k = 1,nr
239     do j = jmin,jmax
240     xx_obcsw0(j,k,bi,bj,iobcs) = tmpfldyz2(j,k,bi,bj)
241     enddo
242     enddo
243     enddo
244     enddo
245     cih
246     call active_read_yz( fnameobcsw, tmpfldyz,
247     & (obcswcount1-1)*nobcs+iobcs,
248     & doglobalread, ladinit, optimcycle,
249     & mythid, xx_obcsw_dummy )
250 mmazloff 1.2 cnma
251     #ifdef CTRL_OBCS_MODES
252 mmazloff 1.1 cih
253     if ( optimcycle .ge. 0) then
254     cih If normal velocity
255     if (iobcs .eq. 3) then
256     cih Begin loop over y-points.
257     do bj = jtlo,jthi
258     do bi = itlo, ithi
259     do j = jmin,jmax
260     cih If open boundary.
261     if ( OB_Iw(j,bi,bj) .ne. 0. ) then
262     i = OB_Iw(j,bi,bj)
263     cih Determine number of open vertical layers.
264     nz = 0
265     do k = 1,Nr
266     nz = nz + maskW(i+ip1,j,k,bi,bj)
267     end do
268     cih Compute absolute velocities from the barotropic-baroclinic modes.
269     do k = 1,Nr
270     if (k.le.nz) then
271     stmp = 0.
272     do nk = 1,nz
273     stmp = stmp +
274     & modesv(k,nk,nz)*tmpfldyz(j,nk,bi,bj)
275     end do
276     tmpz(k,bi,bj) = stmp
277     else
278     tmpz(k,bi,bj) = 0.
279     end if
280     end do
281     do k = 1,Nr
282 mmazloff 1.2 tmpfldyz(j,k,bi,bj) = tmpz(k,bi,bj)
283     & *recip_hFacW(i+ip1,j,k,bi,bj)
284 mmazloff 1.1 end do
285     cih End open boundary.
286     end if
287     cih End loop over x-points
288     end do
289     end do
290     end do
291     cih End if iobcs = 3.
292     end if
293     cih
294     cih If tangential velocity
295     if (iobcs .eq. 4) then
296     cih Begin loop over y-points.
297     do bj = jtlo,jthi
298     do bi = itlo, ithi
299     do j = jmin,jmax
300     cih If open boundary.
301     if ( OB_Iw(j,bi,bj) .ne. 0. ) then
302     i = OB_Iw(j,bi,bj)
303     cih Determine number of open vertical layers.
304     nz = 0
305     do k = 1,Nr
306     nz = nz + maskS(i,j,k,bi,bj)
307     end do
308     cih Compute absolute velocities from the barotropic-baroclinic modes.
309     do k = 1,Nr
310     if (k.le.nz) then
311     stmp = 0.
312     do nk = 1,nz
313     stmp = stmp +
314     & modesv(k,nk,nz)*tmpfldyz(j,nk,bi,bj)
315     end do
316     tmpz(k,bi,bj) = stmp
317     else
318     tmpz(k,bi,bj) = 0.
319     end if
320     end do
321     do k = 1,Nr
322 mmazloff 1.2 tmpfldyz(j,k,bi,bj) = tmpz(k,bi,bj)
323     & *recip_hFacS(i,j,k,bi,bj)
324 mmazloff 1.1 end do
325     cih End open boundary.
326     end if
327     cih End loop over x-points
328     end do
329     end do
330     end do
331     cih End if iobcs = 4.
332     end if
333     cih End if optimcycle > 0 .
334     end if
335 mmazloff 1.2 cnma
336     #endif
337 mmazloff 1.1 cih
338     do bj = jtlo,jthi
339     do bi = itlo,ithi
340     do k = 1,nr
341     do j = jmin,jmax
342     xx_obcsw1 (j,k,bi,bj,iobcs) = tmpfldyz (j,k,bi,bj)
343     enddo
344     enddo
345     enddo
346     enddo
347     cih End if obcswfirst .or. obcswchanged
348     endif
349    
350     c-- Add control to model variable.
351     do bj = jtlo, jthi
352     do bi = itlo, ithi
353     c-- Calculate mask for tracer cells (0 => land, 1 => water).
354     do k = 1,nr
355     do j = 1,sny
356     i = OB_Iw(j,bi,bj)
357     if (iobcs .EQ. 1) then
358     OBWt(j,k,bi,bj) = OBWt (j,k,bi,bj)
359     & + obcswfac *xx_obcsw0(j,k,bi,bj,iobcs)
360     & + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)
361     OBWt(j,k,bi,bj) = OBWt(j,k,bi,bj)
362     & *maskW(i+ip1,j,k,bi,bj)
363     else if (iobcs .EQ. 2) then
364     OBWs(j,k,bi,bj) = OBWs (j,k,bi,bj)
365     & + obcswfac *xx_obcsw0(j,k,bi,bj,iobcs)
366     & + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)
367     OBWs(j,k,bi,bj) = OBWs(j,k,bi,bj)
368     & *maskW(i+ip1,j,k,bi,bj)
369     else if (iobcs .EQ. 3) then
370     OBWu(j,k,bi,bj) = OBWu (j,k,bi,bj)
371     & + obcswfac *xx_obcsw0(j,k,bi,bj,iobcs)
372     & + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)
373     OBWu(j,k,bi,bj) = OBWu(j,k,bi,bj)
374     & *maskW(i+ip1,j,k,bi,bj)
375     else if (iobcs .EQ. 4) then
376     OBWv(j,k,bi,bj) = OBWv (j,k,bi,bj)
377     & + obcswfac *xx_obcsw0(j,k,bi,bj,iobcs)
378     & + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)
379     OBWv(j,k,bi,bj) = OBWv(j,k,bi,bj)
380     & *maskS(i,j,k,bi,bj)
381     endif
382     enddo
383     enddo
384     enddo
385     enddo
386    
387     C-- End over iobcs loop
388     enddo
389    
390     #else /* ALLOW_OBCSW_CONTROL undefined */
391    
392     c == routine arguments ==
393    
394     _RL mytime
395     integer myiter
396     integer mythid
397    
398     c-- CPP flag ALLOW_OBCSW_CONTROL undefined.
399    
400     #endif /* ALLOW_OBCSW_CONTROL */
401    
402     end
403    

  ViewVC Help
Powered by ViewVC 1.1.22