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

Annotation 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 - (hide 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 jahn 1.8 C $Header: /u/gcmpack/MITgcm_contrib/darwin2/pkg/darwin/darwin_readparams.F,v 1.7 2013/04/16 20:21:21 jahn Exp $
2 jahn 1.2 C $Name: $
3 jahn 1.1
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 jahn 1.2 & darwin_PARunits, darwin_W_to_uEins,
68 jahn 1.3 & darwin_ironUnits,
69     & darwin_dustFrac,
70 jahn 1.1 #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 jahn 1.4 & ,darwin_wavebands
80 jahn 1.1 #ifdef DAR_CALC_ACDOM
81 jahn 1.8 & ,darwin_Sdom, darwin_lambda_aCDOM, darwin_aCDOM_fac
82 jahn 1.1 #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 jahn 1.2 darwin_PARunits='Ein/m2/d '
139     darwin_W_to_uEins=1. _d 0/0.2174 _d 0
140 jahn 1.3 darwin_ironUnits='mol Fe/m2/s '
141     darwin_dustFrac=0.035 _d 0
142 jahn 1.1 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 jahn 1.4 DO ilam=1,tlam+1
170 jahn 1.7 darwin_wavebands(ilam) = -1
171 jahn 1.4 ENDDO
172 jahn 1.1 #ifdef DAR_CALC_ACDOM
173     darwin_Sdom = 0.014 _d 0
174 jahn 1.8 darwin_lambda_aCDOM = 450 _d 0
175     darwin_aCDOM_fac = 0.2 _d 0
176 jahn 1.1 #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 jahn 1.5 darwin_radtrans_niter = -2
195 jahn 1.1 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 jahn 1.2 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 jahn 1.1
226 jahn 1.3 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 jahn 1.1 #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 jahn 1.4 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 jahn 1.1 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 jahn 1.4 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 jahn 1.1 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 jahn 1.6 CALL WRITE_1D_RL( darwin_waves, tlam, 0,
374 jahn 1.1 & '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 jahn 1.5 &' /* iterations/method for solveing 3-stream equations */')
422 jahn 1.1 #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