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

Contents of /MITgcm_contrib/SOSE/BoxAdj/code_ad/ctrl_getobcsn.F

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


Revision 1.2 - (show annotations) (download)
Wed Jan 19 12:39:13 2011 UTC (14 years, 11 months ago) by mmazloff
Branch: MAIN
Changes since 1.1: +17 -5 lines
account for partial cells

1
2 #include "CTRL_CPPOPTIONS.h"
3 #ifdef ALLOW_OBCS
4 # include "OBCS_OPTIONS.h"
5 #endif
6
7
8 subroutine ctrl_getobcsn(
9 I mytime,
10 I myiter,
11 I mythid
12 & )
13
14 c ==================================================================
15 c SUBROUTINE ctrl_getobcsn
16 c ==================================================================
17 c
18 c o Get northern 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_getobcsn
26 c ==================================================================
27
28 implicit none
29
30 #ifdef ALLOW_OBCSN_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 ilobcsn
59 integer iobcs
60 integer jp1
61
62 _RL dummy
63 _RL obcsnfac
64 logical obcsnfirst
65 logical obcsnchanged
66 integer obcsncount0
67 integer obcsncount1
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) fnameobcsn
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
93 jmin = 1
94 jmax = sny
95 imin = 1
96 imax = snx
97 jp1 = 0
98
99 c-- Now, read the control vector.
100 doglobalread = .false.
101 ladinit = .false.
102
103 if (optimcycle .ge. 0) then
104 ilobcsn=ilnblnk( xx_obcsn_file )
105 write(fnameobcsn(1:80),'(2a,i10.10)')
106 & xx_obcsn_file(1:ilobcsn), '.', optimcycle
107 endif
108
109 c-- Get the counters, flags, and the interpolation factor.
110 call ctrl_get_gen_rec(
111 I xx_obcsnstartdate, xx_obcsnperiod,
112 O obcsnfac, obcsnfirst, obcsnchanged,
113 O obcsncount0,obcsncount1,
114 I mytime, myiter, mythid )
115 cih
116 do iobcs = 1,nobcs
117 if ( obcsnfirst ) then
118 call active_read_xz( fnameobcsn, tmpfldxz,
119 & (obcsncount0-1)*nobcs+iobcs,
120 & doglobalread, ladinit, optimcycle,
121 & mythid, xx_obcsn_dummy )
122 cnma
123 #ifdef CTRL_OBCS_MODES
124 cih
125 if ( optimcycle .ge. 0) then
126 cih If normal velocity
127 if (iobcs .eq. 3) then
128 cih Begin loop over x-points.
129 do bj = jtlo,jthi
130 do bi = itlo, ithi
131 do i = imin,imax
132 cih If open boundary.
133 if ( OB_Jn(i,bi,bj) .ne. 0. ) then
134 j = OB_Jn(i,bi,bj)
135 cih Balance net transport.
136 cih Barotropic velocities are stored in the first level.
137 print*,'Correct vbaro: Apply shiftvel(1)'
138 print*,shiftvel(1),i,j
139 tmpfldxz(i,1,bi,bj) =
140 & tmpfldxz(i,1,bi,bj) + shiftvel(1)
141 cih Determine number of open vertical layers.
142 nz = 0
143 do k = 1,Nr
144 nz = nz + maskS(i,j+jp1,k,bi,bj)
145 end do
146 cih Compute absolute velocities from the barotropic-baroclinic modes.
147 do k = 1,Nr
148 if (k.le.nz) then
149 stmp = 0.
150 do nk = 1,nz
151 stmp = stmp +
152 & modesv(k,nk,nz)*tmpfldxz(i,nk,bi,bj)
153 end do
154 tmpz(k,bi,bj) = stmp
155 else
156 tmpz(k,bi,bj) = 0.
157 end if
158 end do
159 do k = 1,Nr
160 tmpfldxz(i,k,bi,bj) = tmpz(k,bi,bj)
161 & *recip_hFacS(i,j+jp1,k,bi,bj)
162 end do
163 cih End if open boundary.
164 end if
165 cih End loop over x-points
166 end do
167 end do
168 end do
169 cih End if iobcs = 3.
170 end if
171 cih
172 cih
173 cih If tangential velocity
174 if (iobcs .eq. 4) then
175 cih Begin loop over x-points.
176 do bj = jtlo,jthi
177 do bi = itlo, ithi
178 do i = imin,imax
179 cih If open boundary.
180 if ( OB_Jn(i,bi,bj) .ne. 0. ) then
181 j = OB_Jn(i,bi,bj)
182 cih Determine number of open vertical layers.
183 nz = 0
184 do k = 1,Nr
185 nz = nz + maskW(i,j,k,bi,bj)
186 end do
187 cih Compute absolute velocities from the barotropic-baroclinic modes.
188 do k = 1,Nr
189 if (k.le.nz) then
190 stmp = 0.
191 do nk = 1,nz
192 stmp = stmp +
193 & modesv(k,nk,nz)*tmpfldxz(i,nk,bi,bj)
194 end do
195 tmpz(k,bi,bj) = stmp
196 else
197 tmpz(k,bi,bj) = 0.
198 end if
199 end do
200 do k = 1,Nr
201 tmpfldxz(i,k,bi,bj) = tmpz(k,bi,bj)
202 & *recip_hFacW(i,j,k,bi,bj)
203 end do
204 cih End open boundary.
205 end if
206 cih End loop over x-points
207 end do
208 end do
209 end do
210 cih End if iobcs = 4.
211 end if
212 cih End if optimcycle > 0 .
213 end if
214 cnma
215 #endif
216 cih
217 do bj = jtlo,jthi
218 do bi = itlo,ithi
219 do k = 1,nr
220 do i = imin,imax
221 xx_obcsn1(i,k,bi,bj,iobcs) = tmpfldxz (i,k,bi,bj)
222 enddo
223 enddo
224 enddo
225 enddo
226 cih End if obcsnfirst
227 endif
228 cih
229 if ( (obcsnfirst) .or. (obcsnchanged) ) then
230
231 do bj = jtlo,jthi
232 do bi = itlo,ithi
233 do k = 1,nr
234 do i = imin,imax
235 tmpfldxz(i,k,bi,bj) = xx_obcsn1(i,k,bi,bj,iobcs)
236 enddo
237 enddo
238 enddo
239 enddo
240 cih
241 call exf_swapffields_xz( tmpfldxz2, tmpfldxz, mythid)
242 cih
243 do bj = jtlo,jthi
244 do bi = itlo,ithi
245 do k = 1,nr
246 do i = imin,imax
247 xx_obcsn0(i,k,bi,bj,iobcs) = tmpfldxz2(i,k,bi,bj)
248 enddo
249 enddo
250 enddo
251 enddo
252 cih
253 call active_read_xz( fnameobcsn, tmpfldxz,
254 & (obcsncount1-1)*nobcs+iobcs,
255 & doglobalread, ladinit, optimcycle,
256 & mythid, xx_obcsn_dummy )
257 cnma
258 #ifdef CTRL_OBCS_MODES
259 cih
260 if ( optimcycle .ge. 0) then
261 cih If normal velocity
262 if (iobcs .eq. 3) then
263 cih Begin loop over x-points.
264 do bj = jtlo,jthi
265 do bi = itlo, ithi
266 do i = imin,imax
267 cih If open boundary.
268 if ( OB_Jn(i,bi,bj) .ne. 0. ) then
269 j = OB_Jn(i,bi,bj)
270 cih Balance net transport.
271 cih Barotropic velocities are stored in the first level.
272 print*,'Correct vbaro: Apply shiftvel(2)'
273 print*,shiftvel(2),i,j
274 tmpfldxz(i,1,bi,bj) =
275 & tmpfldxz(i,1,bi,bj) + shiftvel(2)
276 cih Determine number of open vertical layers.
277 nz = 0
278 do k = 1,Nr
279 nz = nz + maskS(i,j+jp1,k,bi,bj)
280 end do
281 cih Compute absolute velocities from the barotropic-baroclinic modes.
282 do k = 1,Nr
283 if (k.le.nz) then
284 stmp = 0.
285 do nk = 1,nz
286 stmp = stmp +
287 & modesv(k,nk,nz)*tmpfldxz(i,nk,bi,bj)
288 end do
289 tmpz(k,bi,bj) = stmp
290 else
291 tmpz(k,bi,bj) = 0.
292 end if
293 end do
294 do k = 1,Nr
295 tmpfldxz(i,k,bi,bj) = tmpz(k,bi,bj)
296 & *recip_hFacS(i,j+jp1,k,bi,bj)
297 end do
298 cih End if open boundary.
299 end if
300 cih End loop over x-points
301 end do
302 end do
303 end do
304 cih End if iobcs = 3.
305 end if
306 cih
307 cih
308 cih If tangential velocity
309 if (iobcs .eq. 4) then
310 cih Begin loop over x-points.
311 do bj = jtlo,jthi
312 do bi = itlo, ithi
313 do i = imin,imax
314 cih If open boundary.
315 if ( OB_Jn(i,bi,bj) .ne. 0. ) then
316 j = OB_Jn(i,bi,bj)
317 cih Determine number of open vertical layers.
318 nz = 0
319 do k = 1,Nr
320 nz = nz + maskW(i,j,k,bi,bj)
321 end do
322 cih Compute absolute velocities from the barotropic-baroclinic modes.
323 do k = 1,Nr
324 if (k.le.nz) then
325 stmp = 0.
326 do nk = 1,nz
327 stmp = stmp +
328 & modesv(k,nk,nz)*tmpfldxz(i,nk,bi,bj)
329 end do
330 tmpz(k,bi,bj) = stmp
331 else
332 tmpz(k,bi,bj) = 0.
333 end if
334 end do
335 do k = 1,Nr
336 tmpfldxz(i,k,bi,bj) = tmpz(k,bi,bj)
337 & *recip_hFacW(i,j,k,bi,bj)
338 end do
339 cih End if open boundary.
340 end if
341 cih End loop over x-points
342 end do
343 end do
344 end do
345 cih End if iobcs = 4.
346 end if
347 cih End if optimcycle > 0 .
348 end if
349 cnma
350 #endif
351 cih
352 do bj = jtlo,jthi
353 do bi = itlo,ithi
354 do k = 1,nr
355 do i = imin,imax
356 xx_obcsn1 (i,k,bi,bj,iobcs) = tmpfldxz (i,k,bi,bj)
357 enddo
358 enddo
359 enddo
360 enddo
361 cih End if obcsnfirst .or. obcschanged
362 endif
363 cih
364 c-- Add control to model variable.
365 do bj = jtlo,jthi
366 do bi = itlo,ithi
367 c-- Calculate mask for tracer cells (0 => land, 1 => water).
368 do k = 1,nr
369 do i = 1,snx
370 j = OB_Jn(I,bi,bj)
371 if (iobcs .EQ. 1) then
372 OBNt(i,k,bi,bj) = OBNt (i,k,bi,bj)
373 & + obcsnfac *xx_obcsn0(i,k,bi,bj,iobcs)
374 & + (1. _d 0 - obcsnfac)*xx_obcsn1(i,k,bi,bj,iobcs)
375 OBNt(i,k,bi,bj) = OBNt(i,k,bi,bj)
376 & *maskS(i,j+jp1,k,bi,bj)
377 cgg & *maskxz(i,k,bi,bj)
378 else if (iobcs .EQ. 2) then
379 OBNs(i,k,bi,bj) = OBNs (i,k,bi,bj)
380 & + obcsnfac *xx_obcsn0(i,k,bi,bj,iobcs)
381 & + (1. _d 0 - obcsnfac)*xx_obcsn1(i,k,bi,bj,iobcs)
382 OBNs(i,k,bi,bj) = OBNs(i,k,bi,bj)
383 & *maskS(i,j+jp1,k,bi,bj)
384 cgg & *maskxz(i,k,bi,bj)
385 else if (iobcs .EQ. 4) then
386 OBNu(i,k,bi,bj) = OBNu (i,k,bi,bj)
387 & + obcsnfac *xx_obcsn0(i,k,bi,bj,iobcs)
388 & + (1. _d 0 - obcsnfac)*xx_obcsn1(i,k,bi,bj,iobcs)
389 OBNu(i,k,bi,bj) = OBNu(i,k,bi,bj)
390 & *maskW(i,j,k,bi,bj)
391 cgg & *maskxz(i,k,bi,bj)
392 else if (iobcs .EQ. 3) then
393 OBNv(i,k,bi,bj) = OBNv (i,k,bi,bj)
394 & + obcsnfac *xx_obcsn0(i,k,bi,bj,iobcs)
395 & + (1. _d 0 - obcsnfac)*xx_obcsn1(i,k,bi,bj,iobcs)
396 OBNv(i,k,bi,bj) = OBNv(i,k,bi,bj)
397 & *maskS(i,j+jp1,k,bi,bj)
398 cgg & *maskxz(i,k,bi,bj)
399 endif
400 enddo
401 enddo
402 enddo
403 enddo
404
405 C-- End over iobcs loop
406 enddo
407
408 #else /* ALLOW_OBCSN_CONTROL undefined */
409
410 c == routine arguments ==
411
412 _RL mytime
413 integer myiter
414 integer mythid
415
416 c-- CPP flag ALLOW_OBCSN_CONTROL undefined.
417
418 #endif /* ALLOW_OBCSN_CONTROL */
419
420 end
421
422
423
424
425

  ViewVC Help
Powered by ViewVC 1.1.22