C $Header: /u/gcmpack/MITgcm/pkg/ctrl/ctrl_getobcsn.F,v 1.17 2014/10/09 00:49:26 gforget Exp $
C $Name: $
#include "CTRL_OPTIONS.h"
#ifdef ALLOW_OBCS
# include "OBCS_OPTIONS.h"
#endif
subroutine CTRL_GETOBCSN(
I mytime,
I myiter,
I mythid
& )
c ==================================================================
c SUBROUTINE ctrl_getobcsn
c ==================================================================
c
c o Get northern obc of the control vector and add it
c to dyn. fields
c
c started: heimbach@mit.edu, 29-Aug-2001
c
c modified: gebbie@mit.edu, 18-Mar-2003
c ==================================================================
c SUBROUTINE ctrl_getobcsn
c ==================================================================
implicit none
c == global variables ==
#ifdef ALLOW_OBCSN_CONTROL
#include "EEPARAMS.h"
#include "SIZE.h"
#include "PARAMS.h"
#include "GRID.h"
c#include "OBCS_PARAMS.h"
#include "OBCS_GRID.h"
#include "OBCS_FIELDS.h"
#include "CTRL_SIZE.h"
#include "ctrl.h"
#include "ctrl_dummy.h"
#include "CTRL_OBCS.h"
#include "optim.h"
#endif /* ALLOW_OBCSN_CONTROL */
c == routine arguments ==
_RL mytime
integer myiter
integer mythid
#ifdef ALLOW_OBCSN_CONTROL
c == local variables ==
integer bi,bj
integer i,j,k
integer itlo,ithi
integer jtlo,jthi
integer jmin,jmax
integer imin,imax
integer ilobcsn
integer iobcs
integer jp1
_RL dummy
_RL obcsnfac
logical obcsnfirst
logical obcsnchanged
integer obcsncount0
integer obcsncount1
cgg _RL maskxz (1-olx:snx+olx,nr,nsx,nsy)
_RL tmpfldxz (1-olx:snx+olx,nr,nsx,nsy)
logical doglobalread
logical ladinit
character*(80) fnameobcsn
#ifdef ALLOW_OBCS_CONTROL_MODES
integer nk,nz
_RL tmpz (nr,nsx,nsy)
_RL stmp
#endif
c == external functions ==
integer ilnblnk
external
c == end of interface ==
jtlo = mybylo(mythid)
jthi = mybyhi(mythid)
itlo = mybxlo(mythid)
ithi = mybxhi(mythid)
cgg jmin = 1-oly
cgg jmax = sny+oly
cgg imin = 1-olx
cgg imax = snx+olx
jmin = 1
jmax = sny
imin = 1
imax = snx
jp1 = 0
c-- Now, read the control vector.
doglobalread = .false.
ladinit = .false.
if (optimcycle .ge. 0) then
ilobcsn=ilnblnk( xx_obcsn_file )
write(fnameobcsn(1:80),'(2a,i10.10)')
& xx_obcsn_file(1:ilobcsn), '.', optimcycle
endif
c-- Get the counters, flags, and the interpolation factor.
call CTRL_GET_GEN_REC(
I xx_obcsnstartdate, xx_obcsnperiod,
O obcsnfac, obcsnfirst, obcsnchanged,
O obcsncount0,obcsncount1,
I mytime, myiter, mythid )
do iobcs = 1,nobcs
if ( obcsnfirst ) then
call ACTIVE_READ_XZ( fnameobcsn, tmpfldxz,
& (obcsncount0-1)*nobcs+iobcs,
& doglobalread, ladinit, optimcycle,
& mythid, xx_obcsn_dummy )
do bj = jtlo,jthi
do bi = itlo,ithi
#ifdef ALLOW_OBCS_CONTROL_MODES
if (iobcs .gt. 2) then
do i = imin,imax
j = OB_Jn(i,bi,bj)
IF ( j.EQ.OB_indexNone ) j = 1
cih Determine number of open vertical layers.
nz = 0
do k = 1,Nr
if (iobcs .eq. 3) then
nz = nz + maskS(i,j+jp1,k,bi,bj)
else
nz = nz + maskW(i,j,k,bi,bj)
endif
end
do
cih Compute absolute velocities from the barotropic-baroclinic modes.
do k = 1,Nr
if (k.le.nz) then
stmp = 0.
do nk = 1,nz
stmp = stmp +
& modesv(k,nk,nz)*tmpfldxz(i,nk,bi,bj)
end
do
tmpz(k,bi,bj) = stmp
else
tmpz(k,bi,bj) = 0.
end
if
end
do
do k = 1,Nr
if (iobcs .eq. 3) then
tmpfldxz(i,k,bi,bj) = tmpz(k,bi,bj)
& *recip_hFacS(i,j+jp1,k,bi,bj)
else
tmpfldxz(i,k,bi,bj) = tmpz(k,bi,bj)
& *recip_hFacW(i,j,k,bi,bj)
endif
end
do
enddo
endif
#endif
do k = 1,nr
do i = imin,imax
xx_obcsn1(i,k,bi,bj,iobcs) = tmpfldxz (i,k,bi,bj)
cgg & * maskxz (i,k,bi,bj)
enddo
enddo
enddo
enddo
endif
if ( (obcsnfirst) .or. (obcsnchanged)) then
do bj = jtlo,jthi
do bi = itlo,ithi
do k = 1,nr
do i = imin,imax
xx_obcsn0(i,k,bi,bj,iobcs) = xx_obcsn1(i,k,bi,bj,iobcs)
tmpfldxz (i,k,bi,bj) = 0. _d 0
enddo
enddo
enddo
enddo
call ACTIVE_READ_XZ( fnameobcsn, tmpfldxz,
& (obcsncount1-1)*nobcs+iobcs,
& doglobalread, ladinit, optimcycle,
& mythid, xx_obcsn_dummy )
do bj = jtlo,jthi
do bi = itlo,ithi
#ifdef ALLOW_OBCS_CONTROL_MODES
if (iobcs .gt. 2) then
do i = imin,imax
j = OB_Jn(i,bi,bj)
IF ( j.EQ.OB_indexNone ) j = 1
cih Determine number of open vertical layers.
nz = 0
do k = 1,Nr
if (iobcs .eq. 3) then
nz = nz + maskS(i,j+jp1,k,bi,bj)
else
nz = nz + maskW(i,j,k,bi,bj)
endif
end
do
cih Compute absolute velocities from the barotropic-baroclinic modes.
do k = 1,Nr
if (k.le.nz) then
stmp = 0.
do nk = 1,nz
stmp = stmp +
& modesv(k,nk,nz)*tmpfldxz(i,nk,bi,bj)
end
do
tmpz(k,bi,bj) = stmp
else
tmpz(k,bi,bj) = 0.
end
if
end
do
do k = 1,Nr
if (iobcs .eq. 3) then
tmpfldxz(i,k,bi,bj) = tmpz(k,bi,bj)
& *recip_hFacS(i,j+jp1,k,bi,bj)
else
tmpfldxz(i,k,bi,bj) = tmpz(k,bi,bj)
& *recip_hFacW(i,j,k,bi,bj)
endif
end
do
enddo
endif
#endif
do k = 1,nr
do i = imin,imax
xx_obcsn1 (i,k,bi,bj,iobcs) = tmpfldxz (i,k,bi,bj)
cgg & * maskxz (i,k,bi,bj)
enddo
enddo
enddo
enddo
endif
c-- Add control to model variable.
do bj = jtlo,jthi
do bi = itlo,ithi
c-- Calculate mask for tracer cells (0 => land, 1 => water).
do k = 1,nr
do i = 1,snx
j = OB_Jn(I,bi,bj)
IF ( j.EQ.OB_indexNone ) j = 1
if (iobcs .EQ. 1) then
OBNt(i,k,bi,bj) = OBNt (i,k,bi,bj)
& + obcsnfac *xx_obcsn0(i,k,bi,bj,iobcs)
& + (1. _d 0 - obcsnfac)*xx_obcsn1(i,k,bi,bj,iobcs)
OBNt(i,k,bi,bj) = OBNt(i,k,bi,bj)
& *maskS(i,j+jp1,k,bi,bj)
else if (iobcs .EQ. 2) then
OBNs(i,k,bi,bj) = OBNs (i,k,bi,bj)
& + obcsnfac *xx_obcsn0(i,k,bi,bj,iobcs)
& + (1. _d 0 - obcsnfac)*xx_obcsn1(i,k,bi,bj,iobcs)
OBNs(i,k,bi,bj) = OBNs(i,k,bi,bj)
& *maskS(i,j+jp1,k,bi,bj)
else if (iobcs .EQ. 4) then
OBNu(i,k,bi,bj) = OBNu (i,k,bi,bj)
& + obcsnfac *xx_obcsn0(i,k,bi,bj,iobcs)
& + (1. _d 0 - obcsnfac)*xx_obcsn1(i,k,bi,bj,iobcs)
OBNu(i,k,bi,bj) = OBNu(i,k,bi,bj)
& *maskW(i,j,k,bi,bj)
else if (iobcs .EQ. 3) then
OBNv(i,k,bi,bj) = OBNv (i,k,bi,bj)
& + obcsnfac *xx_obcsn0(i,k,bi,bj,iobcs)
& + (1. _d 0 - obcsnfac)*xx_obcsn1(i,k,bi,bj,iobcs)
OBNv(i,k,bi,bj) = OBNv(i,k,bi,bj)
& *maskS(i,j+jp1,k,bi,bj)
endif
enddo
enddo
enddo
enddo
C-- End over iobcs loop
enddo
#endif /* ALLOW_OBCSN_CONTROL */
return
end