pkg/ctrl/ctrl_map_ini.F

C $Header: /u/gcmpack/MITgcm/pkg/ctrl/ctrl_map_ini.F,v 1.26 2007/10/30 20:19:13 heimbach Exp $
C $Name:  $

#include "CTRL_CPPOPTIONS.h"

CBOP
C     !ROUTINE: ctrl_map_ini
C     !INTERFACE:
      subroutine ctrl_map_ini( mythid )

C     !DESCRIPTION: \bv
c     *=================================================================
c     | SUBROUTINE ctrl_map_ini
c     | Add the temperature, salinity, and diffusivity parts of the
c     | control vector to the model state and update the tile halos.
c     | The control vector is defined in the header file "ctrl.h".
c     *=================================================================
C     \ev

C     !USES:
      implicit none

c     == global variables ==
#include "SIZE.h"
#include "EEPARAMS.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_PTRACERS
# include "PTRACERS_SIZE.h"
c#include "PTRACERS_PARAMS.h"
# include "PTRACERS_FIELDS.h"
#endif
#ifdef ALLOW_ECCO
# include "ecco_cost.h"
#endif

C     !INPUT/OUTPUT PARAMETERS:
c     == routine arguments ==
      integer mythid

C     !LOCAL VARIABLES:
c     == local variables ==

      integer bi,bj
      integer i,j,k
      integer itlo,ithi
      integer jtlo,jthi
      integer jmin,jmax
      integer imin,imax
      integer il

      logical equal
      logical doglobalread
      logical ladinit

      character*( 80)   fnamegeneric

      _RL     fac
      _RL tmptest

c     == external ==
      integer  ilnblnk
      external ilnblnk

c     == end of interface ==
CEOP

      jtlo = mybylo(mythid)
      jthi = mybyhi(mythid)
      itlo = mybxlo(mythid)
      ithi = mybxhi(mythid)
      jmin = 1
      jmax = sny
      imin = 1
      imax = snx

      doglobalread = .false.
      ladinit      = .false.

      equal = .true.

      if ( equal ) then
        fac = 1. _d 0
      else
        fac = 0. _d 0
      endif

#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( fnamegeneric, tmpfld3d, 1,
     &                      doglobalread, ladinit, optimcycle,
     &                      mythid, xx_theta_dummy )

      do bj = jtlo,jthi
        do bi = itlo,ithi
          do k = 1,nr
            do j = jmin,jmax
              do i = imin,imax
#ifdef ALLOW_ECCO
               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))
#endif
#ifdef ALLOW_AUTODIFF_OPENAD
                theta(i,j,k,bi,bj) = theta(i,j,k,bi,bj) +
     &                               fac*xx_theta(i,j,k,bi,bj) +
     &                               fac*tmpfld3d(i,j,k,bi,bj)
#else
                theta(i,j,k,bi,bj) = theta(i,j,k,bi,bj) +
     &                               fac*tmpfld3d(i,j,k,bi,bj)
#endif
#ifndef DISABLE_CTRL_THETA_LIMIT
                if(theta(i,j,k,bi,bj).lt.-2.0)
     &               theta(i,j,k,bi,bj)= -2.0
#endif
              enddo
            enddo
          enddo
       enddo
      enddo

#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( 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
#ifdef ALLOW_ECCO
               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))
#endif
#ifdef ALLOW_AUTODIFF_OPENAD
                salt(i,j,k,bi,bj) = salt(i,j,k,bi,bj) +
     &                               fac*xx_salt(i,j,k,bi,bj) +
     &                               fac*tmpfld3d(i,j,k,bi,bj)
#else
                salt(i,j,k,bi,bj) = salt(i,j,k,bi,bj) +
     &                               fac*tmpfld3d(i,j,k,bi,bj)
#endif

              enddo
            enddo
          enddo
       enddo
      enddO
#endif

#ifdef ALLOW_TR10_CONTROL
#ifdef ALLOW_PTRACERS
c--   Temperature field.
      il=ilnblnk( xx_tr1_file )
      write(fnamegeneric(1:80),'(2a,i10.10)')
     &     xx_tr1_file(1:il),'.',optimcycle
      call active_read_xyz( fnamegeneric, tmpfld3d, 1,
     &                      doglobalread, ladinit, optimcycle,
     &                      mythid, xx_tr1_dummy )

      do bj = jtlo,jthi
        do bi = itlo,ithi
          do k = 1,nr
            do j = jmin,jmax
              do i = imin,imax
                ptracer(i,j,k,bi,bj,1) = ptracer(i,j,k,bi,bj,1) +
     &                               fac*tmpfld3d(i,j,k,bi,bj)
              enddo
            enddo
          enddo
       enddo
      enddo
#endif
#endif

#ifdef ALLOW_SST0_CONTROL
c--   sst0.
      il=ilnblnk( xx_sst_file )
      write(fnamegeneric(1:80),'(2a,i10.10)')
     &     xx_sst_file(1:il),'.',optimcycle
      call active_read_xy ( fnamegeneric, tmpfld2d, 1,
     &                      doglobalread, ladinit, optimcycle,
     &                      mythid, xx_sst_dummy )
      do bj = jtlo,jthi
        do bi = itlo,ithi
          do j = jmin,jmax
            do i = imin,imax
cph              sst(i,j,bi,bj) = sst(i,j,bi,bj) + tmpfld2d(i,j,bi,bj)
              theta(i,j,1,bi,bj) = theta(i,j,1,bi,bj)
     &                             + tmpfld2d(i,j,bi,bj)
            enddo
          enddo
        enddo
      enddo
#endif

#ifdef ALLOW_SSS0_CONTROL
c--   sss0.
      il=ilnblnk( xx_sss_file )
      write(fnamegeneric(1:80),'(2a,i10.10)')
     &     xx_sss_file(1:il),'.',optimcycle
      call active_read_xy ( fnamegeneric, tmpfld2d, 1,
     &                      doglobalread, ladinit, optimcycle,
     &                      mythid, xx_sss_dummy )
      do bj = jtlo,jthi
        do bi = itlo,ithi
          do j = jmin,jmax
            do i = imin,imax
cph              sss(i,j,bi,bj) = sss(i,j,bi,bj) + tmpfld2d(i,j,bi,bj)
              salt(i,j,1,bi,bj) = salt(i,j,1,bi,bj)
     &                             + tmpfld2d(i,j,bi,bj)
            enddo
          enddo
        enddo
      enddo
#endif

#ifdef ALLOW_DIFFKR_CONTROL
c--   diffkr.
      il=ilnblnk( xx_diffkr_file )
      write(fnamegeneric(1:80),'(2a,i10.10)')
     &     xx_diffkr_file(1:il),'.',optimcycle
      call active_read_xyz( fnamegeneric, tmpfld3d, 1,
     &                      doglobalread, ladinit, optimcycle,
     &                      mythid, xx_diffkr_dummy )
      do bj = jtlo,jthi
        do bi = itlo,ithi
          do k = 1,nr
            do j = jmin,jmax
              do i = imin,imax
                diffkr(i,j,k,bi,bj) = diffkr(i,j,k,bi,bj) +
     &                                tmpfld3d(i,j,k,bi,bj)
              enddo
            enddo
          enddo
       enddo
      enddo
#endif

#ifdef ALLOW_KAPGM_CONTROL
c--   kapgm.
      il=ilnblnk( xx_kapgm_file )
      write(fnamegeneric(1:80),'(2a,i10.10)')
     &     xx_kapgm_file(1:il),'.',optimcycle
      call active_read_xyz( fnamegeneric, tmpfld3d, 1,
     &                      doglobalread, ladinit, optimcycle,
     &                      mythid, xx_kapgm_dummy )
      do bj = jtlo,jthi
        do bi = itlo,ithi
          do k = 1,nr
            do j = jmin,jmax
              do i = imin,imax
                kapgm(i,j,k,bi,bj) = kapgm(i,j,k,bi,bj) +
     &                               tmpfld3d(i,j,k,bi,bj)
              enddo
            enddo
          enddo
       enddo
      enddo
#endif

#ifdef ALLOW_EFLUXY0_CONTROL
c--   y-component EP-flux field.
      il=ilnblnk( xx_efluxy_file )
      write(fnamegeneric(1:80),'(2a,i10.10)')
     &     xx_efluxy_file(1:il),'.',optimcycle
      call active_read_xyz( fnamegeneric, tmpfld3d, 1,
     &                      doglobalread, ladinit, optimcycle,
     &                      mythid, xx_efluxy_dummy )

      do bj = jtlo,jthi
        do bi = itlo,ithi
          do k = 1,nr
            do j = jmin,jmax
              do i = imin,imax
                EfluxY(i,j,k,bi,bj) = EfluxY(i,j,k,bi,bj)
     &                                - fac*tmpfld3d(i,j,k,bi,bj)
     &                                  *maskS(i,j,k,bi,bj)
cph                EfluxY(i,j,k,bi,bj) = EfluxY(i,j,k,bi,bj)
cph     &                                - rSphere*cosFacU(J,bi,bj)
cph     &                                  *fac*tmpfld3d(i,j,k,bi,bj)
              enddo
            enddo
          enddo
       enddo
      enddo
#endif

#ifdef ALLOW_EFLUXP0_CONTROL
c--   p-component EP-flux field.
      il=ilnblnk( xx_efluxp_file )
      write(fnamegeneric(1:80),'(2a,i10.10)')
     &     xx_efluxp_file(1:il),'.',optimcycle
      call active_read_xyz( fnamegeneric, tmpfld3d, 1,
     &                      doglobalread, ladinit, optimcycle,
     &                      mythid, xx_efluxp_dummy )

      do bj = jtlo,jthi
        do bi = itlo,ithi
          do k = 1,nr
            do j = jmin,jmax
              do i = imin,imax
                EfluxP(i,j,k,bi,bj) = EfluxP(i,j,k,bi,bj)
     &                                + fCori(i,j,bi,bj)
     &                                  *fac*tmpfld3d(i,j,k,bi,bj)
     &                                  *hFacV(i,j,k,bi,bj)
cph                EfluxP(i,j,k,bi,bj) = EfluxP(i,j,k,bi,bj)
cph     &                                + fCori(i,j,bi,bj)
cph     &                                  *rSphere*cosFacU(J,bi,bj)
cph     &                                  *fac*tmpfld3d(i,j,k,bi,bj)
              enddo
            enddo
          enddo
       enddo
      enddo
#endif

#ifdef ALLOW_BOTTOMDRAG_CONTROL
c--   bottom drag
      il=ilnblnk( xx_bottomdrag_file )
      write(fnamegeneric(1:80),'(2a,i10.10)')
     &     xx_bottomdrag_file(1:il),'.',optimcycle
      call active_read_xy ( fnamegeneric, tmpfld2d, 1,
     &                      doglobalread, ladinit, optimcycle,
     &                      mythid, xx_bottomdrag_dummy )
      do bj = jtlo,jthi
        do bi = itlo,ithi
          do j = jmin,jmax
            do i = imin,imax
              bottomdragfld(i,j,bi,bj) = bottomdragfld(i,j,bi,bj)
     &                                   + tmpfld2d(i,j,bi,bj)
            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) +
     &            fCori(i,j,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) +
     &            fCoriG(i,j,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
#ifdef ALLOW_AUTODIFF_OPENAD
                uVel(i,j,k,bi,bj) = uVel(i,j,k,bi,bj) +
     &                                  fac*xx_uvel(i,j,k,bi,bj)
#else
                uVel(i,j,k,bi,bj) = uVel(i,j,k,bi,bj) +
     &                                  fac*tmpfld3d(i,j,k,bi,bj)
#endif
              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
#ifdef ALLOW_AUTODIFF_OPENAD
                vVel(i,j,k,bi,bj) = vVel(i,j,k,bi,bj) +
     &                                  fac*xx_vvel(i,j,k,bi,bj)
#else
                vVel(i,j,k,bi,bj) = vVel(i,j,k,bi,bj) +
     &                                  fac*tmpfld3d(i,j,k,bi,bj)
#endif
              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 ( 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
#ifdef ALLOW_AUTODIFF_OPENAD
              etaN(i,j,bi,bj) = etaN(i,j,bi,bj) +
     &                              fac*xx_etan(i,j,bi,bj)
#else
              etaN(i,j,bi,bj) = etaN(i,j,bi,bj) +
     &                              fac*tmpfld2d(i,j,bi,bj)
#endif
            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 ( 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 ( 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

#ifdef ALLOW_SEAICE
      call seaice_ctrl_map_ini( mythid )
#endif

c--   Update the tile edges.

#if (defined (ALLOW_THETA0_CONTROL) || defined (ALLOW_SST0_CONTROL))
      _EXCH_XYZ_R8( theta, mythid )
#endif
#if (defined (ALLOW_SALT0_CONTROL) || defined (ALLOW_SSS0_CONTROL))
      _EXCH_XYZ_R8(  salt, mythid )
#endif
#ifdef ALLOW_TR10_CONTROL
#ifdef ALLOW_PTRACERS
      _EXCH_XYZ_R8(pTracer(1-Olx,1-Oly,1,1,1,1),myThid)
#endif
#endif
#ifdef ALLOW_DIFFKR_CONTROL
      _EXCH_XYZ_R8( diffkr, mythid)
#endif
#ifdef ALLOW_KAPGM_CONTROL
      _EXCH_XYZ_R8( kapgm, mythid)
#endif
#ifdef ALLOW_EFLUXY0_CONTROL
      _EXCH_XYZ_R8( EfluxY, mythid )
#endif
#ifdef ALLOW_EFLUXP0_CONTROL
      _EXCH_XYZ_R8( EfluxP, mythid )
#endif
#ifdef ALLOW_BOTTOMDRAG_CONTROL
      _EXCH_XY_R8( bottomdragfld, 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

      return
      end

1	jmc	1.27	C $Header: /u/gcmpack/MITgcm/pkg/ctrl/ctrl_map_ini.F,v 1.26 2007/10/30 20:19:13 heimbach Exp $
2	jmc	1.25	C $Name: $
3	heimbach	1.1
4			#include "CTRL_CPPOPTIONS.h"
5
6	heimbach	1.5	CBOP
7			C !ROUTINE: ctrl_map_ini
8			C !INTERFACE:
9			subroutine ctrl_map_ini( mythid )
10
11			C !DESCRIPTION: \bv
12			c *=================================================================
13			c \| SUBROUTINE ctrl_map_ini
14	jmc	1.25	c \| Add the temperature, salinity, and diffusivity parts of the
15			c \| control vector to the model state and update the tile halos.
16	heimbach	1.5	c \| The control vector is defined in the header file "ctrl.h".
17			c *=================================================================
18			C \ev
19	heimbach	1.1
20	heimbach	1.5	C !USES:
21	heimbach	1.1	implicit none
22
23			c == global variables ==
24	heimbach	1.6	#include "SIZE.h"
25	heimbach	1.1	#include "EEPARAMS.h"
26	heimbach	1.6	#include "PARAMS.h"
27	heimbach	1.15	#include "GRID.h"
28	heimbach	1.1	#include "DYNVARS.h"
29	heimbach	1.15	#include "FFIELDS.h"
30	heimbach	1.1	#include "ctrl.h"
31			#include "ctrl_dummy.h"
32	heimbach	1.2	#include "optim.h"
33	heimbach	1.13	#ifdef ALLOW_PTRACERS
34			# include "PTRACERS_SIZE.h"
35	jmc	1.27	c#include "PTRACERS_PARAMS.h"
36			# include "PTRACERS_FIELDS.h"
37	heimbach	1.13	#endif
38	heimbach	1.12	#ifdef ALLOW_ECCO
39			# include "ecco_cost.h"
40			#endif
41	heimbach	1.1
42	heimbach	1.5	C !INPUT/OUTPUT PARAMETERS:
43	heimbach	1.1	c == routine arguments ==
44			integer mythid
45
46	heimbach	1.5	C !LOCAL VARIABLES:
47	heimbach	1.1	c == local variables ==
48
49			integer bi,bj
50			integer i,j,k
51			integer itlo,ithi
52			integer jtlo,jthi
53			integer jmin,jmax
54			integer imin,imax
55			integer il
56
57			logical equal
58			logical doglobalread
59			logical ladinit
60
61	heimbach	1.17	character*( 80) fnamegeneric
62	heimbach	1.1
63	heimbach	1.5	_RL fac
64	heimbach	1.12	_RL tmptest
65	heimbach	1.5
66	heimbach	1.1	c == external ==
67			integer ilnblnk
68			external ilnblnk
69
70			c == end of interface ==
71	heimbach	1.5	CEOP
72	heimbach	1.1
73			jtlo = mybylo(mythid)
74			jthi = mybyhi(mythid)
75			itlo = mybxlo(mythid)
76			ithi = mybxhi(mythid)
77	heimbach	1.8	jmin = 1
78			jmax = sny
79			imin = 1
80			imax = snx
81	heimbach	1.1
82			doglobalread = .false.
83			ladinit = .false.
84
85			equal = .true.
86
87			if ( equal ) then
88			fac = 1. _d 0
89			else
90			fac = 0. _d 0
91			endif
92
93			#ifdef ALLOW_THETA0_CONTROL
94			c-- Temperature field.
95			il=ilnblnk( xx_theta_file )
96	heimbach	1.17	write(fnamegeneric(1:80),'(2a,i10.10)')
97	heimbach	1.1	& xx_theta_file(1:il),'.',optimcycle
98	heimbach	1.22	call active_read_xyz( fnamegeneric, tmpfld3d, 1,
99	heimbach	1.1	& doglobalread, ladinit, optimcycle,
100			& mythid, xx_theta_dummy )
101
102			do bj = jtlo,jthi
103			do bi = itlo,ithi
104			do k = 1,nr
105			do j = jmin,jmax
106			do i = imin,imax
107	heimbach	1.12	#ifdef ALLOW_ECCO
108			IF (abs(tmpfld3d(i,j,k,bi,bj)).gt.
109			$ 2.0/sqrt(wtheta(k,bi,bj)))
110			$ tmpfld3d(i,j,k,bi,bj)=
111			$ sign(2.0/sqrt(wtheta(k,bi,bj)),tmpfld3d(i,j,k,bi,bj))
112			#endif
113	utke	1.23	#ifdef ALLOW_AUTODIFF_OPENAD
114	heimbach	1.18	theta(i,j,k,bi,bj) = theta(i,j,k,bi,bj) +
115	utke	1.21	& fac*xx_theta(i,j,k,bi,bj) +
116			& fac*tmpfld3d(i,j,k,bi,bj)
117	heimbach	1.18	#else
118	heimbach	1.1	theta(i,j,k,bi,bj) = theta(i,j,k,bi,bj) +
119			& fac*tmpfld3d(i,j,k,bi,bj)
120	heimbach	1.18	#endif
121	heimbach	1.26	#ifndef DISABLE_CTRL_THETA_LIMIT
122	jmc	1.25	if(theta(i,j,k,bi,bj).lt.-2.0)
123			& theta(i,j,k,bi,bj)= -2.0
124	heimbach	1.26	#endif
125	heimbach	1.1	enddo
126			enddo
127			enddo
128			enddo
129			enddo
130	heimbach	1.12
131	heimbach	1.1	#endif
132
133			#ifdef ALLOW_SALT0_CONTROL
134			c-- Temperature field.
135			il=ilnblnk( xx_salt_file )
136	heimbach	1.17	write(fnamegeneric(1:80),'(2a,i10.10)')
137	heimbach	1.1	& xx_salt_file(1:il),'.',optimcycle
138	heimbach	1.22	call active_read_xyz( fnamegeneric, tmpfld3d, 1,
139	heimbach	1.1	& doglobalread, ladinit, optimcycle,
140			& mythid, xx_salt_dummy )
141
142			do bj = jtlo,jthi
143			do bi = itlo,ithi
144			do k = 1,nr
145			do j = jmin,jmax
146			do i = imin,imax
147	heimbach	1.12	#ifdef ALLOW_ECCO
148			IF (abs(tmpfld3d(i,j,k,bi,bj)).gt.
149			$ 2.0/sqrt(wsalt(k,bi,bj)))
150			$ tmpfld3d(i,j,k,bi,bj)=
151			$ sign(2.0/sqrt(wsalt(k,bi,bj)),tmpfld3d(i,j,k,bi,bj))
152			#endif
153	utke	1.23	#ifdef ALLOW_AUTODIFF_OPENAD
154	heimbach	1.18	salt(i,j,k,bi,bj) = salt(i,j,k,bi,bj) +
155	utke	1.21	& fac*xx_salt(i,j,k,bi,bj) +
156			& fac*tmpfld3d(i,j,k,bi,bj)
157	heimbach	1.18	#else
158	heimbach	1.1	salt(i,j,k,bi,bj) = salt(i,j,k,bi,bj) +
159			& fac*tmpfld3d(i,j,k,bi,bj)
160	heimbach	1.18	#endif
161	heimbach	1.12
162	heimbach	1.1	enddo
163			enddo
164			enddo
165			enddo
166	heimbach	1.12	enddO
167	heimbach	1.1	#endif
168
169	heimbach	1.2	#ifdef ALLOW_TR10_CONTROL
170	heimbach	1.13	#ifdef ALLOW_PTRACERS
171	heimbach	1.2	c-- Temperature field.
172			il=ilnblnk( xx_tr1_file )
173	heimbach	1.17	write(fnamegeneric(1:80),'(2a,i10.10)')
174	heimbach	1.2	& xx_tr1_file(1:il),'.',optimcycle
175	heimbach	1.22	call active_read_xyz( fnamegeneric, tmpfld3d, 1,
176	heimbach	1.2	& doglobalread, ladinit, optimcycle,
177			& mythid, xx_tr1_dummy )
178
179			do bj = jtlo,jthi
180			do bi = itlo,ithi
181			do k = 1,nr
182			do j = jmin,jmax
183			do i = imin,imax
184	heimbach	1.13	ptracer(i,j,k,bi,bj,1) = ptracer(i,j,k,bi,bj,1) +
185	heimbach	1.2	& fac*tmpfld3d(i,j,k,bi,bj)
186			enddo
187			enddo
188			enddo
189			enddo
190			enddo
191			#endif
192	heimbach	1.13	#endif
193	heimbach	1.2
194	heimbach	1.14	#ifdef ALLOW_SST0_CONTROL
195			c-- sst0.
196			il=ilnblnk( xx_sst_file )
197	heimbach	1.17	write(fnamegeneric(1:80),'(2a,i10.10)')
198	heimbach	1.14	& xx_sst_file(1:il),'.',optimcycle
199	heimbach	1.22	call active_read_xy ( fnamegeneric, tmpfld2d, 1,
200	heimbach	1.14	& doglobalread, ladinit, optimcycle,
201			& mythid, xx_sst_dummy )
202			do bj = jtlo,jthi
203			do bi = itlo,ithi
204			do j = jmin,jmax
205			do i = imin,imax
206			cph sst(i,j,bi,bj) = sst(i,j,bi,bj) + tmpfld2d(i,j,bi,bj)
207	jmc	1.25	theta(i,j,1,bi,bj) = theta(i,j,1,bi,bj)
208	heimbach	1.14	& + tmpfld2d(i,j,bi,bj)
209			enddo
210			enddo
211			enddo
212			enddo
213			#endif
214
215			#ifdef ALLOW_SSS0_CONTROL
216			c-- sss0.
217			il=ilnblnk( xx_sss_file )
218	heimbach	1.17	write(fnamegeneric(1:80),'(2a,i10.10)')
219	heimbach	1.14	& xx_sss_file(1:il),'.',optimcycle
220	heimbach	1.22	call active_read_xy ( fnamegeneric, tmpfld2d, 1,
221	heimbach	1.14	& doglobalread, ladinit, optimcycle,
222			& mythid, xx_sss_dummy )
223			do bj = jtlo,jthi
224			do bi = itlo,ithi
225			do j = jmin,jmax
226			do i = imin,imax
227			cph sss(i,j,bi,bj) = sss(i,j,bi,bj) + tmpfld2d(i,j,bi,bj)
228			salt(i,j,1,bi,bj) = salt(i,j,1,bi,bj)
229			& + tmpfld2d(i,j,bi,bj)
230			enddo
231			enddo
232			enddo
233			enddo
234			#endif
235
236	heimbach	1.3	#ifdef ALLOW_DIFFKR_CONTROL
237			c-- diffkr.
238			il=ilnblnk( xx_diffkr_file )
239	heimbach	1.17	write(fnamegeneric(1:80),'(2a,i10.10)')
240	heimbach	1.3	& xx_diffkr_file(1:il),'.',optimcycle
241	heimbach	1.22	call active_read_xyz( fnamegeneric, tmpfld3d, 1,
242	heimbach	1.3	& doglobalread, ladinit, optimcycle,
243			& mythid, xx_diffkr_dummy )
244			do bj = jtlo,jthi
245			do bi = itlo,ithi
246			do k = 1,nr
247			do j = jmin,jmax
248			do i = imin,imax
249			diffkr(i,j,k,bi,bj) = diffkr(i,j,k,bi,bj) +
250			& tmpfld3d(i,j,k,bi,bj)
251			enddo
252			enddo
253			enddo
254			enddo
255			enddo
256			#endif
257
258			#ifdef ALLOW_KAPGM_CONTROL
259			c-- kapgm.
260			il=ilnblnk( xx_kapgm_file )
261	heimbach	1.17	write(fnamegeneric(1:80),'(2a,i10.10)')
262	heimbach	1.3	& xx_kapgm_file(1:il),'.',optimcycle
263	heimbach	1.22	call active_read_xyz( fnamegeneric, tmpfld3d, 1,
264	heimbach	1.3	& doglobalread, ladinit, optimcycle,
265			& mythid, xx_kapgm_dummy )
266			do bj = jtlo,jthi
267			do bi = itlo,ithi
268			do k = 1,nr
269			do j = jmin,jmax
270			do i = imin,imax
271			kapgm(i,j,k,bi,bj) = kapgm(i,j,k,bi,bj) +
272			& tmpfld3d(i,j,k,bi,bj)
273			enddo
274			enddo
275			enddo
276			enddo
277			enddo
278			#endif
279
280	heimbach	1.6	#ifdef ALLOW_EFLUXY0_CONTROL
281			c-- y-component EP-flux field.
282			il=ilnblnk( xx_efluxy_file )
283	heimbach	1.17	write(fnamegeneric(1:80),'(2a,i10.10)')
284	heimbach	1.6	& xx_efluxy_file(1:il),'.',optimcycle
285	heimbach	1.22	call active_read_xyz( fnamegeneric, tmpfld3d, 1,
286	heimbach	1.6	& doglobalread, ladinit, optimcycle,
287			& mythid, xx_efluxy_dummy )
288
289			do bj = jtlo,jthi
290			do bi = itlo,ithi
291			do k = 1,nr
292			do j = jmin,jmax
293			do i = imin,imax
294			EfluxY(i,j,k,bi,bj) = EfluxY(i,j,k,bi,bj)
295			& - fac*tmpfld3d(i,j,k,bi,bj)
296			& *maskS(i,j,k,bi,bj)
297			cph EfluxY(i,j,k,bi,bj) = EfluxY(i,j,k,bi,bj)
298			cph & - rSphere*cosFacU(J,bi,bj)
299			cph & factmpfld3d(i,j,k,bi,bj)
300			enddo
301			enddo
302			enddo
303			enddo
304			enddo
305			#endif
306
307			#ifdef ALLOW_EFLUXP0_CONTROL
308			c-- p-component EP-flux field.
309			il=ilnblnk( xx_efluxp_file )
310	heimbach	1.17	write(fnamegeneric(1:80),'(2a,i10.10)')
311	heimbach	1.6	& xx_efluxp_file(1:il),'.',optimcycle
312	heimbach	1.22	call active_read_xyz( fnamegeneric, tmpfld3d, 1,
313	heimbach	1.6	& doglobalread, ladinit, optimcycle,
314			& mythid, xx_efluxp_dummy )
315
316			do bj = jtlo,jthi
317			do bi = itlo,ithi
318			do k = 1,nr
319			do j = jmin,jmax
320			do i = imin,imax
321			EfluxP(i,j,k,bi,bj) = EfluxP(i,j,k,bi,bj)
322			& + fCori(i,j,bi,bj)
323			& factmpfld3d(i,j,k,bi,bj)
324			& *hFacV(i,j,k,bi,bj)
325			cph EfluxP(i,j,k,bi,bj) = EfluxP(i,j,k,bi,bj)
326			cph & + fCori(i,j,bi,bj)
327			cph & rSpherecosFacU(J,bi,bj)
328			cph & factmpfld3d(i,j,k,bi,bj)
329			enddo
330			enddo
331			enddo
332			enddo
333			enddo
334			#endif
335
336	heimbach	1.7	#ifdef ALLOW_BOTTOMDRAG_CONTROL
337			c-- bottom drag
338			il=ilnblnk( xx_bottomdrag_file )
339	heimbach	1.17	write(fnamegeneric(1:80),'(2a,i10.10)')
340	heimbach	1.7	& xx_bottomdrag_file(1:il),'.',optimcycle
341	heimbach	1.22	call active_read_xy ( fnamegeneric, tmpfld2d, 1,
342	heimbach	1.7	& doglobalread, ladinit, optimcycle,
343			& mythid, xx_bottomdrag_dummy )
344			do bj = jtlo,jthi
345			do bi = itlo,ithi
346			do j = jmin,jmax
347			do i = imin,imax
348	jmc	1.25	bottomdragfld(i,j,bi,bj) = bottomdragfld(i,j,bi,bj)
349	heimbach	1.7	& + tmpfld2d(i,j,bi,bj)
350			enddo
351			enddo
352			enddo
353			enddo
354			#endif
355
356	heimbach	1.16	#ifdef ALLOW_EDTAUX_CONTROL
357	heimbach	1.15	c-- zonal eddy stress : edtaux
358			il=ilnblnk( xx_edtaux_file )
359	heimbach	1.17	write(fnamegeneric(1:80),'(2a,i10.10)')
360	heimbach	1.15	& xx_edtaux_file(1:il),'.',optimcycle
361	heimbach	1.17	call active_read_xyz( fnamegeneric, tmpfld3d, 1,
362	heimbach	1.15	& doglobalread, ladinit, optimcycle,
363			& mythid, xx_edtaux_dummy )
364			do bj = jtlo,jthi
365			do bi = itlo,ithi
366			do k = 1,nr
367			do j = jmin,jmax
368			do i = imin,imax
369	heimbach	1.20	eddyTauX(i,j,k,bi,bj) = eddyTauX(i,j,k,bi,bj) +
370			& fCori(i,j,bi,bj)*tmpfld3d(i,j,k,bi,bj)
371	heimbach	1.15	enddo
372			enddo
373			enddo
374			enddo
375			enddo
376			#endif
377
378			#ifdef ALLOW_EDTAUY_CONTROL
379			c-- meridional eddy stress : edtauy
380			il=ilnblnk( xx_edtauy_file )
381	heimbach	1.17	write(fnamegeneric(1:80),'(2a,i10.10)')
382	heimbach	1.15	& xx_edtauy_file(1:il),'.',optimcycle
383	heimbach	1.17	call active_read_xyz( fnamegeneric, tmpfld3d, 1,
384	heimbach	1.15	& doglobalread, ladinit, optimcycle,
385			& mythid, xx_edtauy_dummy )
386			do bj = jtlo,jthi
387			do bi = itlo,ithi
388			do k = 1,nr
389			do j = jmin,jmax
390			do i = imin,imax
391	heimbach	1.20	eddyTauY(i,j,k,bi,bj) = eddyTauY(i,j,k,bi,bj) +
392			& fCoriG(i,j,bi,bj)*tmpfld3d(i,j,k,bi,bj)
393	heimbach	1.15	enddo
394			enddo
395			enddo
396			enddo
397			enddo
398			#endif
399	heimbach	1.1
400	heimbach	1.17	#ifdef ALLOW_UVEL0_CONTROL
401			c-- initial zonal velocity
402			il=ilnblnk( xx_uvel_file )
403			write(fnamegeneric(1:80),'(2a,i10.10)')
404			& xx_uvel_file(1:il),'.',optimcycle
405			call active_read_xyz( fnamegeneric, tmpfld3d, 1,
406			& doglobalread, ladinit, optimcycle,
407			& mythid, xx_uvel_dummy )
408			do bj = jtlo,jthi
409			do bi = itlo,ithi
410			do k = 1,nr
411			do j = jmin,jmax
412			do i = imin,imax
413	utke	1.23	#ifdef ALLOW_AUTODIFF_OPENAD
414	heimbach	1.19	uVel(i,j,k,bi,bj) = uVel(i,j,k,bi,bj) +
415			& fac*xx_uvel(i,j,k,bi,bj)
416			#else
417	heimbach	1.17	uVel(i,j,k,bi,bj) = uVel(i,j,k,bi,bj) +
418	heimbach	1.19	& fac*tmpfld3d(i,j,k,bi,bj)
419			#endif
420	heimbach	1.17	enddo
421			enddo
422			enddo
423			enddo
424			enddo
425			#endif
426
427			#ifdef ALLOW_VVEL0_CONTROL
428			c-- initial merid. velocity
429			il=ilnblnk( xx_vvel_file )
430			write(fnamegeneric(1:80),'(2a,i10.10)')
431			& xx_vvel_file(1:il),'.',optimcycle
432			call active_read_xyz( fnamegeneric, tmpfld3d, 1,
433			& doglobalread, ladinit, optimcycle,
434			& mythid, xx_vvel_dummy )
435			do bj = jtlo,jthi
436			do bi = itlo,ithi
437			do k = 1,nr
438			do j = jmin,jmax
439			do i = imin,imax
440	utke	1.23	#ifdef ALLOW_AUTODIFF_OPENAD
441	heimbach	1.17	vVel(i,j,k,bi,bj) = vVel(i,j,k,bi,bj) +
442	heimbach	1.19	& fac*xx_vvel(i,j,k,bi,bj)
443			#else
444			vVel(i,j,k,bi,bj) = vVel(i,j,k,bi,bj) +
445			& fac*tmpfld3d(i,j,k,bi,bj)
446			#endif
447	heimbach	1.17	enddo
448			enddo
449			enddo
450			enddo
451			enddo
452			#endif
453
454			#ifdef ALLOW_ETAN0_CONTROL
455			c-- initial Eta.
456			il=ilnblnk( xx_etan_file )
457			write(fnamegeneric(1:80),'(2a,i10.10)')
458			& xx_etan_file(1:il),'.',optimcycle
459	heimbach	1.22	call active_read_xy ( fnamegeneric, tmpfld2d, 1,
460	heimbach	1.17	& doglobalread, ladinit, optimcycle,
461			& mythid, xx_etan_dummy )
462			do bj = jtlo,jthi
463			do bi = itlo,ithi
464			do j = jmin,jmax
465			do i = imin,imax
466	utke	1.23	#ifdef ALLOW_AUTODIFF_OPENAD
467	heimbach	1.19	etaN(i,j,bi,bj) = etaN(i,j,bi,bj) +
468			& fac*xx_etan(i,j,bi,bj)
469			#else
470	jmc	1.25	etaN(i,j,bi,bj) = etaN(i,j,bi,bj) +
471	heimbach	1.19	& fac*tmpfld2d(i,j,bi,bj)
472			#endif
473	heimbach	1.17	enddo
474			enddo
475			enddo
476			enddo
477			#endif
478
479			#ifdef ALLOW_RELAXSST_CONTROL
480			c-- SST relaxation coefficient.
481			il=ilnblnk( xx_relaxsst_file )
482			write(fnamegeneric(1:80),'(2a,i10.10)')
483			& xx_relaxsst_file(1:il),'.',optimcycle
484	heimbach	1.22	call active_read_xy ( fnamegeneric, tmpfld2d, 1,
485	heimbach	1.17	& doglobalread, ladinit, optimcycle,
486			& mythid, xx_relaxsst_dummy )
487			do bj = jtlo,jthi
488			do bi = itlo,ithi
489			do j = jmin,jmax
490			do i = imin,imax
491	jmc	1.25	lambdaThetaClimRelax(i,j,bi,bj) =
492			& lambdaThetaClimRelax(i,j,bi,bj)
493	heimbach	1.17	& + tmpfld2d(i,j,bi,bj)
494			enddo
495			enddo
496			enddo
497			enddo
498			#endif
499
500			#ifdef ALLOW_RELAXSSS_CONTROL
501			c-- SSS relaxation coefficient.
502			il=ilnblnk( xx_relaxsss_file )
503			write(fnamegeneric(1:80),'(2a,i10.10)')
504			& xx_relaxsss_file(1:il),'.',optimcycle
505	heimbach	1.22	call active_read_xy ( fnamegeneric, tmpfld2d, 1,
506	heimbach	1.17	& doglobalread, ladinit, optimcycle,
507			& mythid, xx_relaxsss_dummy )
508			do bj = jtlo,jthi
509			do bi = itlo,ithi
510			do j = jmin,jmax
511			do i = imin,imax
512	jmc	1.25	lambdaSaltClimRelax(i,j,bi,bj) =
513			& lambdaSaltClimRelax(i,j,bi,bj)
514	heimbach	1.17	& + tmpfld2d(i,j,bi,bj)
515			enddo
516			enddo
517			enddo
518			enddo
519			#endif
520
521	heimbach	1.24	#ifdef ALLOW_SEAICE
522			call seaice_ctrl_map_ini( mythid )
523			#endif
524
525	heimbach	1.1	c-- Update the tile edges.
526
527	heimbach	1.14	#if (defined (ALLOW_THETA0_CONTROL) \|\| defined (ALLOW_SST0_CONTROL))
528	heimbach	1.1	_EXCH_XYZ_R8( theta, mythid )
529			#endif
530	heimbach	1.14	#if (defined (ALLOW_SALT0_CONTROL) \|\| defined (ALLOW_SSS0_CONTROL))
531	heimbach	1.1	_EXCH_XYZ_R8( salt, mythid )
532	heimbach	1.2	#endif
533			#ifdef ALLOW_TR10_CONTROL
534	heimbach	1.13	#ifdef ALLOW_PTRACERS
535			_EXCH_XYZ_R8(pTracer(1-Olx,1-Oly,1,1,1,1),myThid)
536			#endif
537	heimbach	1.1	#endif
538	heimbach	1.3	#ifdef ALLOW_DIFFKR_CONTROL
539			_EXCH_XYZ_R8( diffkr, mythid)
540			#endif
541			#ifdef ALLOW_KAPGM_CONTROL
542			_EXCH_XYZ_R8( kapgm, mythid)
543	heimbach	1.6	#endif
544			#ifdef ALLOW_EFLUXY0_CONTROL
545			_EXCH_XYZ_R8( EfluxY, mythid )
546			#endif
547			#ifdef ALLOW_EFLUXP0_CONTROL
548			_EXCH_XYZ_R8( EfluxP, mythid )
549	heimbach	1.7	#endif
550			#ifdef ALLOW_BOTTOMDRAG_CONTROL
551			_EXCH_XY_R8( bottomdragfld, mythid )
552	heimbach	1.3	#endif
553
554	heimbach	1.15	#if (defined (ALLOW_EDTAUX_CONTROL) && defined (ALLOW_EDTAUY_CONTROL))
555	heimbach	1.20	CALL EXCH_UV_XYZ_RS(eddyTauX,eddyTauY,.TRUE.,myThid)
556	heimbach	1.15	#elif (defined (ALLOW_EDTAUX_CONTROL) \|\| defined (ALLOW_EDTAUY_CONTROL))
557			STOP 'ctrl_map_forcing: need BOTH ALLOW_EDTAU[X,Y]_CONTROL'
558			#endif
559	heimbach	1.1
560	heimbach	1.17	#ifdef ALLOW_UVEL0_CONTROL
561			_EXCH_XYZ_R8( uVel, mythid)
562			#endif
563
564			#ifdef ALLOW_VVEL0_CONTROL
565			_EXCH_XYZ_R8( vVel, mythid)
566			#endif
567
568			#ifdef ALLOW_ETAN0_CONTROL
569			_EXCH_XY_R8( etaN, mythid )
570			#endif
571
572			#ifdef ALLOW_RELAXSST_CONTROL
573			_EXCH_XY_R4( lambdaThetaClimRelax, mythid )
574			#endif
575
576			#ifdef ALLOW_RELAXSSS_CONTROL
577			_EXCH_XY_R4( lambdaThetaClimRelax, mythid )
578			#endif
579
580	heimbach	1.1	return
581			end
582