1 |
|
C $Header$ |
2 |
|
C $Name$ |
3 |
|
|
4 |
#include "CTRL_CPPOPTIONS.h" |
#include "CTRL_CPPOPTIONS.h" |
5 |
#ifdef ALLOW_OBCS |
#ifdef ALLOW_OBCS |
69 |
integer obcswcount1 |
integer obcswcount1 |
70 |
|
|
71 |
cgg _RL maskyz (1-oly:sny+oly,nr,nsx,nsy) |
cgg _RL maskyz (1-oly:sny+oly,nr,nsx,nsy) |
72 |
|
_RL tmpfldyz (1-oly:sny+oly,nr,nsx,nsy) |
73 |
|
|
74 |
logical doglobalread |
logical doglobalread |
75 |
logical ladinit |
logical ladinit |
110 |
|
|
111 |
if (optimcycle .ge. 0) then |
if (optimcycle .ge. 0) then |
112 |
ilobcsw=ilnblnk( xx_obcsw_file ) |
ilobcsw=ilnblnk( xx_obcsw_file ) |
113 |
write(fnameobcsw(1:80),'(2a,i10.10)') |
write(fnameobcsw(1:80),'(2a,i10.10)') |
114 |
& xx_obcsw_file(1:ilobcsw), '.', optimcycle |
& xx_obcsw_file(1:ilobcsw), '.', optimcycle |
115 |
endif |
endif |
116 |
|
|
121 |
O obcswcount0,obcswcount1, |
O obcswcount0,obcswcount1, |
122 |
I mytime, myiter, mythid ) |
I mytime, myiter, mythid ) |
123 |
|
|
124 |
|
CML print *,'ml-getobcs: ',myIter,obcswfirst,obcswchanged, |
125 |
|
CML & obcswcount0,obcswcount1,obcswfac |
126 |
do iobcs = 1,nobcs |
do iobcs = 1,nobcs |
127 |
if ( obcswfirst ) then |
if ( obcswfirst ) then |
128 |
call active_read_yz_loc( fnameobcsw, tmpfldyz, |
call active_read_yz( fnameobcsw, tmpfldyz, |
129 |
& (obcswcount0-1)*nobcs+iobcs, |
& (obcswcount0-1)*nobcs+iobcs, |
130 |
& doglobalread, ladinit, optimcycle, |
& doglobalread, ladinit, optimcycle, |
131 |
& mythid, xx_obcsw_dummy ) |
& mythid, xx_obcsw_dummy ) |
132 |
|
|
133 |
#ifdef ALLOW_CTRL_OBCS_BALANCE |
#ifdef ALLOW_CTRL_OBCS_BALANCE |
134 |
|
|
135 |
if ( optimcycle .gt. 0) then |
if ( optimcycle .gt. 0) then |
136 |
if (iobcs .eq. 3) then |
if (iobcs .eq. 3) then |
137 |
cgg Special attention is needed for the normal velocity. |
cgg Special attention is needed for the normal velocity. |
138 |
cgg For the north, this is the v velocity, iobcs = 4. |
cgg For the north, this is the v velocity, iobcs = 4. |
149 |
|
|
150 |
do k = 1,Nr |
do k = 1,Nr |
151 |
cgg If cells are not full, this should be modified with hFac. |
cgg If cells are not full, this should be modified with hFac. |
152 |
cgg |
cgg |
153 |
cgg The xx field (tmpfldxz) does not contain the velocity at the |
cgg The xx field (tmpfldyz) does not contain the velocity at the |
154 |
cgg surface level. This velocity is not independent; it must |
cgg surface level. This velocity is not independent; it must |
155 |
cgg exactly balance the volume flux, since we are dealing with |
cgg exactly balance the volume flux, since we are dealing with |
156 |
cgg the baroclinic velocity structure.. |
cgg the baroclinic velocity structure.. |
157 |
utop = tmpfldyz(j,k,bi,bj)* |
utop = tmpfldyz(j,k,bi,bj)* |
158 |
& maskW(i+ip1,j,k,bi,bj) * delR(k) + utop |
& maskW(i+ip1,j,k,bi,bj) * delR(k) + utop |
159 |
cgg Add the barotropic velocity component. |
cgg Add the barotropic velocity component. |
160 |
if (maskW(i+ip1,j,k,bi,bj) .ne. 0.) then |
if (maskW(i+ip1,j,k,bi,bj) .ne. 0.) then |
161 |
tmpfldyz(j,k,bi,bj) = tmpfldyz(j,k,bi,bj)+ ubaro |
tmpfldyz(j,k,bi,bj) = tmpfldyz(j,k,bi,bj)+ ubaro |
162 |
endif |
endif |
163 |
enddo |
enddo |
164 |
cgg Compute the baroclinic velocity at level 1. Should balance flux. |
cgg Compute the baroclinic velocity at level 1. Should balance flux. |
165 |
tmpfldyz(j,1,bi,bj) = tmpfldyz(j,1,bi,bj) |
tmpfldyz(j,1,bi,bj) = tmpfldyz(j,1,bi,bj) |
166 |
& - utop / delR(1) |
& - utop / delR(1) |
167 |
enddo |
enddo |
168 |
enddo |
enddo |
184 |
|
|
185 |
do k = 1,Nr |
do k = 1,Nr |
186 |
cgg If cells are not full, this should be modified with hFac. |
cgg If cells are not full, this should be modified with hFac. |
187 |
cgg |
cgg |
188 |
cgg The xx field (tmpfldxz) does not contain the velocity at the |
cgg The xx field (tmpfldyz) does not contain the velocity at the |
189 |
cgg surface level. This velocity is not independent; it must |
cgg surface level. This velocity is not independent; it must |
190 |
cgg exactly balance the volume flux, since we are dealing with |
cgg exactly balance the volume flux, since we are dealing with |
191 |
cgg the baroclinic velocity structure.. |
cgg the baroclinic velocity structure.. |
192 |
utop = tmpfldyz(j,k,bi,bj)* |
utop = tmpfldyz(j,k,bi,bj)* |
193 |
& maskS(i,j,k,bi,bj) * delR(k) + utop |
& maskS(i,j,k,bi,bj) * delR(k) + utop |
194 |
cgg Add the barotropic velocity component. |
cgg Add the barotropic velocity component. |
195 |
if (maskS(i,j,k,bi,bj) .ne. 0.) then |
if (maskS(i,j,k,bi,bj) .ne. 0.) then |
196 |
tmpfldyz(j,k,bi,bj) = tmpfldyz(j,k,bi,bj)+ ubaro |
tmpfldyz(j,k,bi,bj) = tmpfldyz(j,k,bi,bj)+ ubaro |
197 |
endif |
endif |
198 |
enddo |
enddo |
199 |
cgg Compute the baroclinic velocity at level 1. Should balance flux. |
cgg Compute the baroclinic velocity at level 1. Should balance flux. |
200 |
tmpfldyz(j,1,bi,bj) = tmpfldyz(j,1,bi,bj) |
tmpfldyz(j,1,bi,bj) = tmpfldyz(j,1,bi,bj) |
201 |
& - utop / delR(1) |
& - utop / delR(1) |
202 |
enddo |
enddo |
203 |
enddo |
enddo |
220 |
endif |
endif |
221 |
|
|
222 |
if ( (obcswfirst) .or. (obcswchanged)) then |
if ( (obcswfirst) .or. (obcswchanged)) then |
|
|
|
|
cgg( This is a terribly long way to do it. However, the dimensions do not exactly |
|
|
cgg match up. I will blame Fortran for the ugliness. |
|
223 |
|
|
224 |
do bj = jtlo,jthi |
do bj = jtlo,jthi |
225 |
do bi = itlo,ithi |
do bi = itlo,ithi |
226 |
do k = 1,nr |
do k = 1,nr |
227 |
do j = jmin,jmax |
do j = jmin,jmax |
228 |
tmpfldyz(j,k,bi,bj) = xx_obcsw1(j,k,bi,bj,iobcs) |
xx_obcsw0(j,k,bi,bj,iobcs) = xx_obcsw1(j,k,bi,bj,iobcs) |
229 |
enddo |
tmpfldyz (j,k,bi,bj) = 0. _d 0 |
230 |
enddo |
enddo |
|
enddo |
|
|
enddo |
|
|
|
|
|
call exf_swapffields_yz( tmpfldyz2, tmpfldyz, mythid) |
|
|
|
|
|
do bj = jtlo,jthi |
|
|
do bi = itlo,ithi |
|
|
do k = 1,nr |
|
|
do j = jmin,jmax |
|
|
xx_obcsw0(j,k,bi,bj,iobcs) = tmpfldyz2(j,k,bi,bj) |
|
|
enddo |
|
|
enddo |
|
231 |
enddo |
enddo |
232 |
|
enddo |
233 |
enddo |
enddo |
234 |
|
|
235 |
call active_read_yz_loc( fnameobcsw, tmpfldyz, |
call active_read_yz( fnameobcsw, tmpfldyz, |
236 |
& (obcswcount1-1)*nobcs+iobcs, |
& (obcswcount1-1)*nobcs+iobcs, |
237 |
& doglobalread, ladinit, optimcycle, |
& doglobalread, ladinit, optimcycle, |
238 |
& mythid, xx_obcsw_dummy ) |
& mythid, xx_obcsw_dummy ) |
239 |
|
|
240 |
#ifdef ALLOW_CTRL_OBCS_BALANCE |
#ifdef ALLOW_CTRL_OBCS_BALANCE |
241 |
|
|
242 |
if ( optimcycle .gt. 0) then |
if ( optimcycle .gt. 0) then |
243 |
if (iobcs .eq. 3) then |
if (iobcs .eq. 3) then |
244 |
cgg Special attention is needed for the normal velocity. |
cgg Special attention is needed for the normal velocity. |
245 |
cgg For the north, this is the v velocity, iobcs = 4. |
cgg For the north, this is the v velocity, iobcs = 4. |
256 |
|
|
257 |
do k = 1,Nr |
do k = 1,Nr |
258 |
cgg If cells are not full, this should be modified with hFac. |
cgg If cells are not full, this should be modified with hFac. |
259 |
cgg |
cgg |
260 |
cgg The xx field (tmpfldxz) does not contain the velocity at the |
cgg The xx field (tmpfldyz) does not contain the velocity at the |
261 |
cgg surface level. This velocity is not independent; it must |
cgg surface level. This velocity is not independent; it must |
262 |
cgg exactly balance the volume flux, since we are dealing with |
cgg exactly balance the volume flux, since we are dealing with |
263 |
cgg the baroclinic velocity structure.. |
cgg the baroclinic velocity structure.. |
264 |
utop = tmpfldyz(j,k,bi,bj)* |
utop = tmpfldyz(j,k,bi,bj)* |
265 |
& maskW(i+ip1,j,k,bi,bj) * delR(k) + utop |
& maskW(i+ip1,j,k,bi,bj) * delR(k) + utop |
266 |
cgg Add the barotropic velocity component. |
cgg Add the barotropic velocity component. |
267 |
if (maskW(i+ip1,j,k,bi,bj) .ne. 0.) then |
if (maskW(i+ip1,j,k,bi,bj) .ne. 0.) then |
268 |
tmpfldyz(j,k,bi,bj) = tmpfldyz(j,k,bi,bj)+ ubaro |
tmpfldyz(j,k,bi,bj) = tmpfldyz(j,k,bi,bj)+ ubaro |
269 |
endif |
endif |
270 |
enddo |
enddo |
271 |
cgg Compute the baroclinic velocity at level 1. Should balance flux. |
cgg Compute the baroclinic velocity at level 1. Should balance flux. |
272 |
tmpfldyz(j,1,bi,bj) = tmpfldyz(j,1,bi,bj) |
tmpfldyz(j,1,bi,bj) = tmpfldyz(j,1,bi,bj) |
273 |
& - utop / delR(1) |
& - utop / delR(1) |
274 |
enddo |
enddo |
275 |
enddo |
enddo |
291 |
|
|
292 |
do k = 1,Nr |
do k = 1,Nr |
293 |
cgg If cells are not full, this should be modified with hFac. |
cgg If cells are not full, this should be modified with hFac. |
294 |
cgg |
cgg |
295 |
cgg The xx field (tmpfldxz) does not contain the velocity at the |
cgg The xx field (tmpfldyz) does not contain the velocity at the |
296 |
cgg surface level. This velocity is not independent; it must |
cgg surface level. This velocity is not independent; it must |
297 |
cgg exactly balance the volume flux, since we are dealing with |
cgg exactly balance the volume flux, since we are dealing with |
298 |
cgg the baroclinic velocity structure.. |
cgg the baroclinic velocity structure.. |
299 |
utop = tmpfldyz(j,k,bi,bj)* |
utop = tmpfldyz(j,k,bi,bj)* |
300 |
& maskS(i,j,k,bi,bj) * delR(k) + utop |
& maskS(i,j,k,bi,bj) * delR(k) + utop |
301 |
cgg Add the barotropic velocity component. |
cgg Add the barotropic velocity component. |
302 |
if (maskS(i,j,k,bi,bj) .ne. 0.) then |
if (maskS(i,j,k,bi,bj) .ne. 0.) then |
303 |
tmpfldyz(j,k,bi,bj) = tmpfldyz(j,k,bi,bj)+ ubaro |
tmpfldyz(j,k,bi,bj) = tmpfldyz(j,k,bi,bj)+ ubaro |
304 |
endif |
endif |
305 |
enddo |
enddo |
306 |
cgg Compute the baroclinic velocity at level 1. Should balance flux. |
cgg Compute the baroclinic velocity at level 1. Should balance flux. |
307 |
tmpfldyz(j,1,bi,bj) = tmpfldyz(j,1,bi,bj) |
tmpfldyz(j,1,bi,bj) = tmpfldyz(j,1,bi,bj) |
308 |
& - utop / delR(1) |
& - utop / delR(1) |
309 |
enddo |
enddo |
310 |
enddo |
enddo |
326 |
enddo |
enddo |
327 |
endif |
endif |
328 |
|
|
329 |
c-- Add control to model variable. |
c-- Add control to model variable. |
330 |
do bj = jtlo, jthi |
do bj = jtlo, jthi |
331 |
do bi = itlo, ithi |
do bi = itlo, ithi |
332 |
c-- Calculate mask for tracer cells (0 => land, 1 => water). |
c-- Calculate mask for tracer cells (0 => land, 1 => water). |
333 |
do k = 1,nr |
do k = 1,nr |
334 |
do j = 1,sny |
do j = 1,sny |
335 |
i = OB_Iw(j,bi,bj) |
i = OB_Iw(j,bi,bj) |
336 |
if (iobcs .EQ. 1) then |
if (iobcs .EQ. 1) then |
337 |
OBWt(j,k,bi,bj) = OBWt (j,k,bi,bj) |
OBWt(j,k,bi,bj) = OBWt (j,k,bi,bj) |
338 |
& + obcswfac *xx_obcsw0(j,k,bi,bj,iobcs) |
& + obcswfac *xx_obcsw0(j,k,bi,bj,iobcs) |
339 |
& + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs) |
& + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs) |
340 |
OBWt(j,k,bi,bj) = OBWt(j,k,bi,bj) |
OBWt(j,k,bi,bj) = OBWt(j,k,bi,bj) |
341 |
& *maskW(i+ip1,j,k,bi,bj) |
& *maskW(i+ip1,j,k,bi,bj) |
342 |
else if (iobcs .EQ. 2) then |
else if (iobcs .EQ. 2) then |
343 |
OBWs(j,k,bi,bj) = OBWs (j,k,bi,bj) |
OBWs(j,k,bi,bj) = OBWs (j,k,bi,bj) |
344 |
& + obcswfac *xx_obcsw0(j,k,bi,bj,iobcs) |
& + obcswfac *xx_obcsw0(j,k,bi,bj,iobcs) |
345 |
& + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs) |
& + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs) |
346 |
OBWs(j,k,bi,bj) = OBWs(j,k,bi,bj) |
OBWs(j,k,bi,bj) = OBWs(j,k,bi,bj) |
347 |
& *maskW(i+ip1,j,k,bi,bj) |
& *maskW(i+ip1,j,k,bi,bj) |
348 |
else if (iobcs .EQ. 3) then |
else if (iobcs .EQ. 3) then |
349 |
OBWu(j,k,bi,bj) = OBWu (j,k,bi,bj) |
OBWu(j,k,bi,bj) = OBWu (j,k,bi,bj) |
350 |
& + obcswfac *xx_obcsw0(j,k,bi,bj,iobcs) |
& + obcswfac *xx_obcsw0(j,k,bi,bj,iobcs) |
351 |
& + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs) |
& + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs) |
352 |
OBWu(j,k,bi,bj) = OBWu(j,k,bi,bj) |
OBWu(j,k,bi,bj) = OBWu(j,k,bi,bj) |
353 |
& *maskW(i+ip1,j,k,bi,bj) |
& *maskW(i+ip1,j,k,bi,bj) |
354 |
else if (iobcs .EQ. 4) then |
else if (iobcs .EQ. 4) then |
355 |
OBWv(j,k,bi,bj) = OBWv (j,k,bi,bj) |
OBWv(j,k,bi,bj) = OBWv (j,k,bi,bj) |
356 |
& + obcswfac *xx_obcsw0(j,k,bi,bj,iobcs) |
& + obcswfac *xx_obcsw0(j,k,bi,bj,iobcs) |
357 |
& + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs) |
& + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs) |
358 |
OBWv(j,k,bi,bj) = OBWv(j,k,bi,bj) |
OBWv(j,k,bi,bj) = OBWv(j,k,bi,bj) |
359 |
& *maskS(i,j,k,bi,bj) |
& *maskS(i,j,k,bi,bj) |
360 |
endif |
endif |
|
enddo |
|
|
enddo |
|
361 |
enddo |
enddo |
362 |
|
enddo |
363 |
|
enddo |
364 |
enddo |
enddo |
365 |
|
|
366 |
C-- End over iobcs loop |
C-- End over iobcs loop |