/[MITgcm]/MITgcm_contrib/bling/pkg/bling_carbonate_init.F
ViewVC logotype

Annotation of /MITgcm_contrib/bling/pkg/bling_carbonate_init.F

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


Revision 1.1 - (hide annotations) (download)
Sun Feb 28 21:53:37 2016 UTC (9 years, 4 months ago) by mmazloff
Branch: MAIN
Update to BLING version 2

1 mmazloff 1.1 C $Header: $
2     C $Name: $
3    
4     #include "BLING_OPTIONS.h"
5     #include "PTRACERS_OPTIONS.h"
6    
7     CBOP
8     subroutine BLING_CARBONATE_INIT( myThid )
9    
10     C ==========================================================
11     C | subroutine bling_surfforcing_init
12     C | o Calculate first guess of pH
13     C | Adapted from pkg/dic/dic_surfforcing_init
14     C ==========================================================
15    
16     implicit none
17    
18     C === Global variables ===
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 "BLING_VARS.h"
26     #include "PTRACERS_SIZE.h"
27     #include "PTRACERS_PARAMS.h"
28     #include "PTRACERS_FIELDS.h"
29     #include "BLING_LOAD.h"
30    
31     C === Routine arguments ===
32     C myThid :: thread Id. number
33     INTEGER myThid
34    
35     #ifdef ALLOW_BLING
36     C === Local variables ===
37     INTEGER i,j, k, it
38     INTEGER intimeP, intime0, intime1
39     _RL aWght, bWght, co3dummy
40     C Number of iterations for pCO2 solvers...
41     C Solubility relation coefficients
42     C local variables for carbon chem
43     INTEGER iMin,iMax,jMin,jMax, bi, bj
44     _RL alktmp(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
45     _RL phostmp(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
46     _RL sitmp(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
47     _RL thetatmp(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
48     _RL salttmp(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
49     _RL dictmp(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
50     INTEGER iprt,jprt
51     LOGICAL pH_isLoaded
52     CEOP
53    
54     IF ( periodicExternalForcing ) THEN
55    
56     c read in silica field
57     CALL LEF_ZERO( silica0,myThid )
58     CALL LEF_ZERO( silica1,myThid )
59    
60     C-- Now calculate whether it is time to update the forcing arrays
61     CALL GET_PERIODIC_INTERVAL(
62     O intimeP, intime0, intime1, bWght, aWght,
63     I externForcingCycle, externForcingPeriod,
64     I deltaTclock, startTime, myThid )
65    
66     _BARRIER
67     _BEGIN_MASTER(myThid)
68     WRITE(standardMessageUnit,'(A,I10,A,2(2I5,A))')
69     & ' BLING_SURFFORCING_INIT, it=', nIter0,
70     & ' : Reading new data, i0,i1=', intime0, intime1
71     _END_MASTER(myThid)
72    
73     IF ( BLING_silicaFile .NE. ' ' ) THEN
74     CALL READ_REC_XY_RS( BLING_silicaFile,silica0,intime0,
75     & nIter0,myThid )
76     CALL READ_REC_XY_RS( BLING_silicaFile,silica1,intime1,
77     & nIter0,myThid )
78     ENDIF
79    
80     #ifdef ALLOW_OFFLINE
81     IF ( useOffLine ) THEN
82     CALL OFFLINE_FIELDS_LOAD( startTime, nIter0, myThid )
83     ENDIF
84     #endif
85    
86     _EXCH_XY_RS(silica0, myThid )
87     _EXCH_XY_RS(silica1, myThid )
88    
89     IF ( BLING_silicaFile .NE. ' ' ) THEN
90     DO bj = myByLo(myThid), myByHi(myThid)
91     DO bi = myBxLo(myThid), myBxHi(myThid)
92     DO j=1-Oly,sNy+Oly
93     DO i=1-Olx,sNx+Olx
94     SILICA(i,j,bi,bj)= bWght*silica0(i,j,bi,bj)
95     & + aWght*silica1(i,j,bi,bj)
96     ENDDO
97     ENDDO
98     ENDDO
99     ENDDO
100     ENDIF
101    
102     c end periodicExternalForcing
103     ENDIF
104    
105     C =================================================================
106    
107     jMin=1
108     jMax=sNy
109     iMin=1
110     iMax=sNx
111    
112     DO k=1,Nr
113     DO bj=myByLo(myThid),myByHi(myThid)
114     DO bi=myBxLo(myThid),myBxHi(myThid)
115     DO j=1-OLy,sNy+OLy
116     DO i=1-Olx,sNx+OLx
117     pH(i,j,k,bi,bj) = 8. _d 0
118     ENDDO
119     ENDDO
120     ENDDO
121     ENDDO
122     ENDDO
123    
124     DO bj=myByLo(myThid),myByHi(myThid)
125     DO bi=myBxLo(myThid),myBxHi(myThid)
126     DO j=1-OLy,sNy+OLy
127     DO i=1-OLx,sNx+OLx
128     ak0(i,j,bi,bj)=0. _d 0
129     ak1(i,j,bi,bj)=0. _d 0
130     ak2(i,j,bi,bj)=0. _d 0
131     akw(i,j,bi,bj)=0. _d 0
132     akb(i,j,bi,bj)=0. _d 0
133     akf(i,j,bi,bj)=0. _d 0
134     ak1p(i,j,bi,bj)=0. _d 0
135     ak2p(i,j,bi,bj)=0. _d 0
136     ak3p(i,j,bi,bj)=0. _d 0
137     aksi(i,j,bi,bj)=0. _d 0
138     fugf(i,j,bi,bj)=0. _d 0
139     ff(i,j,bi,bj)=0. _d 0
140     ft(i,j,bi,bj)=0. _d 0
141     st(i,j,bi,bj)=0. _d 0
142     bt(i,j,bi,bj)=0. _d 0
143     ENDDO
144     ENDDO
145     ENDDO
146     ENDDO
147    
148     pH_isLoaded = .FALSE.
149     IF ( nIter0.GT.PTRACERS_Iter0 .OR.
150     & (nIter0.EQ.PTRACERS_Iter0 .AND. pickupSuff.NE.' ')
151     & ) THEN
152     C Read pH from a pickup file if needed
153     CALL BLING_READ_PICKUP(
154     O pH_isLoaded,
155     I nIter0, myThid )
156     ENDIF
157    
158     DO bj=myByLo(myThid),myByHi(myThid)
159     DO bi=myBxLo(myThid),myBxHi(myThid)
160    
161     C determine inorganic carbon chem coefficients
162    
163     IF ( .NOT.pH_isLoaded ) THEN
164     C set guess of pH for first step here
165    
166     C$TAF LOOP = parallel
167     DO k=1,Nr
168     DO j=jMin,jMax
169     DO i=iMin,iMax
170     alktmp(i,j) = PTRACER(i,j,k,bi,bj,2)
171     & * maskC(i,j,k,bi,bj)
172     phostmp(i,j)= PTRACER(i,j,k,bi,bj,5)
173     & * maskC(i,j,k,bi,bj)
174     C FOR NON-INTERACTIVE Si
175     IF ( k.eq.1 ) THEN
176     sitmp(i,j) = Silica(i,j,bi,bj) * maskC(i,j,k,bi,bj)
177     ELSE
178     sitmp(i,j) = 0.03 * maskC(i,j,k,bi,bj)
179     ENDIF
180     dictmp(i,j) = PTRACER(i,j,k,bi,bj,1)
181     & * maskC(i,j,k,bi,bj)
182     thetatmp(i,j) = theta(i,j,k,bi,bj)
183     salttmp(i,j) = salt(i,j,k,bi,bj)
184     ENDDO
185     ENDDO
186    
187     CALL CARBON_COEFFS_PRESSURE_DEP(
188     I thetatmp,salttmp,
189     I bi,bj,iMin,iMax,jMin,jMax,k,myThid)
190    
191     C====================================================================
192    
193     c first approximation
194    
195     C$TAF LOOP = parallel
196     DO j=jMin,jMax
197     C$TAF LOOP = parallel
198     DO i=iMin,iMax
199     IF ( maskC(i,j,k,bi,bj) .NE. 0. _d 0) THEN
200     C$TAF init dic_surf = static, 10
201     DO it=1,10
202     C$TAF STORE pH(i,j,k,bi,bj), dictmp(i,j) = dic_surf
203     C$TAF STORE alktmp(i,j), phostmp(i,j), sitmp(i,j) = dic_surf
204     CALL CALC_PCO2_APPROX(
205     I thetatmp(i,j),salttmp(i,j),
206     I dictmp(i,j), phostmp(i,j),
207     I sitmp(i,j),alktmp(i,j),
208     I ak1(i,j,bi,bj),ak2(i,j,bi,bj),
209     I ak1p(i,j,bi,bj),ak2p(i,j,bi,bj),ak3p(i,j,bi,bj),
210     I aks(i,j,bi,bj),akb(i,j,bi,bj),akw(i,j,bi,bj),
211     I aksi(i,j,bi,bj),akf(i,j,bi,bj),
212     I ak0(i,j,bi,bj), fugf(i,j,bi,bj),
213     I ff(i,j,bi,bj),
214     I bt(i,j,bi,bj),st(i,j,bi,bj),ft(i,j,bi,bj),
215     U pH(i,j,k,bi,bj),pCO2(i,j,bi,bj),co3dummy,
216     I i,j,k,bi,bj, it , myThid )
217     ENDDO
218     ENDIF
219     ENDDO
220     ENDDO
221     ENDDO
222     ENDIF
223    
224     C end bi,bj loops
225     ENDDO
226     ENDDO
227    
228     #endif /* ALLOW_BLING */
229    
230     RETURN
231     END

  ViewVC Help
Powered by ViewVC 1.1.22