pkg/exf/exf_mapfields.F

c $Header: /u/gcmpack/MITgcm/pkg/exf/exf_mapfields.F,v 1.12 2005/06/29 07:11:18 heimbach Exp $

#include "EXF_OPTIONS.h"

      subroutine exf_mapfields( mythid )

c     ==================================================================
c     SUBROUTINE exf_mapfields
c     ==================================================================
c
c     o Map external forcing fields (ustress, vstress, hflux, sflux,
c       swflux, apressure, climsss, climsst, etc.) onto ocean model
c       arrays (fu, fv, Qnet, EmPmR, Qsw, pload, sss, sst, etc.).
c       This routine is included to separate the ocean state estimation
c       tool as much as possible from the ocean model.  Unit and sign
c       conventions can be customized using variables exf_outscal_*,
c       which are set in exf_readparms.F.  See the header files
c       exf_fields.h and FFIELDS.h for definitions of the various input
c       and output fields and for default unit and sign convetions.
c
c     started: Christian Eckert eckert@mit.edu  09-Aug-1999
c
c     changed: Christian Eckert eckert@mit.edu  11-Jan-2000
c              - Restructured the code in order to create a package
c                for the MITgcmUV.
c
c              Christian Eckert eckert@mit.edu  12-Feb-2000
c              - Changed Routine names (package prefix: exf_)
c
c              Patrick Heimbach, heimbach@mit.edu  06-May-2000
c              - added and changed CPP flag structure for
c                ALLOW_BULKFORMULAE, ALLOW_ATM_TEMP
c
c              Patrick Heimbach, heimbach@mit.edu  23-May-2000
c              - sign change of ustress/vstress incorporated into
c                scaling factors exf_outscal_ust, exf_outscal_vst
c
c     mods for pkg/seaice: menemenlis@jpl.nasa.gov 20-Dec-2002
c
c     ==================================================================
c     SUBROUTINE exf_mapfields
c     ==================================================================

      implicit none

c     == global variables ==

#include "EEPARAMS.h"
#include "SIZE.h"
#include "FFIELDS.h"
#include "GRID.h"

#include "exf_param.h"
#include "exf_constants.h"
#include "exf_fields.h"
#include "exf_clim_param.h"
#include "exf_clim_fields.h"
#ifdef ALLOW_AUTODIFF_TAMC
# include "tamc.h"
# include "tamc_keys.h"
#endif
c     == routine arguments ==

c     mythid - thread number for this instance of the routine.

      integer mythid

c     == local variables ==

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

c     == end of interface ==

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

      do bj = jtlo,jthi
        do bi = itlo,ithi

#ifdef ALLOW_AUTODIFF_TAMC
          act1 = bi - myBxLo(myThid)
          max1 = myBxHi(myThid) - myBxLo(myThid) + 1
          act2 = bj - myByLo(myThid)
          max2 = myByHi(myThid) - myByLo(myThid) + 1
          act3 = myThid - 1
          max3 = nTx*nTy
          act4 = ikey_dynamics - 1
          ikey = (act1 + 1) + act2*max1
     &                      + act3*max1*max2
     &                      + act4*max1*max2*max3
#endif /* ALLOW_AUTODIFF_TAMC */

          do j = jmin,jmax
            do i = imin,imax
c            Heat flux.
             qnet(i,j,bi,bj) = exf_outscal_hflux*hflux(i,j,bi,bj)
            enddo
          enddo


          do j = jmin,jmax
            do i = imin,imax
c            Salt flux.
             empmr(i,j,bi,bj)= exf_outscal_sflux*sflux(i,j,bi,bj)
            enddo
          enddo

#ifdef ALLOW_AUTODIFF_TAMC
CADJ STORE ustress(:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte
#endif
          do j = jmin,jmax
            do i = imin,imax
c             Zonal wind stress.
              if (ustress(i,j,bi,bj).gt.windstressmax) then
                ustress(i,j,bi,bj)=windstressmax
              endif
            enddo
          enddo
#ifdef ALLOW_AUTODIFF_TAMC
CADJ STORE ustress(:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte
#endif
          do j = jmin,jmax
            do i = imin,imax
              if (ustress(i,j,bi,bj).lt.-windstressmax) then
                ustress(i,j,bi,bj)=-windstressmax
              endif
            enddo
          enddo
          do j = jmin,jmax
            do i = imin+1,imax
#if (defined (ALLOW_BULKFORMULAE) || defined (USE_EXF_INTERPOLATION))
c     Shift wind stresses calculated at C-points to W/S points
              fu(i,j,bi,bj) = exf_outscal_ustress*
     &              (ustress(i,j,bi,bj)+ustress(i-1,j,bi,bj))/2.*
     &              maskW(i,j,1,bi,bj)
#else
              fu(i,j,bi,bj) = exf_outscal_ustress*ustress(i,j,bi,bj)
#endif
            enddo
          enddo

#ifdef ALLOW_AUTODIFF_TAMC
CADJ STORE vstress(:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte
#endif
          do j = jmin,jmax
            do i = imin,imax
c             Meridional wind stress.
              if (vstress(i,j,bi,bj).gt.windstressmax) then
                vstress(i,j,bi,bj)=windstressmax
              endif
            enddo
          enddo
#ifdef ALLOW_AUTODIFF_TAMC
CADJ STORE vstress(:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte
#endif
          do j = jmin,jmax
            do i = imin,imax
              if (vstress(i,j,bi,bj).lt.-windstressmax) then
                vstress(i,j,bi,bj)=-windstressmax
              endif
            enddo
          enddo
          do j = jmin+1,jmax
            do i = imin,imax
#if (defined (ALLOW_BULKFORMULAE) || defined (USE_EXF_INTERPOLATION))
c     Shift wind stresses calculated at C-points to W/S points
              fv(i,j,bi,bj) = exf_outscal_vstress*
     &              (vstress(i,j,bi,bj)+vstress(i,j-1,bi,bj))/2.*
     &              maskS(i,j,1,bi,bj)
#else
              fv(i,j,bi,bj) = exf_outscal_vstress*vstress(i,j,bi,bj)
#endif
            enddo
          enddo

#ifdef SHORTWAVE_HEATING
c             Short wave radiative flux.
          do j = jmin,jmax
            do i = imin,imax
             qsw(i,j,bi,bj)  = exf_outscal_swflux*swflux(i,j,bi,bj)
            enddo
          enddo
#endif

#ifdef ALLOW_CLIMSST_RELAXATION
          do j = jmin,jmax
            do i = imin,imax
             sst(i,j,bi,bj)  = exf_outscal_sst*climsst(i,j,bi,bj)
            enddo
          enddo
#endif

#ifdef ALLOW_CLIMSSS_RELAXATION
          do j = jmin,jmax
            do i = imin,imax
             sss(i,j,bi,bj)  = exf_outscal_sss*climsss(i,j,bi,bj)
            enddo
          enddo
#endif

#ifdef ALLOW_CLIMTEMP_RELAXATION
          if ( climtempfile .NE. ' ' ) then
             do k = 1, Nr
                do j = jmin,jmax
                   do i = imin,imax
                      thetaStar(i,j,k,bi,bj)  = climtemp(i,j,k,bi,bj)
                   enddo
                enddo
             enddo
          endif
#endif

#ifdef ALLOW_CLIMSALT_RELAXATION
          if ( climsaltfile .NE. ' ' ) then
             do k = 1, Nr
                do j = jmin,jmax
                   do i = imin,imax
                      saltStar(i,j,k,bi,bj)  = climsalt(i,j,k,bi,bj)
                   enddo
                enddo
             enddo
          endif
#endif

#ifdef ATMOSPHERIC_LOADING
          do j = jmin,jmax
            do i = imin,imax
             pload(i,j,bi,bj)=exf_outscal_apressure*apressure(i,j,bi,bj)
            enddo
          enddo
#endif

        enddo
      enddo

c     Update the tile edges.

      _EXCH_XY_R4(  qnet, mythid )
      _EXCH_XY_R4( empmr, mythid )
c      _EXCH_XY_R4(    fu, mythid )
c      _EXCH_XY_R4(    fv, mythid )
       CALL EXCH_UV_XY_RS(fu, fv, .TRUE., myThid)
#ifdef SHORTWAVE_HEATING
      _EXCH_XY_R4(   qsw, mythid )
#endif
#ifdef ALLOW_CLIMSST_RELAXATION
      _EXCH_XY_R4(   sst, mythid )
#endif
#ifdef ALLOW_CLIMSSS_RELAXATION
      _EXCH_XY_R4(   sss, mythid )
#endif
#ifdef ALLOW_CLIMTEMP_RELAXATION
      _EXCH_XYZ_R4( thetaStar, mythid )
#endif
#ifdef ALLOW_CLIMSALT_RELAXATION
      _EXCH_XYZ_R4( saltStar, mythid )
#endif
#ifdef ATMOSPHERIC_LOADING
      _EXCH_XY_R4( pload, mythid )
#endif

      end
1	heimbach	1.13	c $Header: /u/gcmpack/MITgcm/pkg/exf/exf_mapfields.F,v 1.12 2005/06/29 07:11:18 heimbach Exp $
2	heimbach	1.1
3	edhill	1.7	#include "EXF_OPTIONS.h"
4	heimbach	1.1
5	heimbach	1.2	subroutine exf_mapfields( mythid )
6	heimbach	1.1
7			c ==================================================================
8	heimbach	1.2	c SUBROUTINE exf_mapfields
9	heimbach	1.1	c ==================================================================
10			c
11	dimitri	1.5	c o Map external forcing fields (ustress, vstress, hflux, sflux,
12			c swflux, apressure, climsss, climsst, etc.) onto ocean model
13			c arrays (fu, fv, Qnet, EmPmR, Qsw, pload, sss, sst, etc.).
14			c This routine is included to separate the ocean state estimation
15			c tool as much as possible from the ocean model. Unit and sign
16			c conventions can be customized using variables exf_outscal_*,
17			c which are set in exf_readparms.F. See the header files
18			c exf_fields.h and FFIELDS.h for definitions of the various input
19			c and output fields and for default unit and sign convetions.
20	heimbach	1.1	c
21			c started: Christian Eckert eckert@mit.edu 09-Aug-1999
22			c
23			c changed: Christian Eckert eckert@mit.edu 11-Jan-2000
24			c - Restructured the code in order to create a package
25			c for the MITgcmUV.
26			c
27			c Christian Eckert eckert@mit.edu 12-Feb-2000
28			c - Changed Routine names (package prefix: exf_)
29			c
30			c Patrick Heimbach, heimbach@mit.edu 06-May-2000
31			c - added and changed CPP flag structure for
32			c ALLOW_BULKFORMULAE, ALLOW_ATM_TEMP
33			c
34			c Patrick Heimbach, heimbach@mit.edu 23-May-2000
35			c - sign change of ustress/vstress incorporated into
36	dimitri	1.5	c scaling factors exf_outscal_ust, exf_outscal_vst
37	dimitri	1.4	c
38	dimitri	1.5	c mods for pkg/seaice: menemenlis@jpl.nasa.gov 20-Dec-2002
39	dimitri	1.4	c
40	heimbach	1.1	c ==================================================================
41	heimbach	1.2	c SUBROUTINE exf_mapfields
42	heimbach	1.1	c ==================================================================
43
44			implicit none
45
46			c == global variables ==
47
48			#include "EEPARAMS.h"
49			#include "SIZE.h"
50			#include "FFIELDS.h"
51	heimbach	1.8	#include "GRID.h"
52
53	heimbach	1.2	#include "exf_param.h"
54	heimbach	1.1	#include "exf_constants.h"
55			#include "exf_fields.h"
56	heimbach	1.13	#include "exf_clim_param.h"
57	heimbach	1.1	#include "exf_clim_fields.h"
58	heimbach	1.2	#ifdef ALLOW_AUTODIFF_TAMC
59			# include "tamc.h"
60			# include "tamc_keys.h"
61			#endif
62	heimbach	1.1	c == routine arguments ==
63
64			c mythid - thread number for this instance of the routine.
65
66			integer mythid
67
68			c == local variables ==
69
70			integer bi,bj
71	heimbach	1.13	integer i,j,k
72	heimbach	1.1	integer jtlo
73			integer jthi
74			integer itlo
75			integer ithi
76			integer jmin
77			integer jmax
78			integer imin
79			integer imax
80
81			c == end of interface ==
82
83			jtlo = mybylo(mythid)
84			jthi = mybyhi(mythid)
85			itlo = mybxlo(mythid)
86			ithi = mybxhi(mythid)
87			jmin = 1-oly
88			jmax = sny+oly
89			imin = 1-olx
90			imax = snx+olx
91
92			do bj = jtlo,jthi
93			do bi = itlo,ithi
94	heimbach	1.2
95			#ifdef ALLOW_AUTODIFF_TAMC
96			act1 = bi - myBxLo(myThid)
97			max1 = myBxHi(myThid) - myBxLo(myThid) + 1
98			act2 = bj - myByLo(myThid)
99			max2 = myByHi(myThid) - myByLo(myThid) + 1
100			act3 = myThid - 1
101			max3 = nTx*nTy
102			act4 = ikey_dynamics - 1
103			ikey = (act1 + 1) + act2*max1
104			& + act3max1max2
105			& + act4max1max2*max3
106			#endif /* ALLOW_AUTODIFF_TAMC */
107
108	heimbach	1.1	do j = jmin,jmax
109			do i = imin,imax
110	dimitri	1.5	c Heat flux.
111	heimbach	1.6	qnet(i,j,bi,bj) = exf_outscal_hflux*hflux(i,j,bi,bj)
112	heimbach	1.2	enddo
113			enddo
114	heimbach	1.1
115	heimbach	1.2
116			do j = jmin,jmax
117			do i = imin,imax
118	dimitri	1.5	c Salt flux.
119	heimbach	1.6	empmr(i,j,bi,bj)= exf_outscal_sflux*sflux(i,j,bi,bj)
120	heimbach	1.2	enddo
121			enddo
122	heimbach	1.1
123	heimbach	1.2	#ifdef ALLOW_AUTODIFF_TAMC
124			CADJ STORE ustress(:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte
125			#endif
126			do j = jmin,jmax
127			do i = imin,imax
128	heimbach	1.1	c Zonal wind stress.
129	mlosch	1.10	if (ustress(i,j,bi,bj).gt.windstressmax) then
130			ustress(i,j,bi,bj)=windstressmax
131	heimbach	1.2	endif
132			enddo
133			enddo
134			#ifdef ALLOW_AUTODIFF_TAMC
135			CADJ STORE ustress(:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte
136			#endif
137			do j = jmin,jmax
138			do i = imin,imax
139	mlosch	1.10	if (ustress(i,j,bi,bj).lt.-windstressmax) then
140			ustress(i,j,bi,bj)=-windstressmax
141	heimbach	1.2	endif
142			enddo
143			enddo
144			do j = jmin,jmax
145	heimbach	1.8	do i = imin+1,imax
146	dimitri	1.9	#if (defined (ALLOW_BULKFORMULAE) \|\| defined (USE_EXF_INTERPOLATION))
147	heimbach	1.8	c Shift wind stresses calculated at C-points to W/S points
148			fu(i,j,bi,bj) = exf_outscal_ustress*
149			& (ustress(i,j,bi,bj)+ustress(i-1,j,bi,bj))/2.*
150			& maskW(i,j,1,bi,bj)
151			#else
152			fu(i,j,bi,bj) = exf_outscal_ustress*ustress(i,j,bi,bj)
153			#endif
154	heimbach	1.2	enddo
155			enddo
156	heimbach	1.1
157	heimbach	1.2	#ifdef ALLOW_AUTODIFF_TAMC
158			CADJ STORE vstress(:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte
159			#endif
160			do j = jmin,jmax
161			do i = imin,imax
162	heimbach	1.1	c Meridional wind stress.
163	mlosch	1.10	if (vstress(i,j,bi,bj).gt.windstressmax) then
164			vstress(i,j,bi,bj)=windstressmax
165	heimbach	1.2	endif
166			enddo
167			enddo
168			#ifdef ALLOW_AUTODIFF_TAMC
169			CADJ STORE vstress(:,:,bi,bj) = comlev1_bibj, key=ikey, byte=isbyte
170			#endif
171			do j = jmin,jmax
172			do i = imin,imax
173	mlosch	1.10	if (vstress(i,j,bi,bj).lt.-windstressmax) then
174			vstress(i,j,bi,bj)=-windstressmax
175	heimbach	1.2	endif
176			enddo
177			enddo
178	heimbach	1.8	do j = jmin+1,jmax
179	heimbach	1.2	do i = imin,imax
180	dimitri	1.9	#if (defined (ALLOW_BULKFORMULAE) \|\| defined (USE_EXF_INTERPOLATION))
181	heimbach	1.8	c Shift wind stresses calculated at C-points to W/S points
182			fv(i,j,bi,bj) = exf_outscal_vstress*
183			& (vstress(i,j,bi,bj)+vstress(i,j-1,bi,bj))/2.*
184			& maskS(i,j,1,bi,bj)
185			#else
186			fv(i,j,bi,bj) = exf_outscal_vstress*vstress(i,j,bi,bj)
187			#endif
188	heimbach	1.2	enddo
189			enddo
190	heimbach	1.1
191	dimitri	1.5	#ifdef SHORTWAVE_HEATING
192	heimbach	1.1	c Short wave radiative flux.
193	heimbach	1.2	do j = jmin,jmax
194			do i = imin,imax
195	heimbach	1.6	qsw(i,j,bi,bj) = exf_outscal_swflux*swflux(i,j,bi,bj)
196	heimbach	1.2	enddo
197			enddo
198	heimbach	1.1	#endif
199
200			#ifdef ALLOW_CLIMSST_RELAXATION
201	heimbach	1.2	do j = jmin,jmax
202			do i = imin,imax
203	dimitri	1.5	sst(i,j,bi,bj) = exf_outscal_sst*climsst(i,j,bi,bj)
204	heimbach	1.2	enddo
205			enddo
206	heimbach	1.1	#endif
207
208			#ifdef ALLOW_CLIMSSS_RELAXATION
209	heimbach	1.2	do j = jmin,jmax
210			do i = imin,imax
211	dimitri	1.5	sss(i,j,bi,bj) = exf_outscal_sss*climsss(i,j,bi,bj)
212	heimbach	1.2	enddo
213			enddo
214	heimbach	1.1	#endif
215
216	heimbach	1.13	#ifdef ALLOW_CLIMTEMP_RELAXATION
217			if ( climtempfile .NE. ' ' ) then
218			do k = 1, Nr
219			do j = jmin,jmax
220			do i = imin,imax
221			thetaStar(i,j,k,bi,bj) = climtemp(i,j,k,bi,bj)
222			enddo
223			enddo
224			enddo
225			endif
226			#endif
227
228			#ifdef ALLOW_CLIMSALT_RELAXATION
229			if ( climsaltfile .NE. ' ' ) then
230			do k = 1, Nr
231			do j = jmin,jmax
232			do i = imin,imax
233			saltStar(i,j,k,bi,bj) = climsalt(i,j,k,bi,bj)
234			enddo
235			enddo
236			enddo
237			endif
238			#endif
239
240	heimbach	1.2	#ifdef ATMOSPHERIC_LOADING
241			do j = jmin,jmax
242			do i = imin,imax
243	dimitri	1.5	pload(i,j,bi,bj)=exf_outscal_apressure*apressure(i,j,bi,bj)
244	heimbach	1.1	enddo
245			enddo
246	heimbach	1.2	#endif
247
248	heimbach	1.1	enddo
249			enddo
250
251			c Update the tile edges.
252
253			_EXCH_XY_R4( qnet, mythid )
254			_EXCH_XY_R4( empmr, mythid )
255	cheisey	1.3	c _EXCH_XY_R4( fu, mythid )
256			c _EXCH_XY_R4( fv, mythid )
257			CALL EXCH_UV_XY_RS(fu, fv, .TRUE., myThid)
258	dimitri	1.5	#ifdef SHORTWAVE_HEATING
259	heimbach	1.1	_EXCH_XY_R4( qsw, mythid )
260			#endif
261			#ifdef ALLOW_CLIMSST_RELAXATION
262			_EXCH_XY_R4( sst, mythid )
263			#endif
264			#ifdef ALLOW_CLIMSSS_RELAXATION
265			_EXCH_XY_R4( sss, mythid )
266	heimbach	1.2	#endif
267	heimbach	1.13	#ifdef ALLOW_CLIMTEMP_RELAXATION
268			_EXCH_XYZ_R4( thetaStar, mythid )
269			#endif
270			#ifdef ALLOW_CLIMSALT_RELAXATION
271			_EXCH_XYZ_R4( saltStar, mythid )
272			#endif
273	heimbach	1.2	#ifdef ATMOSPHERIC_LOADING
274			_EXCH_XY_R4( pload, mythid )
275	heimbach	1.1	#endif
276
277			end