model/src/ini_forcing.F

C $Header: /u/gcmpack/MITgcm/model/src/ini_forcing.F,v 1.39 2005/11/07 18:26:02 cnh Exp $
C $Name:  $

#include "PACKAGES_CONFIG.h"
#include "CPP_OPTIONS.h"

CBOP
C     !ROUTINE: INI_FORCING
C     !INTERFACE:
      SUBROUTINE INI_FORCING( myThid )

C     !DESCRIPTION: \bv
C     *==========================================================*
C     | SUBROUTINE INI_FORCING                                    
C     | o Set model initial forcing fields.                       
C     *==========================================================*
C     \ev

C     !USES:
      IMPLICIT NONE
C     === Global variables ===
#include "SIZE.h"
#include "EEPARAMS.h"
#include "PARAMS.h"
#include "GRID.h"
#include "FFIELDS.h"

C     !INPUT/OUTPUT PARAMETERS:
C     == Routine arguments ==
C     myThid -  Number of this instance of INI_FORCING
      INTEGER myThid

C     !LOCAL VARIABLES:
C     == Local variables ==
C     bi,bj  - Loop counters
C     I,J
      INTEGER bi, bj
      INTEGER  I, J, K
CEOP

      _BARRIER

      DO bj = myByLo(myThid), myByHi(myThid)
       DO bi = myBxLo(myThid), myBxHi(myThid)
        DO j=1-OLy,sNy+OLy
         DO i=1-OLx,sNx+OLx
          fu              (i,j,bi,bj) = 0. _d 0
          fv              (i,j,bi,bj) = 0. _d 0
          Qnet            (i,j,bi,bj) = 0. _d 0
          EmPmR           (i,j,bi,bj) = 0. _d 0
          saltFlux        (i,j,bi,bj) = 0. _d 0
          SST             (i,j,bi,bj) = 0. _d 0
          SSS             (i,j,bi,bj) = 0. _d 0
          Qsw             (i,j,bi,bj) = 0. _d 0
#ifdef ATMOSPHERIC_LOADING
          pload           (i,j,bi,bj) = 0. _d 0
          sIceLoad        (i,j,bi,bj) = 0. _d 0
#endif
          surfaceForcingU(i,j,bi,bj) = 0. _d 0
          surfaceForcingV(i,j,bi,bj) = 0. _d 0
          surfaceForcingT(i,j,bi,bj) = 0. _d 0
          surfaceForcingS(i,j,bi,bj) = 0. _d 0
          surfaceForcingTice(i,j,bi,bj) = 0. _d 0
#ifndef ALLOW_EXF
          taux0           (i,j,bi,bj) = 0. _d 0
          taux1           (i,j,bi,bj) = 0. _d 0
          tauy0           (i,j,bi,bj) = 0. _d 0
          tauy1           (i,j,bi,bj) = 0. _d 0
          Qnet0           (i,j,bi,bj) = 0. _d 0
          Qnet1           (i,j,bi,bj) = 0. _d 0
          EmPmR0          (i,j,bi,bj) = 0. _d 0
          EmPmR1          (i,j,bi,bj) = 0. _d 0
          saltFlux0       (i,j,bi,bj) = 0. _d 0
          saltFlux1       (i,j,bi,bj) = 0. _d 0
          SST0            (i,j,bi,bj) = 0. _d 0
          SST1            (i,j,bi,bj) = 0. _d 0
          SSS0            (i,j,bi,bj) = 0. _d 0
          SSS1            (i,j,bi,bj) = 0. _d 0
#ifdef SHORTWAVE_HEATING          
          Qsw0            (i,j,bi,bj) = 0. _d 0
          Qsw1            (i,j,bi,bj) = 0. _d 0
#endif
#ifdef ATMOSPHERIC_LOADING
          pload0          (i,j,bi,bj) = 0. _d 0
          pload1          (i,j,bi,bj) = 0. _d 0
#endif
#endif
         ENDDO
        ENDDO
        DO k=1,Nr
         DO j=1-OLy,sNy+OLy
          DO i=1-OLx,sNx+OLx
#ifdef ALLOW_CLIMTEMP_RELAXATION
           thetaStar (i,j,k,bi,bj) = 0. _d 0
#endif
#ifdef ALLOW_CLIMSALT_RELAXATION
           saltStar  (i,j,k,bi,bj) = 0. _d 0
#endif
          ENDDO
         ENDDO
        ENDDO
       ENDDO
      ENDDO
C
      DO bj = myByLo(myThid), myByHi(myThid)
       DO bi = myBxLo(myThid), myBxHi(myThid)
        DO J=1-Oly,sNy+Oly
         DO I=1-Olx,sNx+Olx
          IF ( doThetaClimRelax .AND. 
     &         abs(yC(i,j,bi,bj)).LE.latBandClimRelax ) THEN
           lambdaThetaClimRelax(I,J,bi,bj) = 1./tauThetaClimRelax
          ELSE
           lambdaThetaClimRelax(I,J,bi,bj) = 0.D0
          ENDIF
          IF ( doSaltClimRelax .AND.
     &         abs(yC(i,j,bi,bj)).LE.latBandClimRelax ) THEN
           lambdaSaltClimRelax(I,J,bi,bj) = 1./tauSaltClimRelax
          ELSE
           lambdaSaltClimRelax(I,J,bi,bj) = 0.D0
          ENDIF
         ENDDO
        ENDDO
       ENDDO
      ENDDO
C
      _BARRIER
      _BEGIN_MASTER(myThid)
      IF ( zonalWindFile .NE. ' '  ) THEN
       CALL READ_FLD_XY_RS( zonalWindFile, ' ', fu, 0, myThid )
      ENDIF
      IF ( meridWindFile .NE. ' '  ) THEN
       CALL READ_FLD_XY_RS( meridWindFile, ' ', fv, 0, myThid )
      ENDIF
      IF ( surfQFile .NE. ' '  ) THEN
       CALL READ_FLD_XY_RS( surfQFile, ' ', Qnet, 0, myThid )
      ELSEIF ( surfQnetFile .NE. ' '  ) THEN
       CALL READ_FLD_XY_RS( surfQnetFile, ' ', Qnet, 0, myThid )
      ENDIF
      IF ( EmPmRfile .NE. ' '  ) THEN
       CALL READ_FLD_XY_RS( EmPmRfile, ' ', EmPmR, 0, myThid )
      ENDIF
      IF ( saltFluxFile .NE. ' '  ) THEN
       CALL READ_FLD_XY_RS( saltFluxFile, ' ', saltFlux, 0, myThid )
      ENDIF
      IF ( thetaClimFile .NE. ' '  ) THEN
       CALL READ_FLD_XY_RS( thetaClimFile, ' ', SST, 0, myThid )
      ENDIF
      IF ( saltClimFile .NE. ' '  ) THEN
       CALL READ_FLD_XY_RS( saltClimFile, ' ', SSS, 0, myThid )
      ENDIF
      IF ( lambdaThetaFile .NE. ' '  ) THEN
       CALL READ_FLD_XY_RS( lambdaThetaFile, ' ', 
     &  lambdaThetaClimRelax, 0, myThid )
      ENDIF
      IF ( lambdaSaltFile .NE. ' '  ) THEN
       CALL READ_FLD_XY_RS( lambdaSaltFile, ' ', 
     &  lambdaSaltClimRelax, 0, myThid )
      ENDIF
#ifdef SHORTWAVE_HEATING
      IF ( surfQswFile .NE. ' '  ) THEN
       CALL READ_FLD_XY_RS( surfQswFile, ' ', Qsw, 0, myThid )
       IF ( surfQFile .NE. ' '  ) THEN
C-     Qnet is now (after c54) the net Heat Flux (including SW)
        DO bj = 1,nSy
         DO bi = 1,nSx
          DO j=1-OLy,sNy+OLy
           DO i=1-OLx,sNx+OLx
            Qnet(i,j,bi,bj) = Qnet(i,j,bi,bj) + Qsw(i,j,bi,bj)
           ENDDO
          ENDDO
         ENDDO
        ENDDO
       ENDIF
      ENDIF
#endif
#ifdef ATMOSPHERIC_LOADING
      IF ( pLoadFile .NE. ' '  ) THEN
       CALL READ_FLD_XY_RS( pLoadFile, ' ', pload, 0, myThid )
      ENDIF
#endif
#ifdef ALLOW_CLIMTEMP_RELAXATION
      IF ( hydrogThetaFile .NE. ' '  ) THEN
       CALL READ_FLD_XYZ_RS( 
     &  hydrogThetaFile, ' ', thetaStar, 0, myThid )
      ENDIF
#endif
#ifdef ALLOW_CLIMSALT_RELAXATION
      IF ( hydrogSaltFile .NE. ' '  ) THEN
       CALL READ_FLD_XYZ_RS( 
     &  hydrogSaltFile, ' ', saltStar, 0, myThid )
      ENDIF
#endif

      _END_MASTER(myThid)
C
      _EXCH_XY_R4(fu   , myThid )
      _EXCH_XY_R4(fv   , myThid )
      _EXCH_XY_R4(Qnet , myThid )
      _EXCH_XY_R4(EmPmR, myThid )
      _EXCH_XY_R4( saltFlux, myThid )
      _EXCH_XY_R4(SST  , myThid )
      _EXCH_XY_R4(SSS  , myThid )
      _EXCH_XY_R4(lambdaThetaClimRelax , myThid )
      _EXCH_XY_R4(lambdaSaltClimRelax , myThid )
#ifdef SHORTWAVE_HEATING
      _EXCH_XY_R4(Qsw  , myThid )
#endif
#ifdef ATMOSPHERIC_LOADING
      _EXCH_XY_R4(pload  , myThid )
C     CALL PLOT_FIELD_XYRS( pload, 'S/R INI_FORCING pload',1,myThid)
#endif
#ifdef ALLOW_CLIMTEMP_RELAXATION
      _EXCH_XYZ_R4(thetaStar , myThid )
#endif
#ifdef ALLOW_CLIMSALT_RELAXATION
      _EXCH_XYZ_R4(saltStar , myThid )
#endif
C     CALL PLOT_FIELD_XYRS( fu, 'S/R INI_FORCING FU',1,myThid)
C     CALL PLOT_FIELD_XYRS( fv, 'S/R INI_FORCING FV',1,myThid)

      RETURN
      END
1	heimbach	1.40	C $Header: /u/gcmpack/MITgcm/model/src/ini_forcing.F,v 1.39 2005/11/07 18:26:02 cnh Exp $
2	cnh	1.22	C $Name: $
3	cnh	1.1
4	edhill	1.26	#include "PACKAGES_CONFIG.h"
5	cnh	1.12	#include "CPP_OPTIONS.h"
6	cnh	1.1
7	cnh	1.22	CBOP
8			C !ROUTINE: INI_FORCING
9			C !INTERFACE:
10	cnh	1.1	SUBROUTINE INI_FORCING( myThid )
11	cnh	1.22
12			C !DESCRIPTION: \bv
13			C ==========================================================
14			C \| SUBROUTINE INI_FORCING
15			C \| o Set model initial forcing fields.
16			C ==========================================================
17			C \ev
18
19			C !USES:
20	adcroft	1.14	IMPLICIT NONE
21	cnh	1.1	C === Global variables ===
22			#include "SIZE.h"
23			#include "EEPARAMS.h"
24			#include "PARAMS.h"
25			#include "GRID.h"
26			#include "FFIELDS.h"
27
28	cnh	1.22	C !INPUT/OUTPUT PARAMETERS:
29	cnh	1.1	C == Routine arguments ==
30			C myThid - Number of this instance of INI_FORCING
31			INTEGER myThid
32
33	cnh	1.22	C !LOCAL VARIABLES:
34	cnh	1.1	C == Local variables ==
35			C bi,bj - Loop counters
36	adcroft	1.17	C I,J
37	cnh	1.1	INTEGER bi, bj
38	heimbach	1.40	INTEGER I, J, K
39	cnh	1.22	CEOP
40	cnh	1.1
41	cnh	1.6	_BARRIER
42
43	adcroft	1.10	DO bj = myByLo(myThid), myByHi(myThid)
44			DO bi = myBxLo(myThid), myBxHi(myThid)
45			DO j=1-OLy,sNy+OLy
46			DO i=1-OLx,sNx+OLx
47	heimbach	1.18	fu (i,j,bi,bj) = 0. _d 0
48			fv (i,j,bi,bj) = 0. _d 0
49			Qnet (i,j,bi,bj) = 0. _d 0
50			EmPmR (i,j,bi,bj) = 0. _d 0
51	jmc	1.27	saltFlux (i,j,bi,bj) = 0. _d 0
52	heimbach	1.18	SST (i,j,bi,bj) = 0. _d 0
53			SSS (i,j,bi,bj) = 0. _d 0
54			Qsw (i,j,bi,bj) = 0. _d 0
55	mlosch	1.23	#ifdef ATMOSPHERIC_LOADING
56			pload (i,j,bi,bj) = 0. _d 0
57	jmc	1.33	sIceLoad (i,j,bi,bj) = 0. _d 0
58	mlosch	1.23	#endif
59	jmc	1.35	surfaceForcingU(i,j,bi,bj) = 0. _d 0
60			surfaceForcingV(i,j,bi,bj) = 0. _d 0
61			surfaceForcingT(i,j,bi,bj) = 0. _d 0
62			surfaceForcingS(i,j,bi,bj) = 0. _d 0
63			surfaceForcingTice(i,j,bi,bj) = 0. _d 0
64	edhill	1.30	#ifndef ALLOW_EXF
65	heimbach	1.24	taux0 (i,j,bi,bj) = 0. _d 0
66			taux1 (i,j,bi,bj) = 0. _d 0
67			tauy0 (i,j,bi,bj) = 0. _d 0
68			tauy1 (i,j,bi,bj) = 0. _d 0
69			Qnet0 (i,j,bi,bj) = 0. _d 0
70			Qnet1 (i,j,bi,bj) = 0. _d 0
71			EmPmR0 (i,j,bi,bj) = 0. _d 0
72			EmPmR1 (i,j,bi,bj) = 0. _d 0
73	jmc	1.37	saltFlux0 (i,j,bi,bj) = 0. _d 0
74			saltFlux1 (i,j,bi,bj) = 0. _d 0
75	heimbach	1.24	SST0 (i,j,bi,bj) = 0. _d 0
76			SST1 (i,j,bi,bj) = 0. _d 0
77			SSS0 (i,j,bi,bj) = 0. _d 0
78			SSS1 (i,j,bi,bj) = 0. _d 0
79			#ifdef SHORTWAVE_HEATING
80			Qsw0 (i,j,bi,bj) = 0. _d 0
81			Qsw1 (i,j,bi,bj) = 0. _d 0
82			#endif
83			#ifdef ATMOSPHERIC_LOADING
84			pload0 (i,j,bi,bj) = 0. _d 0
85			pload1 (i,j,bi,bj) = 0. _d 0
86			#endif
87			#endif
88	adcroft	1.10	ENDDO
89			ENDDO
90	heimbach	1.40	DO k=1,Nr
91			DO j=1-OLy,sNy+OLy
92			DO i=1-OLx,sNx+OLx
93			#ifdef ALLOW_CLIMTEMP_RELAXATION
94			thetaStar (i,j,k,bi,bj) = 0. _d 0
95			#endif
96			#ifdef ALLOW_CLIMSALT_RELAXATION
97			saltStar (i,j,k,bi,bj) = 0. _d 0
98			#endif
99			ENDDO
100			ENDDO
101			ENDDO
102	adcroft	1.10	ENDDO
103			ENDDO
104	adcroft	1.15	C
105	heimbach	1.38	DO bj = myByLo(myThid), myByHi(myThid)
106			DO bi = myBxLo(myThid), myBxHi(myThid)
107			DO J=1-Oly,sNy+Oly
108			DO I=1-Olx,sNx+Olx
109			IF ( doThetaClimRelax .AND.
110			& abs(yC(i,j,bi,bj)).LE.latBandClimRelax ) THEN
111			lambdaThetaClimRelax(I,J,bi,bj) = 1./tauThetaClimRelax
112			ELSE
113			lambdaThetaClimRelax(I,J,bi,bj) = 0.D0
114			ENDIF
115			IF ( doSaltClimRelax .AND.
116			& abs(yC(i,j,bi,bj)).LE.latBandClimRelax ) THEN
117			lambdaSaltClimRelax(I,J,bi,bj) = 1./tauSaltClimRelax
118			ELSE
119			lambdaSaltClimRelax(I,J,bi,bj) = 0.D0
120			ENDIF
121			ENDDO
122			ENDDO
123			ENDDO
124			ENDDO
125			C
126	cnh	1.39	_BARRIER
127	adcroft	1.16	_BEGIN_MASTER(myThid)
128			IF ( zonalWindFile .NE. ' ' ) THEN
129			CALL READ_FLD_XY_RS( zonalWindFile, ' ', fu, 0, myThid )
130			ENDIF
131			IF ( meridWindFile .NE. ' ' ) THEN
132			CALL READ_FLD_XY_RS( meridWindFile, ' ', fv, 0, myThid )
133			ENDIF
134			IF ( surfQFile .NE. ' ' ) THEN
135			CALL READ_FLD_XY_RS( surfQFile, ' ', Qnet, 0, myThid )
136	jmc	1.34	ELSEIF ( surfQnetFile .NE. ' ' ) THEN
137			CALL READ_FLD_XY_RS( surfQnetFile, ' ', Qnet, 0, myThid )
138	adcroft	1.16	ENDIF
139			IF ( EmPmRfile .NE. ' ' ) THEN
140			CALL READ_FLD_XY_RS( EmPmRfile, ' ', EmPmR, 0, myThid )
141			ENDIF
142	jmc	1.37	IF ( saltFluxFile .NE. ' ' ) THEN
143			CALL READ_FLD_XY_RS( saltFluxFile, ' ', saltFlux, 0, myThid )
144			ENDIF
145	adcroft	1.16	IF ( thetaClimFile .NE. ' ' ) THEN
146			CALL READ_FLD_XY_RS( thetaClimFile, ' ', SST, 0, myThid )
147			ENDIF
148			IF ( saltClimFile .NE. ' ' ) THEN
149			CALL READ_FLD_XY_RS( saltClimFile, ' ', SSS, 0, myThid )
150			ENDIF
151	heimbach	1.38	IF ( lambdaThetaFile .NE. ' ' ) THEN
152			CALL READ_FLD_XY_RS( lambdaThetaFile, ' ',
153			& lambdaThetaClimRelax, 0, myThid )
154			ENDIF
155			IF ( lambdaSaltFile .NE. ' ' ) THEN
156			CALL READ_FLD_XY_RS( lambdaSaltFile, ' ',
157			& lambdaSaltClimRelax, 0, myThid )
158			ENDIF
159	heimbach	1.18	#ifdef SHORTWAVE_HEATING
160			IF ( surfQswFile .NE. ' ' ) THEN
161			CALL READ_FLD_XY_RS( surfQswFile, ' ', Qsw, 0, myThid )
162	jmc	1.34	IF ( surfQFile .NE. ' ' ) THEN
163			C- Qnet is now (after c54) the net Heat Flux (including SW)
164	jmc	1.36	DO bj = 1,nSy
165			DO bi = 1,nSx
166			DO j=1-OLy,sNy+OLy
167			DO i=1-OLx,sNx+OLx
168	jmc	1.34	Qnet(i,j,bi,bj) = Qnet(i,j,bi,bj) + Qsw(i,j,bi,bj)
169	jmc	1.36	ENDDO
170	jmc	1.34	ENDDO
171			ENDDO
172	jmc	1.36	ENDDO
173	jmc	1.34	ENDIF
174	heimbach	1.18	ENDIF
175			#endif
176	mlosch	1.23	#ifdef ATMOSPHERIC_LOADING
177			IF ( pLoadFile .NE. ' ' ) THEN
178			CALL READ_FLD_XY_RS( pLoadFile, ' ', pload, 0, myThid )
179			ENDIF
180			#endif
181	heimbach	1.40	#ifdef ALLOW_CLIMTEMP_RELAXATION
182			IF ( hydrogThetaFile .NE. ' ' ) THEN
183			CALL READ_FLD_XYZ_RS(
184			& hydrogThetaFile, ' ', thetaStar, 0, myThid )
185			ENDIF
186			#endif
187			#ifdef ALLOW_CLIMSALT_RELAXATION
188			IF ( hydrogSaltFile .NE. ' ' ) THEN
189			CALL READ_FLD_XYZ_RS(
190			& hydrogSaltFile, ' ', saltStar, 0, myThid )
191			ENDIF
192			#endif
193
194	adcroft	1.16	_END_MASTER(myThid)
195			C
196	adcroft	1.15	_EXCH_XY_R4(fu , myThid )
197			_EXCH_XY_R4(fv , myThid )
198			_EXCH_XY_R4(Qnet , myThid )
199			_EXCH_XY_R4(EmPmR, myThid )
200	jmc	1.37	_EXCH_XY_R4( saltFlux, myThid )
201	adcroft	1.15	_EXCH_XY_R4(SST , myThid )
202			_EXCH_XY_R4(SSS , myThid )
203	heimbach	1.38	_EXCH_XY_R4(lambdaThetaClimRelax , myThid )
204			_EXCH_XY_R4(lambdaSaltClimRelax , myThid )
205	mlosch	1.23	#ifdef SHORTWAVE_HEATING
206	adcroft	1.16	_EXCH_XY_R4(Qsw , myThid )
207	mlosch	1.23	#endif
208			#ifdef ATMOSPHERIC_LOADING
209			_EXCH_XY_R4(pload , myThid )
210			C CALL PLOT_FIELD_XYRS( pload, 'S/R INI_FORCING pload',1,myThid)
211			#endif
212	heimbach	1.40	#ifdef ALLOW_CLIMTEMP_RELAXATION
213			_EXCH_XYZ_R4(thetaStar , myThid )
214			#endif
215			#ifdef ALLOW_CLIMSALT_RELAXATION
216			_EXCH_XYZ_R4(saltStar , myThid )
217			#endif
218	adcroft	1.15	C CALL PLOT_FIELD_XYRS( fu, 'S/R INI_FORCING FU',1,myThid)
219			C CALL PLOT_FIELD_XYRS( fv, 'S/R INI_FORCING FV',1,myThid)
220	cnh	1.1
221			RETURN
222			END