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

Contents of /MITgcm/pkg/dic/dic_surfforcing_init.F

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


Revision 1.33 - (show annotations) (download)
Fri Oct 7 21:36:39 2011 UTC (12 years, 8 months ago) by dfer
Branch: MAIN
CVS Tags: checkpoint63p, checkpoint63q, checkpoint63r, checkpoint63l, checkpoint63m, checkpoint63n, checkpoint63o, checkpoint63h, checkpoint63i, checkpoint63j, checkpoint63k, checkpoint63d, checkpoint63e, checkpoint63f, checkpoint63g
Changes since 1.32: +3 -3 lines
Remove subroutine CALC_PCO2_APPROX_CO3 from carbon_chem.F
and add carbonate computation/output to CALC_PCO2_APPROX

1 C $Header: /u/gcmpack/MITgcm/pkg/dic/dic_surfforcing_init.F,v 1.32 2011/05/05 22:23:27 stephd Exp $
2 C $Name: $
3
4 #include "DIC_OPTIONS.h"
5 #include "PTRACERS_OPTIONS.h"
6
7 CBOP
8 C !ROUTINE: DIC_SURFFORCING_INIT
9
10 C !INTERFACE: ==========================================================
11 SUBROUTINE DIC_SURFFORCING_INIT(
12 I myThid)
13
14 C !DESCRIPTION:
15 C Calculate first guess of pH
16
17 C !USES: ===============================================================
18 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_VARS.h"
26 #include "PTRACERS_SIZE.h"
27 #include "PTRACERS_PARAMS.h"
28 #include "PTRACERS_FIELDS.h"
29 #include "DIC_LOAD.h"
30
31 C !INPUT PARAMETERS: ===================================================
32 C myThid :: thread number
33 INTEGER myThid
34
35 #ifdef ALLOW_DIC
36
37 C !LOCAL VARIABLES: ====================================================
38 INTEGER i,j, kLev, it
39 INTEGER intimeP, intime0, intime1
40 _RL aWght, bWght, co3dummy
41 C Number of iterations for pCO2 solvers...
42 C Solubility relation coefficients
43 C local variables for carbon chem
44 _RL PTR_CO2(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
45 INTEGER iMin,iMax,jMin,jMax, bi, bj
46 _RL surfalk(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
47 _RL surfphos(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
48 _RL surfsi(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
49 _RL surftemp(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
50 _RL surfsalt(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
51 _RL surfdic(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
52 INTEGER iprt,jprt
53 LOGICAL pH_isLoaded
54 CEOP
55
56 cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
57
58 kLev=1
59
60 CALL DIC_ATMOS(0, startTime, nIter0, myThid )
61
62 ccccccccccccccccccccccccccccccccccccccccc
63 IF ( periodicExternalForcing ) THEN
64
65 c read in silica field
66 CALL LEF_ZERO( silica0,myThid )
67 CALL LEF_ZERO( silica1,myThid )
68
69 C-- Now calculate whether it is time to update the forcing arrays
70 CALL GET_PERIODIC_INTERVAL(
71 O intimeP, intime0, intime1, bWght, aWght,
72 I externForcingCycle, externForcingPeriod,
73 I deltaTclock, startTime, myThid )
74
75 _BARRIER
76 _BEGIN_MASTER(myThid)
77 WRITE(standardMessageUnit,'(A,I10,A,2(2I5,A))')
78 & ' DIC_SURFFORCING_INIT, it=', nIter0,
79 & ' : Reading new data, i0,i1=', intime0, intime1
80 _END_MASTER(myThid)
81
82 IF ( DIC_silicaFile .NE. ' ' ) THEN
83 CALL READ_REC_XY_RS( DIC_silicaFile,silica0,intime0,
84 & nIter0,myThid )
85 CALL READ_REC_XY_RS( DIC_silicaFile,silica1,intime1,
86 & nIter0,myThid )
87 ENDIF
88
89 #ifdef ALLOW_OFFLINE
90 IF ( useOffLine ) THEN
91 CALL OFFLINE_FIELDS_LOAD( startTime, nIter0, myThid )
92 ENDIF
93 #endif
94
95 _EXCH_XY_RS(silica0, myThid )
96 _EXCH_XY_RS(silica1, myThid )
97
98 IF ( DIC_silicaFile .NE. ' ' ) THEN
99 DO bj = myByLo(myThid), myByHi(myThid)
100 DO bi = myBxLo(myThid), myBxHi(myThid)
101 DO j=1-Oly,sNy+Oly
102 DO i=1-Olx,sNx+Olx
103 SILICA(i,j,bi,bj)= bWght*silica0(i,j,bi,bj)
104 & + aWght*silica1(i,j,bi,bj)
105 ENDDO
106 ENDDO
107 ENDDO
108 ENDDO
109 ENDIF
110
111 c end periodicExternalForcing
112 ENDIF
113
114 C =================================================================
115
116 jMin=1
117 jMax=sNy
118 iMin=1
119 iMax=sNx
120
121 DO bj=myByLo(myThid),myByHi(myThid)
122 DO bi=myBxLo(myThid),myBxHi(myThid)
123 DO j=1-OLy,sNy+OLy
124 DO i=1-Olx,sNx+OLx
125 pH(i,j,bi,bj) = 8. _d 0
126 ENDDO
127 ENDDO
128 ENDDO
129 ENDDO
130
131 DO bj=myByLo(myThid),myByHi(myThid)
132 DO bi=myBxLo(myThid),myBxHi(myThid)
133 DO j=1-OLy,sNy+OLy
134 DO i=1-OLx,sNx+OLx
135 ak0(i,j,bi,bj)=0. _d 0
136 ak1(i,j,bi,bj)=0. _d 0
137 ak2(i,j,bi,bj)=0. _d 0
138 akw(i,j,bi,bj)=0. _d 0
139 akb(i,j,bi,bj)=0. _d 0
140 akf(i,j,bi,bj)=0. _d 0
141 ak1p(i,j,bi,bj)=0. _d 0
142 ak2p(i,j,bi,bj)=0. _d 0
143 ak3p(i,j,bi,bj)=0. _d 0
144 aksi(i,j,bi,bj)=0. _d 0
145 fugf(i,j,bi,bj)=0. _d 0
146 ff(i,j,bi,bj)=0. _d 0
147 ft(i,j,bi,bj)=0. _d 0
148 st(i,j,bi,bj)=0. _d 0
149 bt(i,j,bi,bj)=0. _d 0
150 ENDDO
151 ENDDO
152 ENDDO
153 ENDDO
154
155 pH_isLoaded = .FALSE.
156 IF ( nIter0.GT.PTRACERS_Iter0 .OR.
157 & (nIter0.EQ.PTRACERS_Iter0 .AND. pickupSuff.NE.' ')
158 & ) THEN
159 C Read pH from a pickup file if needed
160 CALL DIC_READ_PICKUP(
161 O pH_isLoaded,
162 I nIter0, myThid )
163 ENDIF
164
165 DO bj=myByLo(myThid),myByHi(myThid)
166 DO bi=myBxLo(myThid),myBxHi(myThid)
167
168 C determine inorganic carbon chem coefficients
169 DO j=jMin,jMax
170 DO i=iMin,iMax
171
172 #ifdef DIC_BIOTIC
173 #ifdef DIC_BOUNDS
174 surfalk(i,j) = max(0.4 _d 0,
175 & min(10. _d 0,PTRACER(i,j,klev,bi,bj,2)))
176 & * maskC(i,j,kLev,bi,bj)
177 surfphos(i,j) = max(1.0 _d -11,
178 & min(1._d -1, PTRACER(i,j,klev,bi,bj,3)))
179 & * maskC(i,j,kLev,bi,bj)
180 #else
181 surfalk(i,j) = PTRACER(i,j,klev,bi,bj,2)
182 & * maskC(i,j,kLev,bi,bj)
183 surfphos(i,j)= PTRACER(i,j,klev,bi,bj,3)
184 & * maskC(i,j,kLev,bi,bj)
185 #endif
186 #else
187 surfalk(i,j) = 2.366595 _d 0 *salt(i,j,kLev,bi,bj)/35. _d 0
188 & * maskC(i,j,kLev,bi,bj)
189 surfphos(i,j)= 5.1225 _d -4 * maskC(i,j,kLev,bi,bj)
190 #endif
191 C FOR NON-INTERACTIVE Si
192 surfsi(i,j) = Silica(i,j,bi,bj) * maskC(i,j,kLev,bi,bj)
193 PTR_CO2(i,j,kLev) = PTRACER(i,j,klev,bi,bj,1)
194 & * maskC(i,j,kLev,bi,bj)
195 #ifdef DIC_BOUNDS
196 surftemp(i,j) = max(-4. _d 0,
197 & min(50. _d 0, theta(i,j,kLev,bi,bj)))
198 surfsalt(i,j) = max(4. _d 0,
199 & min(50. _d 0, salt(i,j,kLev,bi,bj)))
200 surfdic(i,j) = max(0.4 _d 0,
201 & min(10. _d 0, PTRACER(i,j,klev,bi,bj,1)))
202 #else
203 surftemp(i,j) = theta(i,j,kLev,bi,bj)
204 surfsalt(i,j) = salt(i,j,kLev,bi,bj)
205 surfdic(i,j) = PTR_CO2(i,j,kLev)
206 #endif
207 ENDDO
208 ENDDO
209
210 CALL CARBON_COEFFS(
211 I surftemp,surfsalt,
212 I bi,bj,iMin,iMax,jMin,jMax,myThid)
213
214 C====================================================================
215
216 IF ( .NOT.pH_isLoaded ) THEN
217 C set guess of pH for first step here
218
219 WRITE(standardMessageUnit,*) 'QQ: pCO2 approximation method'
220 c first approximation
221 C$TAF LOOP = parallel
222 DO j=jMin,jMax
223 C$TAF LOOP = parallel
224 DO i=iMin,iMax
225 IF ( maskC(i,j,kLev,bi,bj) .NE. 0. _d 0) THEN
226 C$TAF init dic_surf = static, 10
227 DO it=1,10
228 C$TAF STORE pH(i,j,bi,bj), PTR_CO2(i,j,kLev) = dic_surf
229 C$TAF STORE surfalk(i,j), surfphos(i,j), surfsi(i,j) = dic_surf
230 CALL CALC_PCO2_APPROX(
231 I surftemp(i,j),surfsalt(i,j),
232 I surfdic(i,j), surfphos(i,j),
233 I surfsi(i,j),surfalk(i,j),
234 I ak1(i,j,bi,bj),ak2(i,j,bi,bj),
235 I ak1p(i,j,bi,bj),ak2p(i,j,bi,bj),ak3p(i,j,bi,bj),
236 I aks(i,j,bi,bj),akb(i,j,bi,bj),akw(i,j,bi,bj),
237 I aksi(i,j,bi,bj),akf(i,j,bi,bj),
238 I ak0(i,j,bi,bj), fugf(i,j,bi,bj),
239 I ff(i,j,bi,bj),
240 I bt(i,j,bi,bj),st(i,j,bi,bj),ft(i,j,bi,bj),
241 U pH(i,j,bi,bj),pCO2(i,j,bi,bj),co3dummy,
242 I i,j,kLev,bi,bj, it , myThid )
243 ENDDO
244 ENDIF
245 ENDDO
246 ENDDO
247 iprt = MIN(20,sNx)
248 jprt = MIN(20,sNy)
249 WRITE(standardMessageUnit,*) 'QQ first guess pH',
250 & pH(iprt,jprt,bi,bj),
251 & theta(iprt,jprt,1,bi,bj), salt(iprt,jprt,1,bi,bj),
252 & PTR_CO2(iprt,jprt,1), surfphos(iprt,jprt),
253 & surfsi(iprt,jprt),surfalk(iprt,jprt)
254
255 ENDIF
256
257 C end bi,bj loops
258 ENDDO
259 ENDDO
260
261 #endif /* ALLOW_DIC */
262 RETURN
263 END

  ViewVC Help
Powered by ViewVC 1.1.22