1 |
C $Header: /u/gcmpack/MITgcm_contrib/ocean_inversion_project/code_bombC14/ini_tr1.F,v 1.2 2006/05/25 17:16:56 dimitri Exp $ |
2 |
C $Name: $ |
3 |
|
4 |
#include "CPP_OPTIONS.h" |
5 |
|
6 |
CBOP |
7 |
C !ROUTINE: INI_TR1 |
8 |
C !INTERFACE: |
9 |
SUBROUTINE INI_TR1( myThid ) |
10 |
C !DESCRIPTION: \bv |
11 |
C *==========================================================* |
12 |
C | SUBROUTINE INI_TR1 |
13 |
C | o Set initial tracer 1 distribution. |
14 |
C *==========================================================* |
15 |
C | Passive tracers 1-N can be initialised so as to dye fluid. |
16 |
C | This routine is hardcoded for Nir's bomb C14 experiment. |
17 |
C *==========================================================* |
18 |
C \ev |
19 |
|
20 |
C !USES: |
21 |
IMPLICIT NONE |
22 |
C === Global variables === |
23 |
#include "SIZE.h" |
24 |
#include "EEPARAMS.h" |
25 |
#include "PARAMS.h" |
26 |
#include "GRID.h" |
27 |
#include "DYNVARS.h" |
28 |
#ifdef ALLOW_PASSIVE_TRACER |
29 |
#include "TR1.h" |
30 |
#endif |
31 |
|
32 |
C !INPUT/OUTPUT PARAMETERS: |
33 |
C == Routine arguments == |
34 |
C myThid - Number of this instance of INI_TR1 |
35 |
INTEGER myThid |
36 |
|
37 |
#ifdef ALLOW_PASSIVE_TRACER |
38 |
|
39 |
C !LOCAL VARIABLES: |
40 |
C == Local variables == |
41 |
C iC, jC - Center of domain |
42 |
C iD, jD - Disitance from domain center. |
43 |
C rad - Radius of initial patch |
44 |
C rD - Radial displacement of point I,J |
45 |
C iG, jG - Global coordinate index |
46 |
C bi,bj - Loop counters |
47 |
C I,J,K |
48 |
INTEGER iC, jC, iD, jD |
49 |
INTEGER iG, jG |
50 |
INTEGER bi, bj |
51 |
INTEGER I, J, K, localWarnings |
52 |
_RL rad, rD |
53 |
CHARACTER*(MAX_LEN_MBUF) msgBuf |
54 |
CEOP |
55 |
|
56 |
_BARRIER |
57 |
|
58 |
C-- Initialise tracer to inline. |
59 |
DO bj = myByLo(myThid), myByHi(myThid) |
60 |
DO bi = myBxLo(myThid), myBxHi(myThid) |
61 |
DO K=1,Nr |
62 |
DO J=1,sNy |
63 |
DO I=1,sNx |
64 |
tr1(I,J,K,bi,bj) = 0. _d 0 |
65 |
ENDDO |
66 |
ENDDO |
67 |
ENDDO |
68 |
ENDDO |
69 |
ENDDO |
70 |
C Set initial tendency terms |
71 |
localWarnings=0 |
72 |
DO bj = myByLo(myThid), myByHi(myThid) |
73 |
DO bi = myBxLo(myThid), myBxHi(myThid) |
74 |
DO K=1,Nr |
75 |
DO J=1,sNy |
76 |
DO I=1,sNx |
77 |
gtr1 (I,J,K,bi,bj) = 0. _d 0 |
78 |
gtr1NM1(I,J,K,bi,bj) = 0. _d 0 |
79 |
ENDDO |
80 |
ENDDO |
81 |
ENDDO |
82 |
ENDDO |
83 |
ENDDO |
84 |
C |
85 |
_EXCH_XYZ_R8(tr1 , myThid ) |
86 |
_EXCH_XYZ_R8(gtr1 , myThid ) |
87 |
_EXCH_XYZ_R8(gtr1NM1 , myThid ) |
88 |
|
89 |
C-- Input file names |
90 |
FiceFile = 'FICE' |
91 |
XkwFile = 'XKW' |
92 |
PatmFile = 'PATM' |
93 |
WOASssFile = 'WOA01_S' |
94 |
WOASstFile = 'WOA01_T' |
95 |
SurfDicPreindFile = 'surf_dic_preind' |
96 |
SurfDicAnthroFile = 'surf_dic_anthro' |
97 |
SurfNatD14CFile = 'surface_naturalD14C' |
98 |
|
99 |
C-- Compute numbers needed to determine DC14 and CO2_atm |
100 |
CALL read_dc14 |
101 |
CALL read_CO2 |
102 |
|
103 |
_BEGIN_MASTER(myThid) |
104 |
#define OCMIP_USE_INTERP |
105 |
C-- OCMIP_USE_INTERP option is useful for high-resolution integrations. |
106 |
C For low-resolution, less that 1-deg, it's best to generate files |
107 |
C separately because of sea-ice fraction. |
108 |
C Caution: OCMIP_USE_INTERP as used is not thread-safe. |
109 |
#ifdef OCMIP_USE_INTERP |
110 |
CALL EXF_INTERP( SurfDicPreindFile,readBinaryPrec, |
111 |
& SurfDicPreind,1, |
112 |
& xC,yC,lon0,lon_inc,lat0,lat_inc,nlon,nlat,1,mythid ) |
113 |
CALL EXF_INTERP( SurfDicAnthroFile,readBinaryPrec, |
114 |
& SurfDicAnthro,1, |
115 |
& xC,yC,lon0,lon_inc,lat0,lat_inc,nlon,nlat,1,mythid ) |
116 |
CALL EXF_INTERP( SurfNatD14CFile,readBinaryPrec, |
117 |
& SurfNatD14C,1, |
118 |
& xC,yC,lon0,lon_inc,lat0,lat_inc,nlon,nlat,1,mythid ) |
119 |
|
120 |
#else |
121 |
CALL MDSREADFIELD ( SurfDicPreindFile, readBinaryPrec, |
122 |
& 'RS', 1, SurfDicPreind, 1,, myThid ) |
123 |
CALL MDSREADFIELD ( SurfDicAnthroFile, readBinaryPrec, |
124 |
& 'RS', 1, SurfDicAnthro, 1,, myThid ) |
125 |
CALL MDSREADFIELD ( SurfNatD14CFile, readBinaryPrec, |
126 |
& 'RS', 1, SurfNatD14C, 1,, myThid ) |
127 |
#endif |
128 |
_END_MASTER(myThid) |
129 |
_EXCH_XY_R4(SurfDicPreind , myThid ) |
130 |
_EXCH_XY_R4(SurfDicAnthro , myThid ) |
131 |
_EXCH_XY_R4(SurfNatD14C , myThid ) |
132 |
|
133 |
#endif /* ALLOW_PASSIVE_TRACER */ |
134 |
|
135 |
RETURN |
136 |
END |
137 |
|
138 |
C==================================================================== |
139 |
|
140 |
SUBROUTINE read_dc14 |
141 |
|
142 |
C -------------------------------------------------------------------- |
143 |
C Reads temporal history of atmospheric DC14 in permil |
144 |
C |
145 |
C DC14_atm_year = year of measurement (mid-year) |
146 |
C DC14_atm_n = DC14 in northern hemisphere |
147 |
C DC14_atm_eq = DC14 at equator |
148 |
C DC14_atm_s = DC14 in southern hemisphere |
149 |
C -------------------------------------------------------------------- |
150 |
|
151 |
IMPLICIT NONE |
152 |
#include "SIZE.h" |
153 |
#include "EEPARAMS.h" |
154 |
#include "TR1.h" |
155 |
|
156 |
INTEGER n, iu |
157 |
PARAMETER (iu=10) |
158 |
|
159 |
CHARACTER*80 filen |
160 |
C |
161 |
C OPEN FILE |
162 |
C --------- |
163 |
filen='DC14_atm' |
164 |
OPEN(UNIT=iu,FILE=filen,STATUS='old') |
165 |
C |
166 |
C READ FILE |
167 |
C --------- |
168 |
DO n = 1, 240 |
169 |
READ(iu,*) |
170 |
+ DC14_atm_year(n), DC14_atm_n(n), |
171 |
+ DC14_atm_eq (n), DC14_atm_s(n) |
172 |
C WRITE(6,100) |
173 |
C + DC14_atm_year(n), DC14_atm_n(n), |
174 |
C + DC14_atm_eq (n), DC14_atm_s(n) |
175 |
END DO |
176 |
|
177 |
100 FORMAT(f7.2, 4f8.2) |
178 |
|
179 |
RETURN |
180 |
END |
181 |
|
182 |
C==================================================================== |
183 |
|
184 |
SUBROUTINE read_CO2 |
185 |
|
186 |
C -------------------------------------------------------------------- |
187 |
C Reads temporal history of atmospheric CO2 in ppm |
188 |
C |
189 |
C CO2_atm_year = year of measurement |
190 |
C CO2_atm = atmospheric CO2 |
191 |
C -------------------------------------------------------------------- |
192 |
|
193 |
IMPLICIT NONE |
194 |
#include "SIZE.h" |
195 |
#include "EEPARAMS.h" |
196 |
#include "TR1.h" |
197 |
|
198 |
INTEGER n, iu |
199 |
PARAMETER (iu=10) |
200 |
|
201 |
CHARACTER*80 filen |
202 |
C |
203 |
C OPEN FILE |
204 |
C --------- |
205 |
filen='CO2_atm' |
206 |
OPEN(UNIT=iu,FILE=filen,STATUS='old') |
207 |
C |
208 |
C READ FILE |
209 |
C --------- |
210 |
DO n = 1, 478 |
211 |
READ(iu,*) CO2_atm_year(n), CO2_atm(n) |
212 |
C WRITE(6,100) CO2_atm_year(n), CO2_atm(n) |
213 |
END DO |
214 |
|
215 |
100 FORMAT(f7.2, f8.2) |
216 |
|
217 |
RETURN |
218 |
END |