/[MITgcm]/MITgcm_contrib/darwin2/pkg/darwin/darwin_readparams.F
ViewVC logotype

Contents of /MITgcm_contrib/darwin2/pkg/darwin/darwin_readparams.F

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


Revision 1.8 - (show annotations) (download)
Tue Apr 16 20:21:57 2013 UTC (12 years, 3 months ago) by jahn
Branch: MAIN
CVS Tags: ctrb_darwin2_ckpt64k_20130723, ctrb_darwin2_ckpt64h_20130528, ctrb_darwin2_ckpt64m_20130820, ctrb_darwin2_ckpt64n_20130826, ctrb_darwin2_ckpt64o_20131024, ctrb_darwin2_ckpt64i_20130622, ctrb_darwin2_ckpt64g_20130503, ctrb_darwin2_ckpt64l_20130806, ctrb_darwin2_ckpt64j_20130704, ctrb_darwin2_ckpt64p_20131118, ctrb_darwin2_ckpt64q_20131118, ctrb_darwin2_ckpt64p_20131024
Changes since 1.7: +4 -2 lines
make reference waveband for aCDOM a runtime parameter

1 C $Header: /u/gcmpack/MITgcm_contrib/darwin2/pkg/darwin/darwin_readparams.F,v 1.7 2013/04/16 20:21:21 jahn Exp $
2 C $Name: $
3
4 #include "DARWIN_OPTIONS.h"
5
6 CBOP
7 C !ROUTINE: DARWIN_READPARMS
8
9 C !INTERFACE: ==========================================================
10 SUBROUTINE DARWIN_READPARMS( myThid )
11
12 C !DESCRIPTION:
13 C Initialize DARWIN parameters, read in data.darwin
14
15 C !USES: ===============================================================
16 IMPLICIT NONE
17 #include "SIZE.h"
18 #include "EEPARAMS.h"
19 #include "PARAMS.h"
20 #include "GCHEM.h"
21 #include "DARWIN_SIZE.h"
22 #include "DARWIN_IO.h"
23 #include "DARWIN_PARAMS.h"
24 #ifdef WAVEBANDS
25 #include "SPECTRAL_SIZE.h"
26 #include "SPECTRAL.h"
27 #include "WAVEBANDS_PARAMS.h"
28 #endif
29
30 C !INPUT PARAMETERS: ===================================================
31 C myThid :: thread number
32 INTEGER myThid
33
34 C !OUTPUT PARAMETERS: ==================================================
35 C none
36
37 #ifdef ALLOW_DARWIN
38
39 C !LOCAL VARIABLES: ====================================================
40 C iUnit :: unit number for I/O
41 C msgBuf :: message buffer
42 INTEGER iUnit,errIO
43 CHARACTER*(MAX_LEN_MBUF) msgBuf
44 #if defined(WAVEBANDS) || defined(OASIM)
45 INTEGER ilam,i
46 _RL planck, c, hc, oavo, hcoavo, rlamm
47 #endif
48 CEOP
49
50 NAMELIST /DARWIN_FORCING/
51 & darwin_iceFile,
52 & darwin_ironFile,
53 & darwin_PARFile,
54 & darwin_nutWVelFile,
55 & darwin_PO4_relaxFile, darwin_NO3_relaxFile,
56 & darwin_FeT_relaxFile, darwin_Si_relaxFile,
57 & darwin_PO4_fluxFile, darwin_NO3_fluxFile,
58 & darwin_FeT_fluxFile, darwin_Si_fluxFile,
59 & darwin_waterabsorbFile, darwin_phytoabsorbFile,
60 & darwin_surfacespecFile, darwin_acdomFile,
61 & darwin_particleabsorbFile,
62 #ifdef OASIM
63 & darwin_oasim_edFile, darwin_oasim_esFile,
64 #endif
65 & darwin_relaxscale,
66 & darwin_ForcingPeriod, darwin_ForcingCycle,
67 & darwin_PARunits, darwin_W_to_uEins,
68 & darwin_ironUnits,
69 & darwin_dustFrac,
70 #ifdef ALLOW_PAR_DAY
71 & darwin_PARavPeriod,
72 #endif
73 & darwin_seed
74
75 #if defined(WAVEBANDS) || defined(OASIM)
76 NAMELIST /DARWIN_SPECTRAL_PARM/
77 & darwin_waves
78 #ifdef WAVEBANDS
79 & ,darwin_wavebands
80 #ifdef DAR_CALC_ACDOM
81 & ,darwin_Sdom, darwin_lambda_aCDOM, darwin_aCDOM_fac
82 #endif
83 #ifdef DAR_DIAG_ACDOM
84 & ,darwin_diag_acdom_ilam
85 #endif
86 #ifdef DAR_RADTRANS
87 & ,darwin_PAR_ilamLo
88 & ,darwin_PAR_ilamHi
89 & ,darwin_radmodThresh
90 & ,darwin_Dmax
91 & ,darwin_rmus
92 & ,darwin_rmuu
93 & ,darwin_bbw
94 & ,darwin_bbphy
95 & ,darwin_bbmin
96 & ,darwin_radtrans_kmax
97 & ,darwin_radtrans_niter
98 & ,darwin_part_size_P
99 #endif
100 #endif /* WAVEBANDS */
101 #endif /* WAVEBANDS || OASIM */
102
103 #ifdef DAR_DIAG_CHL
104 NAMELIST /DARWIN_CHL/
105 & Geider_Bigalphachl, Geider_smallalphachl,
106 & Geider_Bigchl2cmax, Geider_smallchl2cmax,
107 & Geider_Bigchl2cmin,
108 & Doney_Bmin, Doney_Bmax, Doney_PARstar,
109 & Cloern_A, Cloern_B, Cloern_C, Cloern_chl2cmin
110 #endif
111
112 #ifdef ALLOW_CARBON
113 NAMELIST /DIC_FORCING/
114 & DIC_windFile, DIC_atmospFile,
115 & dic_pCO2, dic_int1, dic_int2, dic_int3, dic_int4
116 #endif
117
118 C Set defaults values for parameters in DARWIN_IO.h
119 darwin_iceFile=' '
120 darwin_ironFile=' '
121 darwin_PARFile=' '
122 darwin_nutWVelFile=' '
123 darwin_PO4_relaxFile=' '
124 darwin_NO3_relaxFile=' '
125 darwin_FeT_relaxFile=' '
126 darwin_Si_relaxFile=' '
127 darwin_PO4_fluxFile=' '
128 darwin_NO3_fluxFile=' '
129 darwin_FeT_fluxFile=' '
130 darwin_Si_fluxFile=' '
131 darwin_waterabsorbFile=' '
132 darwin_phytoabsorbFile=' '
133 darwin_particleabsorbFile=' '
134 darwin_surfacespecFile=' '
135 darwin_acdomFile=' '
136 darwin_oasim_edFile=' '
137 darwin_oasim_esFile=' '
138 darwin_PARunits='Ein/m2/d '
139 darwin_W_to_uEins=1. _d 0/0.2174 _d 0
140 darwin_ironUnits='mol Fe/m2/s '
141 darwin_dustFrac=0.035 _d 0
142 darwin_relaxscale=0. _d 0
143 darwin_seed=0
144 c default periodic forcing to same as for GCHEM
145 darwin_forcingPeriod=gchem_ForcingPeriod
146 darwin_forcingCycle=gchem_ForcingCycle
147
148 #ifdef ALLOW_CARBON
149 DIC_windFile = ' '
150 DIC_atmospFile= ' '
151 dic_int1 = 0
152 dic_int2 = 0
153 dic_int3 = 0
154 dic_int4 = 0
155 dic_pCO2 = 0. _d 0
156 #endif
157
158 #ifdef ALLOW_PAR_DAY
159 darwin_PARavPeriod=86400. _d 0
160 #endif
161
162 #if defined(WAVEBANDS) || defined(OASIM)
163 DO ilam=1,tlam
164 darwin_waves(ilam) = 0
165 ENDDO
166 #endif
167
168 #ifdef WAVEBANDS
169 DO ilam=1,tlam+1
170 darwin_wavebands(ilam) = -1
171 ENDDO
172 #ifdef DAR_CALC_ACDOM
173 darwin_Sdom = 0.014 _d 0
174 darwin_lambda_aCDOM = 450 _d 0
175 darwin_aCDOM_fac = 0.2 _d 0
176 #endif
177 #ifdef DAR_DIAG_ACDOM
178 c value >= 100 will be converted to index in wavebands_init_fixed
179 darwin_diag_acdom_ilam = 450
180 #endif
181 #ifdef DAR_RADTRANS
182 darwin_PAR_ilamLo = 1
183 darwin_PAR_ilamHi = tlam
184 darwin_radmodThresh = 1 _d -4
185 darwin_Dmax = 500 _d 0
186 darwin_rmus = 1.0/0.83 _d 0
187 darwin_rmuu = 1.0/0.4 _d 0
188 darwin_bbmin = 0.0002 _d 0
189 darwin_bbw = 0.5 _d 0
190 do i=1,tnabp
191 darwin_bbphy(i) = 0 _d 0
192 enddo
193 darwin_radtrans_kmax = Nr
194 darwin_radtrans_niter = -2
195 darwin_part_size_P = 1 _d -15 ! mmol P per particle
196 #endif
197 #endif /* WAVEBANDS */
198
199 C Open and read the data.darwin file
200 _BEGIN_MASTER(myThid)
201 WRITE(msgBuf,'(A)') ' DARWIN_READPARMS: opening data.darwin'
202 CALL PRINT_MESSAGE(msgBuf, standardMessageUnit,
203 & SQUEEZE_RIGHT , 1)
204 CALL OPEN_COPY_DATA_FILE(
205 I 'data.darwin', 'DARWIN_READPARAMS',
206 O iUnit,
207 I myThid )
208 READ(UNIT=iUnit,NML=DARWIN_FORCING)
209 #ifdef ALLOW_PAR_DAY
210 darwin_PARnav = NINT(darwin_PARavPeriod*nsubtime/dTtracerLev(1))
211 #endif
212 C factor for conversion to uEin/m2/s
213 IF ( darwin_PARunits(1:16) .EQ. 'uEin/m2/s ' ) THEN
214 darwin_PARFileConv = 1. _d 0
215 ELSEIF ( darwin_PARunits(1:16) .EQ. 'Ein/m2/d ' ) THEN
216 darwin_PARFileConv = 1. _d 6/86400. _d 0
217 ELSEIF ( darwin_PARunits(1:16) .EQ. 'W/m2 ' ) THEN
218 darwin_PARFileConv = darwin_W_to_uEins
219 ELSE
220 WRITE(msgBuf,'(2A)') 'S/R DARWIN_READPARMS:',
221 & 'darwin_PARunits must be one of Ein/m2/d, uEin/m2/s, W/m2'
222 CALL PRINT_ERROR( msgBuf , 1)
223 STOP 'unknown darwin_PARunits'
224 ENDIF
225
226 C factor for conversion to mmol Fe/m2/s (bioavailable fraction is done later)
227 IF ( darwin_ironUnits(1:16) .EQ. 'mol Fe/m2/s ' ) THEN
228 darwin_ironFileConv = 1000. _d 0
229 ELSEIF ( darwin_ironUnits(1:16) .EQ. 'kg dust/m2/s ' ) THEN
230 darwin_ironFileConv = 1000.*darwin_dustFrac*1000./55.845 _d 0
231 ELSE
232 WRITE(msgBuf,'(2A)') 'S/R DARWIN_READPARMS:',
233 & 'darwin_ironUnits must be one of mol Fe/m2/s, kg dust/m2/s'
234 CALL PRINT_ERROR( msgBuf , 1)
235 STOP 'unknown darwin_Feunits'
236 ENDIF
237
238 #ifdef DAR_DIAG_CHL
239 C default values
240 C Geider: chl:c = max(chl2cmin, chl2cmax/(1+(chl2cmax*alphachl*PARday)/(2*Pcm)))
241 C Pcm = mu*limit*phytoTempFunction
242 Geider_smallalphachl = 2. _d -6 ! mmol C (mg Chl)-1 m2 (uEin)-1
243 Geider_Bigalphachl = 1. _d -6 ! mmol C (mg Chl)-1 m2 (uEin)-1
244 Geider_smallchl2cmax = 0.35 _d 0 ! mg Chl (mmol C)-1
245 Geider_Bigchl2cmax = 0.65 _d 0 ! mg Chl (mmol C)-1
246 Geider_smallchl2cmin = 0.003 _d 0 * 12. _d 0 ! mg Chl a/mmol C
247 Geider_Bigchl2cmin = 0.003 _d 0 * 12. _d 0 ! mg Chl a/mmol C
248
249 C Doney: chl:c = (Bmax - (Bmax-Bmin)*MIN(1,PARday/PARstar))*limit
250 Doney_Bmax = 12. _d 0 / 37. _d 0 ! mg Chl a/mmol C
251 Doney_Bmin = 12. _d 0 / 90. _d 0 ! mg Chl a/mmol C
252 Doney_PARstar = 90. _d 0 / 0.2174 _d 0 ! uEin/m2/s
253
254 C Cloern: chl:c = chl2cmin + A*exp(B*T)*exp(-C*PARday)*limit
255 Cloern_chl2cmin = 0.003 _d 0 * 12. _d 0 ! mg Chl a/mmol C
256 Cloern_A = 0.0154 _d 0 * 12. _d 0 ! mg Chl a/mmol C
257 Cloern_B = 0.050 _d 0 ! (degree C)^{-1}
258 Cloern_C = 0.059 _d 0 * 86400. _d 0 / 1. _d 6 ! m^2 s/uEin
259
260 READ(UNIT=iUnit,NML=DARWIN_CHL)
261 #endif /* DAR_DIAG_CHL */
262
263 #ifdef ALLOW_CARBON
264 READ(UNIT=iUnit,NML=DIC_FORCING)
265 #endif
266
267 #ifdef DAR_RADTRANS
268 #ifndef DAR_NONSPECTRAL_BACKSCATTERING_RATIO
269 DO i=1,tnabp
270 IF ( darwin_bbphy(i) .NE. 0 _d 0 ) THEN
271 WRITE(msgBuf,'(2A)') 'S/R DARWIN_READPARMS:',
272 & 'darwin_bbphy is obsolete.'
273 CALL PRINT_ERROR( msgBuf , 1)
274 WRITE(msgBuf,'(2A)') 'S/R DARWIN_READPARMS:',
275 & 'Backscattering coefficients are now read from'
276 CALL PRINT_ERROR( msgBuf , 1)
277 WRITE(msgBuf,'(2A)') 'S/R DARWIN_READPARMS:',
278 & 'darwin_phytoabsorbFile.'
279 CALL PRINT_ERROR( msgBuf , 1)
280 ENDIF
281 ENDDO
282 #endif
283 #endif
284
285 #if defined(WAVEBANDS) || defined(OASIM)
286 READ(UNIT=iUnit,NML=DARWIN_SPECTRAL_PARM,IOSTAT=errIO)
287 IF ( errIO .LT. 0 ) THEN
288 WRITE(msgBuf,'(A)')
289 & 'S/R DARWIN_READPARMS'
290 CALL PRINT_ERROR( msgBuf , 1)
291 WRITE(msgBuf,'(A)')
292 & 'Error reading darwin package'
293 CALL PRINT_ERROR( msgBuf , 1)
294 WRITE(msgBuf,'(A)')
295 & 'parameter file "data.darwin"'
296 CALL PRINT_ERROR( msgBuf , 1)
297 WRITE(msgBuf,'(A)')
298 & 'Problem in namelist DARWIN_SPECTRAL_PARM'
299 CALL PRINT_ERROR( msgBuf , 1)
300 STOP 'ABNORMAL END: S/R DARWIN_READPARMS'
301 ENDIF
302 #endif
303
304 WRITE(msgBuf,'(A)')
305 & ' DARWIN_READPARMS: finished reading data.darwin'
306 CALL PRINT_MESSAGE(msgBuf, standardMessageUnit,
307 & SQUEEZE_RIGHT , 1)
308
309 C Close the open data file
310 CLOSE(iUnit)
311
312 #if defined(WAVEBANDS) || defined(OASIM)
313 IF ( darwin_waves(1).EQ.0.0 .AND. tlam.EQ.13 ) THEN
314 darwin_waves(1) = 400
315 darwin_waves(2) = 425
316 darwin_waves(3) = 450
317 darwin_waves(4) = 475
318 darwin_waves(5) = 500
319 darwin_waves(6) = 525
320 darwin_waves(7) = 550
321 darwin_waves(8) = 575
322 darwin_waves(9) = 600
323 darwin_waves(10) = 625
324 darwin_waves(11) = 650
325 darwin_waves(12) = 675
326 darwin_waves(13) = 700
327 IF ( darwin_wavebands(1).NE.0.0 ) THEN
328 WRITE(msgBuf,'(2A)') 'DARWIN_READPARMS: ',
329 & 'darwin_wavebands given without darwin_waves.'
330 CALL PRINT_ERROR( msgBuf, myThid )
331 WRITE(msgBuf,'(2A)') 'DARWIN_READPARMS: ',
332 & 'using standard waveband centers 400, 425, ...'
333 CALL PRINT_ERROR( msgBuf, myThid )
334 WRITE(msgBuf,'(2A)') 'DARWIN_READPARMS: ',
335 & 'provide darwin_waves if this is not desired.'
336 CALL PRINT_ERROR( msgBuf, myThid )
337 ENDIF
338 ENDIF
339 c Quanta conversion
340 planck = 6.6256 _d -34 !Plancks constant J sec
341 c = 2.998 _d 8 !speed of light m/sec
342 hc = 1.0/(planck*c)
343 oavo = 1.0/6.023 _d 23 ! 1/Avogadros number
344 hcoavo = hc*oavo
345 do ilam = 1,tlam
346 rlamm = darwin_waves(ilam)*1 _d -9 !lambda in m
347 WtouEins(ilam) = 1 _d 6*rlamm*hcoavo !Watts to uEin/s conversion
348 IF ( rlamm .EQ. 0.0 ) THEN
349 WRITE(msgBuf,'(2A)') 'DARWIN_READPARMS: ',
350 & 'please provide wavelengths in darwin_waves.'
351 CALL PRINT_ERROR( msgBuf, myThid )
352 STOP 'ABNORMAL END: S/R DARWIN_READPARMS'
353 ENDIF
354 enddo
355 #endif
356
357 C-- Print a summary of parameter values:
358 iUnit = standardMessageUnit
359 WRITE(msgBuf,'(A)') '// ==================================='
360 CALL PRINT_MESSAGE( msgBuf, iUnit, SQUEEZE_RIGHT , myThid )
361 WRITE(msgBuf,'(A)') '// darwin parameters '
362 CALL PRINT_MESSAGE( msgBuf, iUnit, SQUEEZE_RIGHT , myThid )
363 WRITE(msgBuf,'(A)') '// ==================================='
364 CALL PRINT_MESSAGE( msgBuf, iUnit, SQUEEZE_RIGHT , myThid )
365 CALL WRITE_0D_I( darwin_seed, INDEX_NONE,
366 & 'darwin_seed =',
367 & ' /* seed for random number generator */')
368
369 WRITE(msgBuf,'(A)') ' -----------------------------------'
370 CALL PRINT_MESSAGE( msgBuf, iUnit, SQUEEZE_RIGHT, myThid )
371
372 #if defined(WAVEBANDS) || defined(OASIM)
373 CALL WRITE_1D_RL( darwin_waves, tlam, 0,
374 & 'darwin_waves =',
375 &' /* "central" wavelengths of wavebands */')
376 #endif
377 #ifdef WAVEBANDS
378 #ifdef DAR_CALC_ACDOM
379 CALL WRITE_0D_RL( darwin_Sdom, INDEX_NONE,
380 & 'darwin_Sdom =',
381 &' /* spectral slope for aCDOM */')
382 #endif
383 #ifdef DAR_DIAG_ACDOM
384 CALL WRITE_0D_I( darwin_diag_acdom_ilam, INDEX_NONE,
385 & 'darwin_diag_acdom_ilam =',
386 &' /* waveband to use for aCDOM diagnostic */')
387 #endif
388 #ifdef DAR_RADTRANS
389 CALL WRITE_0D_I( darwin_PAR_ilamLo, INDEX_NONE,
390 & 'darwin_PAR_ilamLo =',
391 &' /* waveband index of PAR lower bound */')
392 CALL WRITE_0D_I( darwin_PAR_ilamHi, INDEX_NONE,
393 & 'darwin_PAR_ilamHi =',
394 &' /* waveband index of PAR upper bound */')
395 CALL WRITE_0D_RL( darwin_radmodThresh, INDEX_NONE,
396 & 'darwin_radmodThresh =',
397 &' /* threshold for calling radmod (W/m2/waveband) */')
398 CALL WRITE_0D_RL( darwin_Dmax, INDEX_NONE,
399 & 'darwin_Dmax =',
400 &' /* depth at which Ed is assumed zero */')
401 CALL WRITE_0D_RL( darwin_rmus, INDEX_NONE,
402 & 'darwin_rmus =',
403 &' /* inverse average cosine of downward diffuse irradiance */')
404 CALL WRITE_0D_RL( darwin_rmuu, INDEX_NONE,
405 & 'darwin_rmuu =',
406 &' /* inverse average cosine of upward diffuse irradiance */')
407 CALL WRITE_0D_RL( darwin_bbw, INDEX_NONE,
408 & 'darwin_bbw =',
409 &' /* backscattering to forward scattering ratio for water */')
410 CALL WRITE_0D_RL( darwin_bbmin, INDEX_NONE,
411 & 'darwin_bbmin =',
412 &' /* minimum backscattering coefficient (1/m) */')
413 CALL WRITE_1D_RL( darwin_bbphy, tnabp, 0,
414 & 'darwin_bbphy =',
415 &' /* backscattering to forward scattering ratio for phyto */')
416 CALL WRITE_0D_I( darwin_radtrans_kmax, INDEX_NONE,
417 & 'darwin_radtrans_kmax =',
418 &' /* deepest level in which to compute irradiances */')
419 CALL WRITE_0D_I( darwin_radtrans_niter, INDEX_NONE,
420 & 'darwin_radtrans_niter =',
421 &' /* iterations/method for solveing 3-stream equations */')
422 #endif /* DAR_RADTRANS */
423 #endif /* WAVEBANDS */
424
425 WRITE(msgBuf,'(A)') ' ==================================='
426 CALL PRINT_MESSAGE( msgBuf, iUnit, SQUEEZE_RIGHT, myThid )
427
428 _END_MASTER(myThid)
429
430 C Everyone else must wait for the parameters to be loaded
431 _BARRIER
432
433 #endif /* ALLOW_DARWIN */
434
435 RETURN
436 END

  ViewVC Help
Powered by ViewVC 1.1.22