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

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

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


Revision 1.4 - (show annotations) (download)
Sun Oct 26 00:58:03 2003 UTC (20 years, 7 months ago) by heimbach
Branch: MAIN
CVS Tags: checkpoint52l_pre, hrcube4, hrcube5, checkpoint52d_pre, checkpoint56b_post, checkpoint52j_pre, checkpoint54d_post, checkpoint54e_post, checkpoint52l_post, checkpoint52k_post, checkpoint55, checkpoint54, checkpoint56, checkpoint53, checkpoint52, checkpoint52f_post, checkpoint54f_post, checkpoint51t_post, checkpoint55i_post, checkpoint52i_pre, hrcube_1, hrcube_2, hrcube_3, checkpoint51s_post, checkpoint55c_post, checkpoint52e_pre, checkpoint52e_post, checkpoint53d_post, checkpoint52b_pre, checkpoint54b_post, checkpoint52m_post, checkpoint55g_post, checkpoint51q_post, checkpoint52b_post, checkpoint52c_post, checkpoint52f_pre, checkpoint55d_post, checkpoint54a_pre, checkpoint53c_post, checkpoint55d_pre, checkpoint55j_post, checkpoint54a_post, checkpoint55h_post, checkpoint51r_post, checkpoint55b_post, checkpoint53a_post, checkpoint55f_post, checkpoint52d_post, checkpoint53g_post, checkpoint52a_pre, checkpoint52i_post, checkpoint52h_pre, checkpoint56a_post, checkpoint53f_post, checkpoint52j_post, branch-netcdf, checkpoint52n_post, checkpoint53b_pre, checkpoint55a_post, checkpoint51o_post, checkpoint53b_post, checkpoint52a_post, ecco_c52_e35, checkpoint53d_pre, checkpoint55e_post, checkpoint54c_post, checkpoint51p_post, checkpoint51u_post
Branch point for: branch-nonh, netcdf-sm0
Changes since 1.3: +8 -8 lines
Replacing delZ by delR in pkg/ctrl/

1
2 #include "CTRL_CPPOPTIONS.h"
3 #ifdef ALLOW_OBCS
4 # include "OBCS_OPTIONS.h"
5 #endif
6
7
8 subroutine ctrl_getobcss(
9 I mytime,
10 I myiter,
11 I mythid
12 & )
13
14 c ==================================================================
15 c SUBROUTINE ctrl_getobcss
16 c ==================================================================
17 c
18 c o Get southern 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 new flags: gebbie@mit.edu, 25 Jan 2003.
24 c
25 c ==================================================================
26 c SUBROUTINE ctrl_getobcss
27 c ==================================================================
28
29 implicit none
30
31 #ifdef ALLOW_OBCSS_CONTROL
32
33 c == global variables ==
34
35 #include "EEPARAMS.h"
36 #include "SIZE.h"
37 #include "PARAMS.h"
38 #include "GRID.h"
39 #include "OBCS.h"
40
41 #include "ctrl.h"
42 #include "ctrl_dummy.h"
43 #include "optim.h"
44
45 c == routine arguments ==
46
47 _RL mytime
48 integer myiter
49 integer mythid
50
51 c == local variables ==
52
53 integer bi,bj
54 integer i,j,k
55 integer itlo,ithi
56 integer jtlo,jthi
57 integer jmin,jmax
58 integer imin,imax
59 integer ilobcss
60 integer iobcs
61
62 _RL dummy
63 _RL obcssfac
64 logical obcssfirst
65 logical obcsschanged
66 integer obcsscount0
67 integer obcsscount1
68 integer jp1
69
70 cgg _RL maskxz (1-olx:snx+olx,nr,nsx,nsy)
71
72 logical doglobalread
73 logical ladinit
74
75 character*(80) fnameobcss
76
77 cgg( Variables for splitting barotropic/baroclinic vels.
78 _RL vbaro
79 _RL vtop
80 cgg)
81
82 c == external functions ==
83
84 integer ilnblnk
85 external ilnblnk
86
87
88 c == end of interface ==
89
90 jtlo = mybylo(mythid)
91 jthi = mybyhi(mythid)
92 itlo = mybxlo(mythid)
93 ithi = mybxhi(mythid)
94 jmin = 1-oly
95 jmax = sny+oly
96 imin = 1-olx
97 imax = snx+olx
98 jp1 = 1
99
100 cgg( Initialize variables for balancing volume flux.
101 vbaro = 0.d0
102 vtop = 0.d0
103 cgg)
104
105 c-- Now, read the control vector.
106 doglobalread = .false.
107 ladinit = .false.
108
109 if (optimcycle .ge. 0) then
110 ilobcss=ilnblnk( xx_obcss_file )
111 write(fnameobcss(1:80),'(2a,i10.10)')
112 & xx_obcss_file(1:ilobcss), '.', optimcycle
113 endif
114
115 c-- Get the counters, flags, and the interpolation factor.
116 call ctrl_get_gen_rec(
117 I xx_obcssstartdate, xx_obcssperiod,
118 O obcssfac, obcssfirst, obcsschanged,
119 O obcsscount0,obcsscount1,
120 I mytime, myiter, mythid )
121
122 do iobcs = 1,nobcs
123 if ( obcssfirst ) then
124 call active_read_xz_loc( fnameobcss, tmpfldxz,
125 & (obcsscount0-1)*nobcs+iobcs,
126 & doglobalread, ladinit, optimcycle,
127 & mythid, xx_obcss_dummy )
128
129 if ( optimcycle .gt. 0) then
130 if (iobcs .eq. 3) then
131 cgg Special attention is needed for the normal velocity.
132 cgg For the north, this is the v velocity, iobcs = 4.
133 cgg This is done on a columnwise basis here.
134 do bj = jtlo,jthi
135 do bi = itlo, ithi
136 do i = imin,imax
137 j = OB_Js(I,bi,bj)
138
139 cgg The barotropic velocity is stored in the level 1.
140 vbaro = tmpfldxz(i,1,bi,bj)
141 tmpfldxz(i,1,bi,bj) = 0.d0
142 vtop = 0.d0
143
144 do k = 1,Nr
145 cgg If cells are not full, this should be modified with hFac.
146 cgg
147 cgg The xx field (tmpfldxz) does not contain the velocity at the
148 cgg surface level. This velocity is not independent; it must
149 cgg exactly balance the volume flux, since we are dealing with
150 cgg the baroclinic velocity structure..
151 vtop = tmpfldxz(i,k,bi,bj)*
152 & maskS(i,j+jp1,k,bi,bj) * delR(k) + vtop
153 cgg Add the barotropic velocity component.
154 if (maskS(i,j+jp1,k,bi,bj) .ne. 0.) then
155 tmpfldxz(i,k,bi,bj) = tmpfldxz(i,k,bi,bj)+ vbaro
156 endif
157 enddo
158 cgg Compute the baroclinic velocity at level 1. Should balance flux.
159 tmpfldxz(i,1,bi,bj) = tmpfldxz(i,1,bi,bj)
160 & - vtop / delR(1)
161 enddo
162 enddo
163 enddo
164 endif
165
166 if (iobcs .eq. 4) then
167 cgg Special attention is needed for the normal velocity.
168 cgg For the north, this is the v velocity, iobcs = 4.
169 cgg This is done on a columnwise basis here.
170 do bj = jtlo,jthi
171 do bi = itlo, ithi
172 do i = imin,imax
173 j = OB_Js(I,bi,bj)
174
175 cgg The barotropic velocity is stored in the level 1.
176 vbaro = tmpfldxz(i,1,bi,bj)
177 tmpfldxz(i,1,bi,bj) = 0.d0
178 vtop = 0.d0
179
180 do k = 1,Nr
181 cgg If cells are not full, this should be modified with hFac.
182 cgg
183 cgg The xx field (tmpfldxz) does not contain the velocity at the
184 cgg surface level. This velocity is not independent; it must
185 cgg exactly balance the volume flux, since we are dealing with
186 cgg the baroclinic velocity structure..
187 vtop = tmpfldxz(i,k,bi,bj)*
188 & maskW(i,j,k,bi,bj) * delR(k) + vtop
189 cgg Add the barotropic velocity component.
190 if (maskW(i,j,k,bi,bj) .ne. 0.) then
191 tmpfldxz(i,k,bi,bj) = tmpfldxz(i,k,bi,bj)+ vbaro
192 endif
193 enddo
194 cgg Compute the baroclinic velocity at level 1. Should balance flux.
195 tmpfldxz(i,1,bi,bj) = tmpfldxz(i,1,bi,bj)
196 & - vtop / delR(1)
197 enddo
198 enddo
199 enddo
200 endif
201 endif
202
203 do bj = jtlo,jthi
204 do bi = itlo,ithi
205 do k = 1,nr
206 do i = imin,imax
207 xx_obcss1(i,k,bi,bj,iobcs) = tmpfldxz (i,k,bi,bj)
208 cgg & * maskxz (i,k,bi,bj)
209 enddo
210 enddo
211 enddo
212 enddo
213 endif
214
215 if ( (obcssfirst) .or. (obcsschanged)) then
216
217 do bj = jtlo,jthi
218 do bi = itlo,ithi
219 do k = 1,nr
220 do i = imin,imax
221 tmpfldxz(i,k,bi,bj) = xx_obcss1(i,k,bi,bj,iobcs)
222 enddo
223 enddo
224 enddo
225 enddo
226
227 call exf_swapffields_xz( tmpfldxz2, tmpfldxz, mythid)
228
229 do bj = jtlo,jthi
230 do bi = itlo,ithi
231 do k = 1,nr
232 do i = imin,imax
233 xx_obcss0(i,k,bi,bj,iobcs) = tmpfldxz2(i,k,bi,bj)
234 enddo
235 enddo
236 enddo
237 enddo
238
239 call active_read_xz_loc( fnameobcss, tmpfldxz,
240 & (obcsscount1-1)*nobcs+iobcs,
241 & doglobalread, ladinit, optimcycle,
242 & mythid, xx_obcss_dummy )
243
244 if ( optimcycle .gt. 0) then
245 if (iobcs .eq. 3) then
246 cgg Special attention is needed for the normal velocity.
247 cgg For the north, this is the v velocity, iobcs = 4.
248 cgg This is done on a columnwise basis here.
249 do bj = jtlo,jthi
250 do bi = itlo, ithi
251 do i = imin,imax
252 j = OB_Js(I,bi,bj)
253
254 cgg The barotropic velocity is stored in the level 1.
255 vbaro = tmpfldxz(i,1,bi,bj)
256 tmpfldxz(i,1,bi,bj) = 0.d0
257 vtop = 0.d0
258
259 do k = 1,Nr
260 cgg If cells are not full, this should be modified with hFac.
261 cgg
262 cgg The xx field (tmpfldxz) does not contain the velocity at the
263 cgg surface level. This velocity is not independent; it must
264 cgg exactly balance the volume flux, since we are dealing with
265 cgg the baroclinic velocity structure..
266 vtop = tmpfldxz(i,k,bi,bj)*
267 & maskS(i,j+jp1,k,bi,bj) * delR(k) + vtop
268 cgg Add the barotropic velocity component.
269 if (maskS(i,j+jp1,k,bi,bj) .ne. 0.) then
270 tmpfldxz(i,k,bi,bj) = tmpfldxz(i,k,bi,bj)+ vbaro
271 endif
272 enddo
273 cgg Compute the baroclinic velocity at level 1. Should balance flux.
274 tmpfldxz(i,1,bi,bj) = tmpfldxz(i,1,bi,bj)
275 & - vtop / delR(1)
276 enddo
277 enddo
278 enddo
279 endif
280
281 if (iobcs .eq. 4) then
282 cgg Special attention is needed for the normal velocity.
283 cgg For the north, this is the v velocity, iobcs = 4.
284 cgg This is done on a columnwise basis here.
285 do bj = jtlo,jthi
286 do bi = itlo, ithi
287 do i = imin,imax
288 j = OB_Js(I,bi,bj)
289
290 cgg The barotropic velocity is stored in the level 1.
291 vbaro = tmpfldxz(i,1,bi,bj)
292 tmpfldxz(i,1,bi,bj) = 0.d0
293 vtop = 0.d0
294
295 do k = 1,Nr
296 cgg If cells are not full, this should be modified with hFac.
297 cgg
298 cgg The xx field (tmpfldxz) does not contain the velocity at the
299 cgg surface level. This velocity is not independent; it must
300 cgg exactly balance the volume flux, since we are dealing with
301 cgg the baroclinic velocity structure..
302 vtop = tmpfldxz(i,k,bi,bj)*
303 & maskW(i,j,k,bi,bj) * delR(k) + vtop
304 cgg Add the barotropic velocity component.
305 if (maskW(i,j,k,bi,bj) .ne. 0.) then
306 tmpfldxz(i,k,bi,bj) = tmpfldxz(i,k,bi,bj)+ vbaro
307 endif
308 enddo
309 cgg Compute the baroclinic velocity at level 1. Should balance flux.
310 tmpfldxz(i,1,bi,bj) = tmpfldxz(i,1,bi,bj)
311 & - vtop / delR(1)
312 enddo
313 enddo
314 enddo
315 endif
316 endif
317
318 do bj = jtlo,jthi
319 do bi = itlo,ithi
320 do k = 1,nr
321 do i = imin,imax
322 xx_obcss1 (i,k,bi,bj,iobcs) = tmpfldxz (i,k,bi,bj)
323 cgg & * maskxz (i,k,bi,bj)
324 enddo
325 enddo
326 enddo
327 enddo
328 endif
329
330 c-- Add control to model variable.
331 do bj = jtlo,jthi
332 do bi = itlo,ithi
333 c-- Calculate mask for tracer cells (0 => land, 1 => water).
334 do k = 1,nr
335 do i = 1,snx
336 j = OB_Js(I,bi,bj)
337 if (iobcs .EQ. 1) then
338 OBSt(i,k,bi,bj) = OBSt (i,k,bi,bj)
339 & + obcssfac *xx_obcss0(i,k,bi,bj,iobcs)
340 & + (1. _d 0 - obcssfac)*xx_obcss1(i,k,bi,bj,iobcs)
341 OBSt(i,k,bi,bj) = OBSt(i,k,bi,bj)
342 & *maskS(i,j+jp1,k,bi,bj)
343 else if (iobcs .EQ. 2) then
344 OBSs(i,k,bi,bj) = OBSs (i,k,bi,bj)
345 & + obcssfac *xx_obcss0(i,k,bi,bj,iobcs)
346 & + (1. _d 0 - obcssfac)*xx_obcss1(i,k,bi,bj,iobcs)
347 OBSs(i,k,bi,bj) = OBSs(i,k,bi,bj)
348 & *maskS(i,j+jp1,k,bi,bj)
349 else if (iobcs .EQ. 4) then
350 OBSu(i,k,bi,bj) = OBSu (i,k,bi,bj)
351 & + obcssfac *xx_obcss0(i,k,bi,bj,iobcs)
352 & + (1. _d 0 - obcssfac)*xx_obcss1(i,k,bi,bj,iobcs)
353 OBSu(i,k,bi,bj) = OBSu(i,k,bi,bj)
354 & *maskW(i,j,k,bi,bj)
355 else if (iobcs .EQ. 3) then
356 OBSv(i,k,bi,bj) = OBSv (i,k,bi,bj)
357 & + obcssfac *xx_obcss0(i,k,bi,bj,iobcs)
358 & + (1. _d 0 - obcssfac)*xx_obcss1(i,k,bi,bj,iobcs)
359 OBSv(i,k,bi,bj) = OBSv(i,k,bi,bj)
360 & *maskS(i,j+jp1,k,bi,bj)
361 endif
362 enddo
363 enddo
364 enddo
365 enddo
366
367 C-- End over iobcs loop
368 enddo
369
370 #else /* ALLOW_OBCSS_CONTROL undefined */
371
372 c == routine arguments ==
373
374 _RL mytime
375 integer myiter
376 integer mythid
377
378 c-- CPP flag ALLOW_OBCSS_CONTROL undefined.
379
380 #endif /* ALLOW_OBCSS_CONTROL */
381
382 end
383

  ViewVC Help
Powered by ViewVC 1.1.22