model/src/ini_forcing.F

C $Header: /u/gcmpack/MITgcm/model/src/ini_forcing.F,v 1.37 2005/01/04 15:57:22 jmc 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
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
       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
      _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
      _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

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.38	C $Header: /u/gcmpack/MITgcm/model/src/ini_forcing.F,v 1.37 2005/01/04 15:57:22 jmc 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	adcroft	1.17	INTEGER I, J
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			ENDDO
91			ENDDO
92	adcroft	1.15	C
93	heimbach	1.38	DO bj = myByLo(myThid), myByHi(myThid)
94			DO bi = myBxLo(myThid), myBxHi(myThid)
95			DO J=1-Oly,sNy+Oly
96			DO I=1-Olx,sNx+Olx
97			IF ( doThetaClimRelax .AND.
98			& abs(yC(i,j,bi,bj)).LE.latBandClimRelax ) THEN
99			lambdaThetaClimRelax(I,J,bi,bj) = 1./tauThetaClimRelax
100			ELSE
101			lambdaThetaClimRelax(I,J,bi,bj) = 0.D0
102			ENDIF
103			IF ( doSaltClimRelax .AND.
104			& abs(yC(i,j,bi,bj)).LE.latBandClimRelax ) THEN
105			lambdaSaltClimRelax(I,J,bi,bj) = 1./tauSaltClimRelax
106			ELSE
107			lambdaSaltClimRelax(I,J,bi,bj) = 0.D0
108			ENDIF
109			ENDDO
110			ENDDO
111			ENDDO
112			ENDDO
113			C
114	adcroft	1.16	_BEGIN_MASTER(myThid)
115			IF ( zonalWindFile .NE. ' ' ) THEN
116			CALL READ_FLD_XY_RS( zonalWindFile, ' ', fu, 0, myThid )
117			ENDIF
118			IF ( meridWindFile .NE. ' ' ) THEN
119			CALL READ_FLD_XY_RS( meridWindFile, ' ', fv, 0, myThid )
120			ENDIF
121			IF ( surfQFile .NE. ' ' ) THEN
122			CALL READ_FLD_XY_RS( surfQFile, ' ', Qnet, 0, myThid )
123	jmc	1.34	ELSEIF ( surfQnetFile .NE. ' ' ) THEN
124			CALL READ_FLD_XY_RS( surfQnetFile, ' ', Qnet, 0, myThid )
125	adcroft	1.16	ENDIF
126			IF ( EmPmRfile .NE. ' ' ) THEN
127			CALL READ_FLD_XY_RS( EmPmRfile, ' ', EmPmR, 0, myThid )
128			ENDIF
129	jmc	1.37	IF ( saltFluxFile .NE. ' ' ) THEN
130			CALL READ_FLD_XY_RS( saltFluxFile, ' ', saltFlux, 0, myThid )
131			ENDIF
132	adcroft	1.16	IF ( thetaClimFile .NE. ' ' ) THEN
133			CALL READ_FLD_XY_RS( thetaClimFile, ' ', SST, 0, myThid )
134			ENDIF
135			IF ( saltClimFile .NE. ' ' ) THEN
136			CALL READ_FLD_XY_RS( saltClimFile, ' ', SSS, 0, myThid )
137			ENDIF
138	heimbach	1.38	IF ( lambdaThetaFile .NE. ' ' ) THEN
139			CALL READ_FLD_XY_RS( lambdaThetaFile, ' ',
140			& lambdaThetaClimRelax, 0, myThid )
141			ENDIF
142			IF ( lambdaSaltFile .NE. ' ' ) THEN
143			CALL READ_FLD_XY_RS( lambdaSaltFile, ' ',
144			& lambdaSaltClimRelax, 0, myThid )
145			ENDIF
146	heimbach	1.18	#ifdef SHORTWAVE_HEATING
147			IF ( surfQswFile .NE. ' ' ) THEN
148			CALL READ_FLD_XY_RS( surfQswFile, ' ', Qsw, 0, myThid )
149	jmc	1.34	IF ( surfQFile .NE. ' ' ) THEN
150			C- Qnet is now (after c54) the net Heat Flux (including SW)
151	jmc	1.36	DO bj = 1,nSy
152			DO bi = 1,nSx
153			DO j=1-OLy,sNy+OLy
154			DO i=1-OLx,sNx+OLx
155	jmc	1.34	Qnet(i,j,bi,bj) = Qnet(i,j,bi,bj) + Qsw(i,j,bi,bj)
156	jmc	1.36	ENDDO
157	jmc	1.34	ENDDO
158			ENDDO
159	jmc	1.36	ENDDO
160	jmc	1.34	ENDIF
161	heimbach	1.18	ENDIF
162			#endif
163	mlosch	1.23	#ifdef ATMOSPHERIC_LOADING
164			IF ( pLoadFile .NE. ' ' ) THEN
165			CALL READ_FLD_XY_RS( pLoadFile, ' ', pload, 0, myThid )
166			ENDIF
167			#endif
168	adcroft	1.16	_END_MASTER(myThid)
169			C
170	adcroft	1.15	_EXCH_XY_R4(fu , myThid )
171			_EXCH_XY_R4(fv , myThid )
172			_EXCH_XY_R4(Qnet , myThid )
173			_EXCH_XY_R4(EmPmR, myThid )
174	jmc	1.37	_EXCH_XY_R4( saltFlux, myThid )
175	adcroft	1.15	_EXCH_XY_R4(SST , myThid )
176			_EXCH_XY_R4(SSS , myThid )
177	heimbach	1.38	_EXCH_XY_R4(lambdaThetaClimRelax , myThid )
178			_EXCH_XY_R4(lambdaSaltClimRelax , myThid )
179	mlosch	1.23	#ifdef SHORTWAVE_HEATING
180	adcroft	1.16	_EXCH_XY_R4(Qsw , myThid )
181	mlosch	1.23	#endif
182			#ifdef ATMOSPHERIC_LOADING
183			_EXCH_XY_R4(pload , myThid )
184			C CALL PLOT_FIELD_XYRS( pload, 'S/R INI_FORCING pload',1,myThid)
185			#endif
186	adcroft	1.10
187	adcroft	1.15	C CALL PLOT_FIELD_XYRS( fu, 'S/R INI_FORCING FU',1,myThid)
188			C CALL PLOT_FIELD_XYRS( fv, 'S/R INI_FORCING FV',1,myThid)
189	cnh	1.1
190			RETURN
191			END