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

Contents of /MITgcm/pkg/ctrl/ctrl_volflux.F

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


Revision 1.6 - (show annotations) (download)
Thu Nov 6 22:05:08 2003 UTC (20 years, 7 months ago) by heimbach
Branch: MAIN
CVS Tags: checkpoint53f_post, checkpoint54a_pre, checkpoint55c_post, checkpoint53b_pre, checkpoint55e_post, checkpoint52l_pre, checkpoint52e_pre, hrcube4, checkpoint52n_post, checkpoint52j_post, checkpoint52e_post, checkpoint52d_pre, checkpoint53c_post, checkpoint53d_post, checkpoint55d_pre, checkpoint55j_post, checkpoint56b_post, checkpoint52j_pre, checkpoint54a_post, branch-netcdf, checkpoint52l_post, checkpoint55h_post, checkpoint52k_post, checkpoint52b_pre, checkpoint54b_post, checkpoint54d_post, checkpoint54e_post, checkpoint55b_post, checkpoint52m_post, checkpoint55, checkpoint53a_post, checkpoint55a_post, checkpoint54, checkpoint54f_post, checkpoint53b_post, checkpoint55g_post, checkpoint55f_post, checkpoint55i_post, checkpoint56, checkpoint53, checkpoint52, checkpoint52d_post, checkpoint52a_post, checkpoint52b_post, checkpoint53g_post, checkpoint52f_post, checkpoint52c_post, ecco_c52_e35, hrcube5, checkpoint52a_pre, checkpoint52i_post, checkpoint53d_pre, checkpoint54c_post, checkpoint52i_pre, checkpoint51u_post, checkpoint52h_pre, checkpoint52f_pre, hrcube_1, hrcube_2, hrcube_3, checkpoint56a_post, checkpoint55d_post
Branch point for: netcdf-sm0
Changes since 1.5: +0 -4 lines
o merging from ecco-branch
o cleaned some cross-dependencies and updated CPP options

1
2 #include "CTRL_CPPOPTIONS.h"
3 #ifdef ALLOW_OBCS
4 # include "OBCS_OPTIONS.h"
5 #endif
6
7 subroutine ctrl_volflux(
8 I obcsncount,
9 O sumarea,
10 O sumflux, mythid
11 & )
12
13 c ==================================================================
14 c SUBROUTINE ctrl_volflux
15 c ==================================================================
16 c
17 c o calculate the o.b. volume flux due to control adjustments.
18 c o Assume the calendar is identical
19 c for all open boundaries. Need to save the barotropic adjustment
20 c velocity so it can be used in all ctrl_getobcs files.
21 c o WARNING: eastern boundary (not defined) filenames have been a
22 c problem in the past.
23 c
24 c - started G. Gebbie, MIT-WHOI, 15-June-2002
25 c ==================================================================
26 c SUBROUTINE ctrl_obcsvol
27 c ==================================================================
28
29 implicit none
30
31 c == global variables ==
32
33 #include "EEPARAMS.h"
34 #include "SIZE.h"
35 #include "PARAMS.h"
36 #include "GRID.h"
37 #include "DYNVARS.h"
38 #ifdef ALLOW_OBCS
39 # include "OBCS.h"
40 #endif
41
42 #include "ctrl.h"
43 #include "ctrl_dummy.h"
44 #include "optim.h"
45
46 c == routine arguments ==
47
48 integer obcsncount
49 _RL sumflux
50 _RL sumarea
51 integer mythid
52
53 c == local variables ==
54
55 integer bi,bj
56 integer i,j,k
57 integer itlo,ithi
58 integer jtlo,jthi
59 integer jmin,jmax
60 integer imin,imax
61 integer irec
62 integer il
63 integer iobcs
64 integer ip1
65 integer jp1
66 integer nrec
67 integer ilfld
68 integer igg
69
70 _RL tmpflux
71 _RL tmparea
72 _RL dummy
73 _RL gg
74 _RL tmpx
75 _RL tmpy
76 _RL obcsnfac
77 character*(80) fnamefldn
78 character*(80) fnameflds
79 character*(80) fnamefldw
80 character*(80) fnameflde
81
82 logical doglobalread
83 logical ladinit
84
85 #ifdef ECCO_VERBOSE
86 character*(MAX_LEN_MBUF) msgbuf
87 #endif
88
89 c == external functions ==
90
91 integer ilnblnk
92 external ilnblnk
93
94 c == end of interface ==
95
96 jtlo = mybylo(mythid)
97 jthi = mybyhi(mythid)
98 itlo = mybxlo(mythid)
99 ithi = mybxhi(mythid)
100 jmin = 1
101 jmax = sny
102 imin = 1
103 imax = snx
104
105 c-- Read tiled data.
106 doglobalread = .false.
107 ladinit = .false.
108
109 cgg Assume the number of records is the same for
110 cgg all boundaries. Needs to be improved someday.
111
112 #if (defined (ALLOW_OBCS_CONTROL) || \
113 defined (ALLOW_OBCS_COST_CONTRIBUTION))
114
115 tmpflux = 0. d 0
116 tmparea = 0. d 0
117 sumarea = 0. d 0
118 sumflux = 0. d 0
119
120 #ifdef ECCO_VERBOSE
121 _BEGIN_MASTER( mythid )
122 write(msgbuf,'(a)') ' '
123 call print_message( msgbuf, standardmessageunit,
124 & SQUEEZE_RIGHT , mythid)
125 write(msgbuf,'(a)') ' '
126 call print_message( msgbuf, standardmessageunit,
127 & SQUEEZE_RIGHT , mythid)
128 write(msgbuf,'(a,i9.8)')
129 & ' ctrl_volflux: number of records to process: ',nrec
130 call print_message( msgbuf, standardmessageunit,
131 & SQUEEZE_RIGHT , mythid)
132 write(msgbuf,'(a)') ' '
133 call print_message( msgbuf, standardmessageunit,
134 & SQUEEZE_RIGHT , mythid)
135 _END_MASTER( mythid )
136 #endif
137
138 if (optimcycle .ge. 0) then
139 c
140 #ifdef ALLOW_OBCSN_CONTROL
141 ilfld=ilnblnk( xx_obcsn_file )
142 write(fnamefldn(1:80),'(2a,i10.10)')
143 & xx_obcsn_file(1:ilfld),'.', optimcycle
144 #endif
145 #ifdef ALLOW_OBCSS_CONTROL
146 ilfld=ilnblnk( xx_obcss_file )
147 write(fnameflds(1:80),'(2a,i10.10)')
148 & xx_obcss_file(1:ilfld),'.',optimcycle
149 #endif
150 #ifdef ALLOW_OBCSW_CONTROL
151 ilfld=ilnblnk( xx_obcsw_file )
152 write(fnamefldw(1:80),'(2a,i10.10)')
153 & xx_obcsw_file(1:ilfld),'.',optimcycle
154 #endif
155 #ifdef ALLOW_OBCSE_CONTROL
156 ilfld=ilnblnk( xx_obcse_file )
157 write(fnameflde(1:80),'(2a,i10.10)')
158 & xx_obcse_file(1:ilfld),'.',optimcycle
159 #endif
160 c
161 endif
162
163 #ifdef ALLOW_OBCSN_CONTROL
164 jp1 = 0
165
166 call active_read_xz_loc(fnamefldn,tmpfldxz,
167 & (obcsncount-1)*nobcs+3, doglobalread,
168 & ladinit, optimcycle, mythid
169 & , xx_obcsn_dummy )
170
171 c-- Loop over this thread's tiles.
172 do bj = jtlo,jthi
173 do bi = itlo,ithi
174
175 tmpflux = 0. d0
176 tmparea = 0. d0
177
178 do k = 1, Nr
179 do i = imin,imax
180 j = Ob_Jn(I,bi,bj)
181 if (j.ne.0) then
182 cgg -- Alternatively I could read the maskobcs file. But this gives the same result.
183 if (maskS(i,j+jp1,k,bi,bj) .ne. 0.) then
184 cgg -- Do not let the corners contribute to the volume flux.
185 if(ob_iw(j,bi,bj).ne.i .and.ob_ie(j,bi,bj).ne.i)then
186 CGG -- Barotropic velocity stored in level 1.
187 tmpx = tmpfldxz(i,1,bi,bj)
188
189 cgg -- Pick the special point where barotropic velocity loses one degree of freedom.
190 cgg -- Add up the cross-sectional area of this column for later calculations.
191 if (ob_iw(j,bi,bj).eq.(i-1) .and.
192 & ob_iw(j,bi,bj).ne. 0) then
193 tmpx = 0.
194 tmparea = tmparea + delR(k) * dxg(i,j+jp1,bi,bj)
195 print*,'tmparea',tmparea
196 endif
197 cgg -- Positive is flux in.
198 tmpflux = tmpflux -tmpx*delR(k)*dxg(i,j+jp1,bi,bj)
199 endif
200 endif
201 endif
202 enddo
203 enddo
204
205 sumarea = sumarea+ tmparea
206 sumflux = sumflux + tmpflux
207 enddo
208 enddo
209 #endif
210
211 #ifdef ALLOW_OBCSS_CONTROL
212 jp1 = 1
213
214 call active_read_xz_loc(fnameflds,tmpfldxz,
215 & (obcsncount-1)*nobcs+3, doglobalread,
216 & ladinit, optimcycle, mythid
217 & , xx_obcss_dummy )
218
219 c-- Loop over this thread's tiles.
220 do bj = jtlo,jthi
221 do bi = itlo,ithi
222
223 tmpflux = 0. d 0
224 #ifndef ALLOW_OBCSN_CONTROL
225 tmparea = 0. d 0
226 #endif
227 do k = 1, Nr
228 do i = imin,imax
229 j = Ob_Js(I,bi,bj)
230 if (j .ne. 0) then
231 if (maskS(i,j+jp1,k,bi,bj) .ne. 0.) then
232 cgg -- Do not let the corners contribute to the volume flux.
233 if (ob_iw(j,bi,bj).ne.i.and.ob_ie(j,bi,bj).ne.i)then
234 tmpx = tmpfldxz(i,1,bi,bj)
235 #ifndef ALLOW_OBCSN_CONTROL
236 cgg -- Pick the special point where barotropic velocity loses one degree of freedom.
237 cgg -- Add up the cross-sectional area of this column for later calculations.
238 cgg -- This is just the backup case where the northern boundary does not exist.
239 cgg -- warning: never been tested.
240 if (ob_iw(j,bi,bj).eq.(i-1).and.
241 & ob_iw(j,bi,bj).ne. 0) then
242 tmpx = 0.
243 tmparea = tmparea + delR(k) * dxg(i,j+jp1,bi,bj)
244 print*,'tmparea',tmparea
245 endif
246 #endif
247 cgg -- Positive is flux in.
248 tmpflux = tmpflux +tmpx*delR(k)*dxg(i,j+jp1,bi,bj)
249 endif
250 endif
251 endif
252 enddo
253 enddo
254 #ifndef ALLOW_OBCSN_CONTROL
255 sumarea = sumarea+ tmparea
256 #endif
257 sumflux = sumflux + tmpflux
258 enddo
259 enddo
260 #endif
261
262 #ifdef ALLOW_OBCSW_CONTROL
263 ip1 = 1
264
265 call active_read_yz_loc( fnamefldw, tmpfldyz,
266 & (obcsncount-1)*nobcs+3, doglobalread,
267 & ladinit, optimcycle, mythid
268 & , xx_obcsw_dummy )
269
270 c-- Loop over this thread's tiles.
271 do bj = jtlo,jthi
272 do bi = itlo,ithi
273
274 c-- Determine the weights to be used.
275 tmpflux = 0. d 0
276 #ifndef ALLOW_OBCSN_CONTROL
277 #ifndef ALLOW_OBCSS_CONTROL
278 tmparea = 0. d 0
279 #endif
280 #endif
281 do k = 1, Nr
282 do j = jmin,jmax
283 i = ob_iw(j,bi,bj)
284 if ( i .ne. 0) then
285 if (maskW(i+ip1,j,k,bi,bj) .ne. 0.) then
286 cgg -- Do not let the corners contribute to the volume flux.
287 if (ob_jn(i,bi,bj).ne.j.and. ob_js(i,bi,bj).ne.j)then
288 tmpy = tmpfldyz(j,1,bi,bj)
289
290 #ifndef ALLOW_OBCSN_CONTROL
291 #ifndef ALLOW_OBCSS_CONTROL
292 cgg -- Pick the special point where barotropic velocity loses one degree of freedom.
293 cgg -- Add up the cross-sectional area of this column for later calculations.
294 cgg -- This is an untested backup case.
295 if (ob_jn(i,bi,bj).eq.(j+1) .and.
296 & ob_jn(i,bi,bj).ne. 0) then
297 tmpy = 0.
298 tmparea = tmparea + delR(k) * dyg(i+ip1,j,bi,bj)
299 print*,'tmparea',tmparea
300 endif
301 #endif
302 #endif
303 cgg -- Positive is flux in.
304 tmpflux = tmpflux + tmpy* delR(k)*dyg(i+ip1,j,bi,bj)
305 endif
306 endif
307 endif
308 enddo
309 enddo
310 #ifndef ALLOW_OBCSN_CONTROL
311 #ifndef ALLOW_OBCSS_CONTROL
312 sumarea =sumarea + tmparea
313 #endif
314 #endif
315 sumflux = sumflux + tmpflux
316 enddo
317 enddo
318 #endif
319
320 #ifdef ALLOW_OBCSE_CONTROL
321 ip1 = 0
322
323 call active_read_yz_loc( fnameflde, tmpfldyz,
324 (obcsncount-1)*nobcs+3, doglobalread,
325 ladinit, optimcycle, mythid
326 & , xx_obcse_dummy )
327
328 c-- Loop over this thread's tiles.
329 do bj = jtlo,jthi
330 do bi = itlo,ithi
331
332 c-- Determine the weights to be used.
333 tmpflux = 0. d 0
334
335 #ifndef ALLOW_OBCSN_CONTROL
336 #ifndef ALLOW_OBCSS_CONTROL
337 #ifndef ALLOW_OBCSW_CONTROL
338 tmparea = 0. d 0
339 #endif
340 #endif
341 #endif
342
343 do k = 1, Nr
344 do j = jmin,jmax
345 i = ob_ie(j,bi,bj)
346 if ( i .ne. 0) then
347 if (maskW(i+ip1,j,k,bi,bj) .ne. 0.) then
348 cgg -- Do not let the corners contribute to the volume flux.
349 if (ob_jn(i,bi,bj).ne.j .and.ob_js(i,bi,bj).ne.j)then
350 tmpy = tmpfldyz(j,1,bi,bj)
351
352 #ifndef ALLOW_OBCSN_CONTROL
353 #ifndef ALLOW_OBCSS_CONTROL
354 #ifndef ALLOW_OBCSW_CONTROL
355 cgg -- Pick the special point where barotropic velocity loses one degree of freedom.
356 cgg -- Add up the cross-sectional area of this column for later calculations.
357 cgg -- This is an untested backup case.
358 if (ob_jn(i,bi,bj).eq.(j+1) .and.
359 & ob_jn(i,bi,bj).ne. 0) then
360 tmpy = 0.
361 tmparea = tmparea + delR(k) * dyg(i+ip1,j,bi,bj)
362 print*,'tmparea',tmparea
363 endif
364 #endif
365 #endif
366 #endif
367
368 cgg -- Positive is flux in.
369 tmpflux = tmpflux - tmpy* delR(k)*dyg(i+ip1,j,bi,bj)
370 #ifndef ALLOW_OBCSN_CONTROL
371 #ifndef ALLOW_OBCSS_CONTROL
372 #ifndef ALLOW_OBCSW_CONTROL
373 tmparea = tmparea + delR(k) *dyg(i+ip1,j,bi,bj)
374 #endif
375 #endif
376 #endif
377 endif
378 endif
379 endif
380 enddo
381 enddo
382
383 #ifndef ALLOW_OBCSN_CONTROL
384 #ifndef ALLOW_OBCSS_CONTROL
385 #ifndef ALLOW_OBCSW_CONTROL
386 sumarea = sumarea+ tmparea
387 #endif
388 #endif
389 #endif
390 sumflux = sumflux + tmpflux
391 enddo
392 enddo
393 #endif
394
395 #endif
396
397 return
398 end
399
400
401
402
403
404
405

  ViewVC Help
Powered by ViewVC 1.1.22