C $Header: /u/gcmpack/MITgcm/pkg/ctrl/ctrl_map_ini_ecco.F,v 1.3 2005/04/07 23:38:43 heimbach Exp $
#include "CTRL_CPPOPTIONS.h"
subroutine CTRL_MAP_INI_ECCO( mythid )
c ==================================================================
c SUBROUTINE ctrl_map_ini_ecco
c ==================================================================
c
c o Add the temperature and salinity parts of the control vector to
c the model state and update the tile edges. The control vector is
c defined in the header file "ctrl.h".
c
c started: Christian Eckert eckert@mit.edu 30-Jun-1999
c
c changed: Christian Eckert eckert@mit.edu 23-Feb-2000
c
c - Restructured the code in order to create a package
c for the MITgcmUV.
c
c ==================================================================
c SUBROUTINE ctrl_map_ini_ecco
c ==================================================================
implicit none
c == global variables ==
#include "EEPARAMS.h"
#include "SIZE.h"
#include "PARAMS.h"
#include "GRID.h"
#include "DYNVARS.h"
#include "FFIELDS.h"
#include "ctrl.h"
#include "ctrl_dummy.h"
#include "optim.h"
#ifdef ALLOW_ECCO
#include "ecco_cost.h"
#endif
#ifdef ALLOW_AUTODIFF_TAMC
#include "tamc.h"
#include "tamc_keys.h"
#endif /* ALLOW_AUTODIFF_TAMC */
c == routine arguments ==
integer mythid
#ifdef ALLOW_ECCO
c == local variables ==
_RL fac
integer bi,bj
integer i,j,k
integer itlo,ithi
integer jtlo,jthi
integer jmin,jmax
integer imin,imax
integer il
logical doglobalread
logical ladinit
character*( 80) fnamegeneric
c == external ==
integer ilnblnk
external
c == end of interface ==
#ifdef ALLOW_AUTODIFF_TAMC
act3 = myThid - 1
max3 = nTx*nTy
act4 = 0
ikey = (act3 + 1) + act4*max3
#endif /* ALLOW_AUTODIFF_TAMC */
jtlo = mybylo(mythid)
jthi = mybyhi(mythid)
itlo = mybxlo(mythid)
ithi = mybxhi(mythid)
jmin = 1-oly
jmax = sny+oly
imin = 1-olx
imax = snx+olx
doglobalread = .false.
ladinit = .false.
fac = 1. _d 0
#ifdef ALLOW_THETA0_CONTROL
c-- Temperature field.
il=ilnblnk( xx_theta_file )
write(fnamegeneric(1:80),'(2a,i10.10)')
& xx_theta_file(1:il),'.',optimcycle
call ACTIVE_READ_XYZ_LOC( fnamegeneric, tmpfld3d, 1,
& doglobalread, ladinit, optimcycle,
& mythid, xx_theta_dummy )
cph(
print *, 'ph-ctrl theta0 a ',
& tmpfld3d(15,15,1,1,1)
cph)
do bj = jtlo,jthi
do bi = itlo,ithi
do k = 1,nr
do j = jmin,jmax
do i = imin,imax
IF (abs(tmpfld3d(i,j,k,bi,bj)).gt.
& 2.0/sqrt(wtheta(k,bi,bj)))
& tmpfld3d(i,j,k,bi,bj)=
& sign(2.0/sqrt(wtheta(k,bi,bj)),
& tmpfld3d(i,j,k,bi,bj))
enddo
enddo
enddo
enddo
enddo
#ifdef ALLOW_AUTODIFF_TAMC
CADJ STORE theta = tapelev_init, key = ikey
#endif
c if(abs(xc(i,j,bi,bj)-123) .gt. 2 .or.
c $ abs(yc(i,j,bi,bj)-84) .gt. 2 )
c if((abs(xc(i,j,bi,bj)-125) .gt. 8 .or.
c $ abs(yc(i,j,bi,bj)-5.5) .gt. 9 )
c $ .and. (k.le.10.or.
c $ (abs(xc(i,j,bi,bj)-164) .gt. 16 .or.
c $ abs(yc(i,j,bi,bj)+60.5) .gt. 6) .and.
c $ (abs(xc(i,j,bi,bj)-359).gt.2 .or.
c $ abs(yc(i,j,bi,bj)-35.5).gt.2) ) .and.
c $ (k.le.17.or.
c $ (abs(xc(i,j,bi,bj)-203) .gt. 26 .or.
c $ abs(yc(i,j,bi,bj)+58.5) .gt. 6)
c $ ))
do bj = jtlo,jthi
do bi = itlo,ithi
do k = 1,nr
do j = jmin,jmax
do i = imin,imax
theta(i,j,k,bi,bj) = theta(i,j,k,bi,bj) +
& fac*tmpfld3d(i,j,k,bi,bj)
if(theta(i,j,k,bi,bj).lt.-2.0)
& theta(i,j,k,bi,bj)= -2.0
enddo
enddo
enddo
enddo
enddo
cph(
print *, 'ph-ctrl theta0 b ',
& tmpfld3d(15,15,1,1,1)
cph)
#endif
#ifdef ALLOW_SALT0_CONTROL
c-- Temperature field.
il=ilnblnk( xx_salt_file )
write(fnamegeneric(1:80),'(2a,i10.10)')
& xx_salt_file(1:il),'.',optimcycle
call ACTIVE_READ_XYZ_LOC( fnamegeneric, tmpfld3d, 1,
& doglobalread, ladinit, optimcycle,
& mythid, xx_salt_dummy )
do bj = jtlo,jthi
do bi = itlo,ithi
do k = 1,nr
do j = jmin,jmax
do i = imin,imax
IF (abs(tmpfld3d(i,j,k,bi,bj)).gt.
& 2.0/sqrt(wsalt(k,bi,bj)))
& tmpfld3d(i,j,k,bi,bj)=
& sign(2.0/sqrt(wsalt(k,bi,bj)),
& tmpfld3d(i,j,k,bi,bj))
c if(abs(xc(i,j,bi,bj)-123) .gt. 2 .or.
c $ abs(yc(i,j,bi,bj)-84) .gt. 2 )
c if((abs(xc(i,j,bi,bj)-125) .gt. 8 .or.
c $ abs(yc(i,j,bi,bj)-5.5) .gt. 9)
c $ .and. (k.le.10.or.
c $ (abs(xc(i,j,bi,bj)-164) .gt. 16 .or.
c $ abs(yc(i,j,bi,bj)+60.5) .gt. 6) )
c $ .and. (k.le.17.or.
c $ (abs(xc(i,j,bi,bj)-203) .gt. 26 .or.
c $ abs(yc(i,j,bi,bj)+58.5) .gt. 6)
c $ ) )
enddo
enddo
enddo
enddo
enddo
#ifdef ALLOW_AUTODIFF_TAMC
CADJ STORE salt = tapelev_init, key = ikey
#endif
do bj = jtlo,jthi
do bi = itlo,ithi
do k = 1,nr
do j = jmin,jmax
do i = imin,imax
salt(i,j,k,bi,bj) = salt(i,j,k,bi,bj) +
& fac*tmpfld3d(i,j,k,bi,bj)
enddo
enddo
enddo
enddo
enddo
#endif
#ifdef ALLOW_EDTAUX_CONTROL
c-- zonal eddy stress : edtaux
il=ilnblnk( xx_edtaux_file )
write(fnamegeneric(1:80),'(2a,i10.10)')
& xx_edtaux_file(1:il),'.',optimcycle
call ACTIVE_READ_XYZ( fnamegeneric, tmpfld3d, 1,
& doglobalread, ladinit, optimcycle,
& mythid, xx_edtaux_dummy )
do bj = jtlo,jthi
do bi = itlo,ithi
do k = 1,nr
do j = jmin,jmax
do i = imin,imax
Eddytaux(i,j,k,bi,bj) = Eddytaux(i,j,k,bi,bj) +
& tmpfld3d(i,j,k,bi,bj)
enddo
enddo
enddo
enddo
enddo
#endif
#ifdef ALLOW_EDTAUY_CONTROL
c-- meridional eddy stress : edtauy
il=ilnblnk( xx_edtauy_file )
write(fnamegeneric(1:80),'(2a,i10.10)')
& xx_edtauy_file(1:il),'.',optimcycle
call ACTIVE_READ_XYZ( fnamegeneric, tmpfld3d, 1,
& doglobalread, ladinit, optimcycle,
& mythid, xx_edtauy_dummy )
do bj = jtlo,jthi
do bi = itlo,ithi
do k = 1,nr
do j = jmin,jmax
do i = imin,imax
Eddytauy(i,j,k,bi,bj) = Eddytauy(i,j,k,bi,bj) +
& tmpfld3d(i,j,k,bi,bj)
enddo
enddo
enddo
enddo
enddo
#endif
#ifdef ALLOW_UVEL0_CONTROL
c-- initial zonal velocity
il=ilnblnk( xx_uvel_file )
write(fnamegeneric(1:80),'(2a,i10.10)')
& xx_uvel_file(1:il),'.',optimcycle
call ACTIVE_READ_XYZ( fnamegeneric, tmpfld3d, 1,
& doglobalread, ladinit, optimcycle,
& mythid, xx_uvel_dummy )
do bj = jtlo,jthi
do bi = itlo,ithi
do k = 1,nr
do j = jmin,jmax
do i = imin,imax
uVel(i,j,k,bi,bj) = uVel(i,j,k,bi,bj) +
& tmpfld3d(i,j,k,bi,bj)
enddo
enddo
enddo
enddo
enddo
#endif
#ifdef ALLOW_VVEL0_CONTROL
c-- initial merid. velocity
il=ilnblnk( xx_vvel_file )
write(fnamegeneric(1:80),'(2a,i10.10)')
& xx_vvel_file(1:il),'.',optimcycle
call ACTIVE_READ_XYZ( fnamegeneric, tmpfld3d, 1,
& doglobalread, ladinit, optimcycle,
& mythid, xx_vvel_dummy )
do bj = jtlo,jthi
do bi = itlo,ithi
do k = 1,nr
do j = jmin,jmax
do i = imin,imax
vVel(i,j,k,bi,bj) = vVel(i,j,k,bi,bj) +
& tmpfld3d(i,j,k,bi,bj)
enddo
enddo
enddo
enddo
enddo
#endif
#ifdef ALLOW_ETAN0_CONTROL
c-- initial Eta.
il=ilnblnk( xx_etan_file )
write(fnamegeneric(1:80),'(2a,i10.10)')
& xx_etan_file(1:il),'.',optimcycle
call ACTIVE_READ_XY_LOC ( fnamegeneric, tmpfld2d, 1,
& doglobalread, ladinit, optimcycle,
& mythid, xx_etan_dummy )
do bj = jtlo,jthi
do bi = itlo,ithi
do j = jmin,jmax
do i = imin,imax
etaN(i,j,bi,bj) = etaN(i,j,bi,bj) + tmpfld2d(i,j,bi,bj)
enddo
enddo
enddo
enddo
#endif
#ifdef ALLOW_RELAXSST_CONTROL
c-- SST relaxation coefficient.
il=ilnblnk( xx_relaxsst_file )
write(fnamegeneric(1:80),'(2a,i10.10)')
& xx_relaxsst_file(1:il),'.',optimcycle
call ACTIVE_READ_XY_LOC ( fnamegeneric, tmpfld2d, 1,
& doglobalread, ladinit, optimcycle,
& mythid, xx_relaxsst_dummy )
do bj = jtlo,jthi
do bi = itlo,ithi
do j = jmin,jmax
do i = imin,imax
lambdaThetaClimRelax(i,j,bi,bj) =
& lambdaThetaClimRelax(i,j,bi,bj)
& + tmpfld2d(i,j,bi,bj)
enddo
enddo
enddo
enddo
#endif
#ifdef ALLOW_RELAXSSS_CONTROL
c-- SSS relaxation coefficient.
il=ilnblnk( xx_relaxsss_file )
write(fnamegeneric(1:80),'(2a,i10.10)')
& xx_relaxsss_file(1:il),'.',optimcycle
call ACTIVE_READ_XY_LOC ( fnamegeneric, tmpfld2d, 1,
& doglobalread, ladinit, optimcycle,
& mythid, xx_relaxsss_dummy )
do bj = jtlo,jthi
do bi = itlo,ithi
do j = jmin,jmax
do i = imin,imax
lambdaSaltClimRelax(i,j,bi,bj) =
& lambdaSaltClimRelax(i,j,bi,bj)
& + tmpfld2d(i,j,bi,bj)
enddo
enddo
enddo
enddo
#endif
c-- Update the tile edges.
#ifdef ALLOW_THETA0_CONTROL
_EXCH_XYZ_R8( theta, mythid )
#endif
#ifdef ALLOW_SALT0_CONTROL
_EXCH_XYZ_R8( salt, mythid )
#endif
#if (defined (ALLOW_EDTAUX_CONTROL) defined (ALLOW_EDTAUY_CONTROL))
CALL EXCH_UV_XYZ_RS(Eddytaux,Eddytauy,.TRUE.,myThid)
#elif (defined (ALLOW_EDTAUX_CONTROL) defined (ALLOW_EDTAUY_CONTROL))
STOP 'ctrl_map_forcing: need BOTH ALLOW_EDTAU[X,Y]_CONTROL'
#endif
#ifdef ALLOW_UVEL0_CONTROL
_EXCH_XYZ_R8( uVel, mythid)
#endif
#ifdef ALLOW_VVEL0_CONTROL
_EXCH_XYZ_R8( vVel, mythid)
#endif
#ifdef ALLOW_ETAN0_CONTROL
_EXCH_XY_R8( etaN, mythid )
#endif
#ifdef ALLOW_RELAXSST_CONTROL
_EXCH_XY_R4( lambdaThetaClimRelax, mythid )
#endif
#ifdef ALLOW_RELAXSSS_CONTROL
_EXCH_XY_R4( lambdaThetaClimRelax, mythid )
#endif
#endif
return
end