/[MITgcm]/MITgcm/pkg/dic/dic_surfforcing.F
ViewVC logotype

Annotation of /MITgcm/pkg/dic/dic_surfforcing.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.5 - (hide annotations) (download)
Thu Feb 12 16:11:46 2004 UTC (20 years, 7 months ago) by stephd
Branch: MAIN
CVS Tags: checkpoint52l_pre, hrcube4, checkpoint52n_post, checkpoint52j_post, checkpoint53d_post, checkpoint54a_pre, checkpoint54a_post, checkpoint53c_post, checkpoint52l_post, checkpoint52k_post, checkpoint54b_post, checkpoint53b_pre, checkpoint52m_post, checkpoint53a_post, checkpoint54, checkpoint53b_post, checkpoint53, checkpoint53g_post, hrcube5, checkpoint52j_pre, checkpoint53f_post, checkpoint53d_pre
Changes since 1.4: +21 -11 lines
o clean up and add extra documentation

1 edhill 1.4 #include "DIC_OPTIONS.h"
2 stephd 1.1 #include "PTRACERS_OPTIONS.h"
3     #include "GCHEM_OPTIONS.h"
4    
5 stephd 1.5 CBOP
6     C !ROUTINE: DIC_SURFFORCING
7    
8     C !INTERFACE: ==========================================================
9 stephd 1.1 SUBROUTINE DIC_SURFFORCING( PTR_CO2 , GDC,
10     I bi,bj,imin,imax,jmin,jmax,
11     I myIter,myTime,myThid)
12    
13 stephd 1.5 C !DESCRIPTION:
14     C Calculate the carbon air-sea flux terms
15     C following external_forcing_dic.F (OCMIP run) from Mick
16    
17     C !USES: ===============================================================
18 stephd 1.1 IMPLICIT NONE
19     #include "SIZE.h"
20     #include "DYNVARS.h"
21     #include "EEPARAMS.h"
22     #include "PARAMS.h"
23     #include "GRID.h"
24     #include "FFIELDS.h"
25     #include "DIC_ABIOTIC.h"
26     #ifdef DIC_BIOTIC
27     #include "PTRACERS.h"
28     #endif
29    
30 stephd 1.5 C !INPUT PARAMETERS: ===================================================
31     C myThid :: thread number
32     C myIter :: current timestep
33     C myTime :: current time
34     c PTR_CO2 :: DIC tracer field
35 stephd 1.1 INTEGER myIter, myThid
36     _RL myTime
37     _RL PTR_CO2(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
38 stephd 1.5 INTEGER iMin,iMax,jMin,jMax, bi, bj
39    
40     C !OUTPUT PARAMETERS: ===================================================
41     c GDC :: tendency term due to air-sea exchange
42 stephd 1.1 _RL GDC(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
43    
44     #ifdef ALLOW_PTRACERS
45 stephd 1.5
46     C !LOCAL VARIABLES: ====================================================
47 stephd 1.2 INTEGER I,J, kLev, it
48 stephd 1.1 C Number of iterations for pCO2 solvers...
49     C Solubility relation coefficients
50     _RL SchmidtNoDIC(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
51     _RL pCO2sat(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
52     _RL Kwexch(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
53     C local variables for carbon chem
54     _RL surfalk(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
55     _RL surfphos(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
56     _RL surfsi(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
57     _RL VirtualFlux(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
58 stephd 1.5 CEOP
59 stephd 1.1
60     cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
61    
62     kLev=1
63    
64     C PRE-INDUSTRIAL STEADY STATE pCO2 = 278.0 ppmv
65     DO j=1-OLy,sNy+OLy
66     DO i=1-OLx,sNx+OLx
67     AtmospCO2(i,j,bi,bj)=278.0d-6
68     ENDDO
69     ENDDO
70    
71    
72     C =================================================================
73     C determine inorganic carbon chem coefficients
74     DO j=1-OLy,sNy+OLy
75     DO i=1-OLx,sNx+OLx
76    
77     #ifdef DIC_BIOTIC
78     cQQQQ check ptracer numbers
79     surfalk(i,j) = PTRACER(i,j,klev,bi,bj,2)
80     & * maskC(i,j,kLev,bi,bj)
81     surfphos(i,j) = PTRACER(i,j,klev,bi,bj,3)
82     & * maskC(i,j,kLev,bi,bj)
83     #else
84     surfalk(i,j) = 2.366595 * salt(i,j,kLev,bi,bj)/gsm_s
85     & * maskC(i,j,kLev,bi,bj)
86     surfphos(i,j) = 5.1225e-4 * maskC(i,j,kLev,bi,bj)
87     #endif
88     C FOR NON-INTERACTIVE Si
89 stephd 1.3 surfsi(i,j) = SILICA(i,j,bi,bj) * maskC(i,j,kLev,bi,bj)
90 stephd 1.1 ENDDO
91     ENDDO
92    
93     CALL CARBON_COEFFS(
94     I theta,salt,
95     I bi,bj,iMin,iMax,jMin,jMax)
96     C====================================================================
97    
98     c pCO2 solver...
99 stephd 1.3 C$TAF LOOP = parallel
100 stephd 1.1 DO j=1-OLy,sNy+OLy
101 stephd 1.3 C$TAF LOOP = parallel
102 stephd 1.1 DO i=1-OLx,sNx+OLx
103    
104     IF(maskC(i,j,kLev,bi,bj) .NE. 0.)THEN
105     CALL CALC_PCO2_APPROX(
106     I theta(i,j,kLev,bi,bj),salt(i,j,kLev,bi,bj),
107     I PTR_CO2(i,j,kLev), surfphos(i,j),
108     I surfsi(i,j),surfalk(i,j),
109     I ak1(i,j,bi,bj),ak2(i,j,bi,bj),
110     I ak1p(i,j,bi,bj),ak2p(i,j,bi,bj),ak3p(i,j,bi,bj),
111     I aks(i,j,bi,bj),akb(i,j,bi,bj),akw(i,j,bi,bj),
112     I aksi(i,j,bi,bj),akf(i,j,bi,bj),ff(i,j,bi,bj),
113     I bt(i,j,bi,bj),st(i,j,bi,bj),ft(i,j,bi,bj),
114     U pH(i,j,bi,bj),pCO2(i,j,bi,bj) )
115     ELSE
116     pCO2(i,j,bi,bj)=0. _d 0
117     END IF
118     ENDDO
119     ENDDO
120    
121     DO j=1-OLy,sNy+OLy
122     DO i=1-OLx,sNx+OLx
123    
124     IF (maskC(i,j,kLev,bi,bj).NE.0.) THEN
125     C calculate SCHMIDT NO. for CO2
126     SchmidtNoDIC(i,j) =
127     & sca1
128     & + sca2 * theta(i,j,kLev,bi,bj)
129     & + sca3 * theta(i,j,kLev,bi,bj)*theta(i,j,kLev,bi,bj)
130     & + sca4 * theta(i,j,kLev,bi,bj)*theta(i,j,kLev,bi,bj)
131     & *theta(i,j,kLev,bi,bj)
132    
133     C Determine surface flux (FDIC)
134     C first correct pCO2at for surface atmos pressure
135     pCO2sat(i,j) =
136     & AtmosP(i,j,bi,bj)*AtmospCO2(i,j,bi,bj)
137     c find exchange coefficient
138     c account for schmidt number and and varible piston velocity
139     Kwexch(i,j) =
140     & pisvel(i,j,bi,bj)
141     & / sqrt(SchmidtNoDIC(i,j)/660.0)
142     c OR use a constant coeff
143     c Kwexch(i,j) = 5e-5
144     c ice influence
145     cQQ Kwexch(i,j) =(1.d0-Fice(i,j,bi,bj))*Kwexch(i,j)
146    
147    
148     C Calculate flux in terms of DIC units using K0, solubility
149     C Flux = Vp * ([CO2sat] - [CO2])
150     C CO2sat = K0*pCO2atmos*P/P0
151     C Converting pCO2 to [CO2] using ff, as in CALC_PCO2
152 stephd 1.2 FluxCO2(i,j,bi,bj) =
153 stephd 1.1 & maskC(i,j,kLev,bi,bj)*Kwexch(i,j)*(
154     & ak0(i,j,bi,bj)*pCO2sat(i,j) -
155     & ff(i,j,bi,bj)*pCO2(i,j,bi,bj)
156     & )
157     ELSE
158 stephd 1.2 FluxCO2(i,j,bi,bj) = 0.
159 stephd 1.1 ENDIF
160     C convert flux (mol kg-1 m s-1) to (mol m-2 s-1)
161 stephd 1.2 FluxCO2(i,j,bi,bj) = FluxCO2(i,j,bi,bj)/permil
162 stephd 1.1
163     IF (maskC(i,j,kLev,bi,bj).NE.0.) THEN
164     c calculate virtual flux
165     c EminusPforV = dS/dt*(1/Sglob)
166     C NOTE: Be very careful with signs here!
167     C Positive EminusPforV => loss of water to atmos and increase
168     C in salinity. Thus, also increase in other surface tracers
169     C (i.e. positive virtual flux into surface layer)
170     C ...so here, VirtualFLux = dC/dt!
171     VirtualFlux(i,j)=gsm_DIC*surfaceTendencyS(i,j,bi,bj)/gsm_s
172     c OR
173     c let virtual flux be zero
174     c VirtualFlux(i,j)=0.d0
175     c
176     ELSE
177     VirtualFlux(i,j)=0. _d 0
178     ENDIF
179     ENDDO
180     ENDDO
181    
182     C update tendency
183     DO j=1-OLy,sNy+OLy
184     DO i=1-OLx,sNx+OLx
185     GDC(i,j)= maskC(i,j,kLev,bi,bj)*(
186 stephd 1.2 & FluxCO2(i,j,bi,bj)*recip_drF(kLev)
187 stephd 1.1 & + VirtualFlux(i,j)
188     & )
189     ENDDO
190     ENDDO
191    
192     #endif
193     RETURN
194     END

  ViewVC Help
Powered by ViewVC 1.1.22