/[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.34 - (show annotations) (download)
Wed Aug 22 21:45:28 2012 UTC (11 years, 9 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint64a, checkpoint64c, checkpoint64b, checkpoint64e, checkpoint64d, checkpoint64g, checkpoint64f, checkpoint63s, checkpoint64
Changes since 1.33: +14 -16 lines
remove local (3D) array PTR_CO2 (not necessary)

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

  ViewVC Help
Powered by ViewVC 1.1.22