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

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

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


Revision 1.16 - (show annotations) (download)
Tue Sep 18 20:21:23 2012 UTC (12 years, 9 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint65b, checkpoint65c, checkpoint65a, checkpoint65d, checkpoint65e, checkpoint64, checkpoint65, checkpoint64y, checkpoint64x, checkpoint64z, checkpoint64q, checkpoint64p, checkpoint64s, checkpoint64r, checkpoint64u, checkpoint64t, checkpoint64w, checkpoint64v, checkpoint64i, checkpoint64h, checkpoint64k, checkpoint64j, checkpoint64m, checkpoint64l, checkpoint64o, checkpoint64n, checkpoint64a, checkpoint64c, checkpoint64b, checkpoint64e, checkpoint64d, checkpoint64g, checkpoint64f
Changes since 1.15: +4 -1 lines
use new parameter OB_indexNone for null index value (instead of hard-coded 0)

1 C $Header: /u/gcmpack/MITgcm/pkg/ctrl/ctrl_getobcsw.F,v 1.15 2012/08/10 19:38:57 jmc Exp $
2 C $Name: $
3
4 #include "CTRL_OPTIONS.h"
5 #ifdef ALLOW_OBCS
6 # include "OBCS_OPTIONS.h"
7 #endif
8
9 subroutine ctrl_getobcsw(
10 I mytime,
11 I myiter,
12 I mythid
13 & )
14
15 c ==================================================================
16 c SUBROUTINE ctrl_getobcsw
17 c ==================================================================
18 c
19 c o Get western obc of the control vector and add it
20 c to dyn. fields
21 c
22 c started: heimbach@mit.edu, 29-Aug-2001
23 c
24 c modified: gebbie@mit.edu, 18-Mar-2003
25 c ==================================================================
26 c SUBROUTINE ctrl_getobcsw
27 c ==================================================================
28
29 implicit none
30
31 c == global variables ==
32 #ifdef ALLOW_OBCSW_CONTROL
33 #include "EEPARAMS.h"
34 #include "SIZE.h"
35 #include "PARAMS.h"
36 #include "GRID.h"
37 c#include "OBCS_PARAMS.h"
38 #include "OBCS_GRID.h"
39 #include "OBCS_FIELDS.h"
40 #include "CTRL_SIZE.h"
41 #include "ctrl.h"
42 #include "ctrl_dummy.h"
43 #include "optim.h"
44 #endif /* ALLOW_OBCSW_CONTROL */
45
46 c == routine arguments ==
47 _RL mytime
48 integer myiter
49 integer mythid
50
51 #ifdef ALLOW_OBCSW_CONTROL
52 c == local variables ==
53
54 integer bi,bj
55 integer i,j,k
56 integer itlo,ithi
57 integer jtlo,jthi
58 integer jmin,jmax
59 integer imin,imax
60 integer ilobcsw
61 integer iobcs
62 integer ip1
63
64 _RL dummy
65 _RL obcswfac
66 logical obcswfirst
67 logical obcswchanged
68 integer obcswcount0
69 integer obcswcount1
70
71 cgg _RL maskyz (1-oly:sny+oly,nr,nsx,nsy)
72 _RL tmpfldyz (1-oly:sny+oly,nr,nsx,nsy)
73
74 logical doglobalread
75 logical ladinit
76
77 character*(80) fnameobcsw
78
79 #ifdef ALLOW_OBCS_CONTROL_MODES
80 integer nk,nz
81 _RL tmpz (nr,nsx,nsy)
82 _RL stmp
83 #endif
84
85 c == external functions ==
86
87 integer ilnblnk
88 external ilnblnk
89
90
91 c == end of interface ==
92
93 jtlo = mybylo(mythid)
94 jthi = mybyhi(mythid)
95 itlo = mybxlo(mythid)
96 ithi = mybxhi(mythid)
97 jmin = 1-oly
98 jmax = sny+oly
99 imin = 1-olx
100 imax = snx+olx
101 ip1 = 1
102
103 c-- Now, read the control vector.
104 doglobalread = .false.
105 ladinit = .false.
106
107 if (optimcycle .ge. 0) then
108 ilobcsw=ilnblnk( xx_obcsw_file )
109 write(fnameobcsw(1:80),'(2a,i10.10)')
110 & xx_obcsw_file(1:ilobcsw), '.', optimcycle
111 endif
112
113 c-- Get the counters, flags, and the interpolation factor.
114 call ctrl_get_gen_rec(
115 I xx_obcswstartdate, xx_obcswperiod,
116 O obcswfac, obcswfirst, obcswchanged,
117 O obcswcount0,obcswcount1,
118 I mytime, myiter, mythid )
119
120 do iobcs = 1,nobcs
121 if ( obcswfirst ) then
122 call active_read_yz( fnameobcsw, tmpfldyz,
123 & (obcswcount0-1)*nobcs+iobcs,
124 & doglobalread, ladinit, optimcycle,
125 & mythid, xx_obcsw_dummy )
126
127 do bj = jtlo,jthi
128 do bi = itlo,ithi
129 #ifdef ALLOW_OBCS_CONTROL_MODES
130 if (iobcs .gt. 2) then
131 do j = jmin,jmax
132 i = OB_Iw(j,bi,bj)
133 IF ( i.EQ.OB_indexNone ) i = 1
134 cih Determine number of open vertical layers.
135 nz = 0
136 do k = 1,Nr
137 if (iobcs .eq. 3) then
138 nz = nz + maskW(i+ip1,j,k,bi,bj)
139 else
140 nz = nz + maskS(i,j,k,bi,bj)
141 endif
142 end do
143 cih Compute absolute velocities from the barotropic-baroclinic modes.
144 do k = 1,Nr
145 if (k.le.nz) then
146 stmp = 0.
147 do nk = 1,nz
148 stmp = stmp +
149 & modesv(k,nk,nz)*tmpfldyz(j,nk,bi,bj)
150 end do
151 tmpz(k,bi,bj) = stmp
152 else
153 tmpz(k,bi,bj) = 0.
154 end if
155 enddo
156 do k = 1,Nr
157 if (iobcs .eq. 3) then
158 tmpfldyz(j,k,bi,bj) = tmpz(k,bi,bj)
159 & *recip_hFacW(i+ip1,j,k,bi,bj)
160 else
161 tmpfldyz(j,k,bi,bj) = tmpz(k,bi,bj)
162 & *recip_hFacS(i,j,k,bi,bj)
163 endif
164 end do
165 enddo
166 endif
167 #endif
168 do k = 1,nr
169 do j = jmin,jmax
170 xx_obcsw1(j,k,bi,bj,iobcs) = tmpfldyz (j,k,bi,bj)
171 cgg & * maskyz (j,k,bi,bj)
172 enddo
173 enddo
174 enddo
175 enddo
176 endif
177
178 if ( (obcswfirst) .or. (obcswchanged)) then
179
180 do bj = jtlo,jthi
181 do bi = itlo,ithi
182 do k = 1,nr
183 do j = jmin,jmax
184 xx_obcsw0(j,k,bi,bj,iobcs) = xx_obcsw1(j,k,bi,bj,iobcs)
185 tmpfldyz (j,k,bi,bj) = 0. _d 0
186 enddo
187 enddo
188 enddo
189 enddo
190
191 call active_read_yz( fnameobcsw, tmpfldyz,
192 & (obcswcount1-1)*nobcs+iobcs,
193 & doglobalread, ladinit, optimcycle,
194 & mythid, xx_obcsw_dummy )
195
196 do bj = jtlo,jthi
197 do bi = itlo,ithi
198 #ifdef ALLOW_OBCS_CONTROL_MODES
199 if (iobcs .gt. 2) then
200 do j = jmin,jmax
201 i = OB_Iw(j,bi,bj)
202 IF ( i.EQ.OB_indexNone ) i = 1
203 cih Determine number of open vertical layers.
204 nz = 0
205 do k = 1,Nr
206 if (iobcs .eq. 3) then
207 nz = nz + maskW(i+ip1,j,k,bi,bj)
208 else
209 nz = nz + maskS(i,j,k,bi,bj)
210 endif
211 end do
212 cih Compute absolute velocities from the barotropic-baroclinic modes.
213 do k = 1,Nr
214 if (k.le.nz) then
215 stmp = 0.
216 do nk = 1,nz
217 stmp = stmp +
218 & modesv(k,nk,nz)*tmpfldyz(j,nk,bi,bj)
219 end do
220 tmpz(k,bi,bj) = stmp
221 else
222 tmpz(k,bi,bj) = 0.
223 end if
224 enddo
225 do k = 1,Nr
226 if (iobcs .eq. 3) then
227 tmpfldyz(j,k,bi,bj) = tmpz(k,bi,bj)
228 & *recip_hFacW(i+ip1,j,k,bi,bj)
229 else
230 tmpfldyz(j,k,bi,bj) = tmpz(k,bi,bj)
231 & *recip_hFacS(i,j,k,bi,bj)
232 endif
233 end do
234 enddo
235 endif
236 #endif
237 do k = 1,nr
238 do j = jmin,jmax
239 xx_obcsw1 (j,k,bi,bj,iobcs) = tmpfldyz (j,k,bi,bj)
240 cgg & * maskyz (j,k,bi,bj)
241 enddo
242 enddo
243 enddo
244 enddo
245 endif
246
247 c-- Add control to model variable.
248 do bj = jtlo, jthi
249 do bi = itlo, ithi
250 c-- Calculate mask for tracer cells (0 => land, 1 => water).
251 do k = 1,nr
252 do j = 1,sny
253 i = OB_Iw(j,bi,bj)
254 IF ( i.EQ.OB_indexNone ) i = 1
255 if (iobcs .EQ. 1) then
256 OBWt(j,k,bi,bj) = OBWt (j,k,bi,bj)
257 & + obcswfac *xx_obcsw0(j,k,bi,bj,iobcs)
258 & + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)
259 OBWt(j,k,bi,bj) = OBWt(j,k,bi,bj)
260 & *maskW(i+ip1,j,k,bi,bj)
261 else if (iobcs .EQ. 2) then
262 OBWs(j,k,bi,bj) = OBWs (j,k,bi,bj)
263 & + obcswfac *xx_obcsw0(j,k,bi,bj,iobcs)
264 & + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)
265 OBWs(j,k,bi,bj) = OBWs(j,k,bi,bj)
266 & *maskW(i+ip1,j,k,bi,bj)
267 else if (iobcs .EQ. 3) then
268 OBWu(j,k,bi,bj) = OBWu (j,k,bi,bj)
269 & + obcswfac *xx_obcsw0(j,k,bi,bj,iobcs)
270 & + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)
271 OBWu(j,k,bi,bj) = OBWu(j,k,bi,bj)
272 & *maskW(i+ip1,j,k,bi,bj)
273 else if (iobcs .EQ. 4) then
274 OBWv(j,k,bi,bj) = OBWv (j,k,bi,bj)
275 & + obcswfac *xx_obcsw0(j,k,bi,bj,iobcs)
276 & + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)
277 OBWv(j,k,bi,bj) = OBWv(j,k,bi,bj)
278 & *maskS(i,j,k,bi,bj)
279 endif
280 enddo
281 enddo
282 enddo
283 enddo
284
285 C-- End over iobcs loop
286 enddo
287
288 #endif /* ALLOW_OBCSW_CONTROL */
289
290 return
291 end

  ViewVC Help
Powered by ViewVC 1.1.22