Line data Source code
1 : !------------------------------------------------------------------------------
2 : ! Harmonized Emissions Component (HEMCO) !
3 : !------------------------------------------------------------------------------
4 : !BOP
5 : !
6 : ! !MODULE: hcox_paranox_mod.F90
7 : !
8 : ! !DESCRIPTION: Module HCOX\_PARANOX\_MOD contains routines to
9 : ! compute ship emissions and associated concentrations of NO, NO2, HNO3
10 : ! and O3 from NO ship emission data. This follows the implementation of
11 : ! the PARANOX ship plume model in GEOS-Chem.
12 : !\\
13 : !\\
14 : ! This module calculates production rates of NO, NO2, HNO3, and O3, as well
15 : ! as loss rates of O3 and HNO3. All fluxes are in kg species/m2/s. The
16 : ! O3 and HNO3 loss fluxes are not converted to a deposition velocity, but
17 : ! rather saved out as mass fluxes (kg/m2/s) into diagnostics
18 : ! 'PARANOX\_O3\_DEPOSITION\_FLUX' and 'PARANOX\_HNO3\_DEPOSITION\_FLUX',
19 : ! respectively. In order to use them, they must be imported explicitly via
20 : ! routine Diagn\_Get (from module hco\_diagn\_mod.F90). This approach avoids
21 : ! problems with uncrealistically high loss rates for loss ambient air
22 : ! concentrations of O3 or HNO3.
23 : !\\
24 : !\\
25 : ! The PARANOx look-up-table can be provided in netCDF or ASCII (txt) format.
26 : ! The latter is particularly useful for running PARANOx in an ESMF environment,
27 : ! where 7-dimensional netCDF files are currently not supported. The input data
28 : ! format can be specified in the HEMCO configuration file (in the PARANOx
29 : ! extensions section).
30 : ! The txt-files can be generated from the previously read netCDF data using
31 : ! subroutine WRITE\_LUT\_TXTFILE.
32 : !\\
33 : !\\
34 : ! References:
35 : ! \begin{itemize}
36 : ! \item Vinken, G. C. M., Boersma, K. F., Jacob, D. J., and Meijer, E. W.:
37 : ! Accounting for non-linear chemistry of ship plumes in the
38 : ! GEOS-Chem global chemistry transport model, Atmos. Chem. Phys., 11,
39 : ! 11707-11722, doi:10.5194/acp-11-11707-2011, 2011.
40 : ! \end{itemize}
41 : !
42 : ! The initial look up tables (LUT) distributed with GEOS-Chem v9-01-03
43 : ! used 7 input variables: Temperature, J(NO2), J(O1D), solar elevation angles
44 : ! at emission time and 5 hours later, and ambient concentrations of NOx
45 : ! and O3. This version was documented by Vinken et al. (2011). Subsequently,
46 : ! we added wind speed as an input variable. We also use J(OH) rather than J(O1D)
47 : ! to index the LUT (C. Holmes, 6 May 2013)
48 : !
49 : ! The LUTs contain 3 quantities:
50 : ! FracNOx : The fraction of NOx emitted from ships that remains as NOx
51 : ! after 5 hours of plume aging. mol/mol
52 : ! OPE : Ozone production efficiency, mol(O3)/mol(HNO3)
53 : ! The net production of O3 per mole of ship NOx oxidized over
54 : ! 5 hours of plume aging. Can be negative!
55 : ! Defined as OPE = [ P(O3) - L(O3) ] / P(HNO3), where each P
56 : ! and L term is an integral over 5 hours. Net O3 production
57 : ! in the plume is E(NOx) * (1-FracNOx) * OPE, where E(NOx) is
58 : ! the emission rate of NOx from the ship (e.g. units: mol/s).
59 : ! MOE : Methane oxidation efficiency, mol(CH4)/mol(NOx)
60 : ! The net oxidation of CH4 per mole of NOx emitted from ships
61 : ! over 5 hours of plume aging.
62 : ! Defined as MOE = L(CH4) / E(NOx).
63 : !\\
64 : !\\
65 : ! The solar elevation angles 5 hours ago are calculated using HEMCO subroutine
66 : ! HCO\_GetSUNCOS. This is the same routine that is used to calculate the solar
67 : ! zenith angles for the current time.
68 : !\\
69 : !\\
70 : ! !INTERFACE:
71 : !
72 : MODULE HCOX_ParaNOx_MOD
73 : !
74 : ! !USES:
75 : !
76 : USE HCO_Error_MOD
77 : USE HCO_Diagn_MOD
78 : USE HCO_State_MOD, ONLY : HCO_State
79 : USE HCOX_State_MOD, ONLY : Ext_State
80 :
81 : IMPLICIT NONE
82 : PRIVATE
83 : !
84 : ! !PUBLIC MEMBER FUNCTIONS:
85 : !
86 : PUBLIC :: HCOX_ParaNOx_Run
87 : PUBLIC :: HCOX_ParaNOx_Init
88 : PUBLIC :: HCOX_ParaNOx_Final
89 : !
90 : ! !PRIVATE MEMBER FUNCTIONS:
91 : !
92 : !
93 : ! !REMARKS:
94 : ! Adapted from the code in GeosCore/paranox_mod.F prior to GEOS-Chem v10-01.
95 : !
96 : ! !REVISION HISTORY:
97 : ! 06 Aug 2013 - C. Keller - Initial version
98 : ! See https://github.com/geoschem/hemco for complete history
99 : !EOP
100 : !------------------------------------------------------------------------------
101 : !BOC
102 : !
103 : ! !MODULE VARIABLES:
104 :
105 : ! Number of values for each variable in the look-up table
106 : INTEGER, PARAMETER :: nT=4, nJ=4, nO3=4, nNOx=5, nSEA=12, nWS=5
107 :
108 : ! Now place all module variables in a lderived type object (for a linked
109 : ! list) so that we can have one instance per node in an MPI environment.
110 : TYPE :: MyInst
111 :
112 : ! Scalars
113 : INTEGER :: Instance
114 : INTEGER :: ExtNr
115 : INTEGER :: IDTNO
116 : INTEGER :: IDTNO2
117 : INTEGER :: IDTHNO3
118 : INTEGER :: IDTO3
119 : REAL*8 :: MW_O3
120 :
121 : REAL*8 :: MW_NO
122 : REAL*8 :: MW_NO2
123 : REAL*8 :: MW_HNO3
124 : REAL*8 :: MW_AIR
125 :
126 : ! Arrays
127 : REAL(hp), POINTER :: ShipNO(:,:,:)
128 :
129 : ! For SunCosMid 5hrs ago
130 : REAL(hp), POINTER :: SC5(:,:)
131 :
132 : ! Deposition fluxes in kg/m2/s
133 : REAL(sp), POINTER :: DEPO3 (:,:)
134 : REAL(sp), POINTER :: DEPHNO3(:,:)
135 :
136 : ! Reference values of variables in the look-up tables
137 : REAL*4 :: Tlev(nT)
138 : REAL*4 :: JNO2lev(nJ)
139 : REAL*4 :: O3lev(nO3)
140 : REAL*4 :: SEA0lev(nSEA)
141 : REAL*4 :: SEA5lev(nSEA)
142 : REAL*4 :: JRATIOlev(nJ)
143 : REAL*4 :: NOXlev(nNOx)
144 : REAL*4 :: WSlev(nWS)
145 :
146 : ! Look-up tables currently used in GEOS-Chem (likely in v10-01)
147 : ! Described by Holmes et al. (2014), now includes effects of wind speed
148 : ! Last two digits in LUT names indicate wind speed in m/s
149 : REAL(sp), POINTER :: FRACNOX_LUT02(:,:,:,:,:,:,:)
150 : REAL(sp), POINTER :: FRACNOX_LUT06(:,:,:,:,:,:,:)
151 : REAL(sp), POINTER :: FRACNOX_LUT10(:,:,:,:,:,:,:)
152 : REAL(sp), POINTER :: FRACNOX_LUT14(:,:,:,:,:,:,:)
153 : REAL(sp), POINTER :: FRACNOX_LUT18(:,:,:,:,:,:,:)
154 : REAL(sp), POINTER :: OPE_LUT02 (:,:,:,:,:,:,:)
155 : REAL(sp), POINTER :: OPE_LUT06 (:,:,:,:,:,:,:)
156 : REAL(sp), POINTER :: OPE_LUT10 (:,:,:,:,:,:,:)
157 : REAL(sp), POINTER :: OPE_LUT14 (:,:,:,:,:,:,:)
158 : REAL(sp), POINTER :: OPE_LUT18 (:,:,:,:,:,:,:)
159 : REAL(sp), POINTER :: MOE_LUT02 (:,:,:,:,:,:,:)
160 : REAL(sp), POINTER :: MOE_LUT06 (:,:,:,:,:,:,:)
161 : REAL(sp), POINTER :: MOE_LUT10 (:,:,:,:,:,:,:)
162 : REAL(sp), POINTER :: MOE_LUT14 (:,:,:,:,:,:,:)
163 : REAL(sp), POINTER :: MOE_LUT18 (:,:,:,:,:,:,:)
164 : REAL(sp), POINTER :: DNOX_LUT02 (:,:,:,:,:,:,:)
165 : REAL(sp), POINTER :: DNOX_LUT06 (:,:,:,:,:,:,:)
166 : REAL(sp), POINTER :: DNOX_LUT10 (:,:,:,:,:,:,:)
167 : REAL(sp), POINTER :: DNOX_LUT14 (:,:,:,:,:,:,:)
168 : REAL(sp), POINTER :: DNOX_LUT18 (:,:,:,:,:,:,:)
169 :
170 : ! Location and type of look up table data
171 : CHARACTER(LEN=255) :: LutDir
172 : LOGICAL :: IsNc
173 :
174 : TYPE(MyInst), POINTER :: NextInst => NULL()
175 : END TYPE MyInst
176 :
177 : ! Pointer to instances
178 : TYPE(MyInst), POINTER :: AllInst => NULL()
179 :
180 : CONTAINS
181 : !EOC
182 : !------------------------------------------------------------------------------
183 : ! Harmonized Emissions Component (HEMCO) !
184 : !------------------------------------------------------------------------------
185 : !BOP
186 : !
187 : ! !IROUTINE: HCOX_ParaNOx_Run
188 : !
189 : ! !DESCRIPTION: Subroutine HCOX\_ParaNOx\_Run is the driver routine to
190 : ! calculate ship NOx emissions for the current time step. Emissions in
191 : ! [kg/m2/s] are added to the emissions array of the passed
192 : !\\
193 : !\\
194 : ! !INTERFACE:
195 : !
196 0 : SUBROUTINE HCOX_ParaNOx_Run( ExtState, HcoState, RC )
197 : !
198 : ! !USES:
199 : !
200 : USE HCO_Calc_Mod, ONLY : HCO_CalcEmis
201 : !
202 : ! !INPUT PARAMETERS:
203 : !
204 : TYPE(Ext_State), POINTER :: ExtState ! External data fields
205 : !
206 : ! !INPUT/OUTPUT PARAMETERS:
207 : !
208 : TYPE(HCO_State), POINTER :: HcoState ! HEMCO State object
209 : INTEGER, INTENT(INOUT) :: RC ! Success or failure?
210 :
211 : ! !REVISION HISTORY:
212 : ! 06 Aug 2013 - C. Keller - Initial Version
213 : ! See https://github.com/geoschem/hemco for complete history
214 : !EOP
215 : !------------------------------------------------------------------------------
216 : !BOC
217 : !
218 : ! !LOCAL VARIABLES:
219 : !
220 : LOGICAL :: DefScaleEmis
221 : CHARACTER(LEN=255) :: MSG, LOC
222 : TYPE(MyInst), POINTER :: Inst
223 :
224 : !=================================================================
225 : ! HCOX_PARANOX_RUN begins here!
226 : !=================================================================
227 0 : LOC = 'HCOX_PARANOX_RUN (HCOX_PARANOX_MOD.F90)'
228 :
229 : ! Return if extension disabled
230 0 : IF ( ExtState%ParaNOx <= 0 ) RETURN
231 :
232 : ! Enter
233 0 : CALL HCO_ENTER(HcoState%Config%Err, LOC, RC)
234 0 : IF ( RC /= HCO_SUCCESS ) THEN
235 0 : CALL HCO_ERROR( 'ERROR 0', RC, THISLOC=LOC )
236 0 : RETURN
237 : ENDIF
238 :
239 : ! Get local instance
240 0 : Inst => NULL()
241 0 : CALL InstGet ( ExtState%ParaNOx, Inst, RC )
242 0 : IF ( RC /= HCO_SUCCESS ) THEN
243 0 : WRITE(MSG,*) 'Cannot find ParaNOx instance Nr. ', ExtState%ParaNOx
244 0 : CALL HCO_ERROR(MSG,RC)
245 0 : RETURN
246 : ENDIF
247 :
248 : ! ----------------------------------------------------------------
249 : ! Use HEMCO core routines to get ship NO emissions
250 : ! ----------------------------------------------------------------
251 :
252 : ! Prepare HEMCO core run (Hco_CalcEmis):
253 : ! --> Set tracer and category range + extension number.
254 : ! Note: Set species min and max to the full range of species.
255 : ! For the ParaNox extension, emission fields of only one species
256 : ! should be defined. Hco_CalcEmis will exit w/ error if this is
257 : ! not the case.
258 0 : HcoState%Options%SpcMin = 1
259 0 : HcoState%Options%SpcMax = -1
260 0 : HcoState%Options%CatMin = 1
261 0 : HcoState%Options%CatMax = -1
262 0 : HcoState%Options%ExtNr = Inst%ExtNr
263 :
264 : ! --> Define array to write emissions into.
265 0 : Inst%ShipNO = 0.0d0
266 0 : HcoState%Options%AutoFillDiagn = .FALSE.
267 0 : HcoState%Options%FillBuffer = .TRUE.
268 0 : HcoState%Buffer3D%Val => Inst%ShipNO
269 :
270 : ! Calculate ship NO emissions and write them into the ShipNO
271 : ! array [kg/m2/s].
272 0 : CALL HCO_CalcEmis( HcoState, .FALSE., RC )
273 0 : IF ( RC /= HCO_SUCCESS ) THEN
274 0 : CALL HCO_ERROR( 'ERROR 1', RC, THISLOC=LOC )
275 0 : RETURN
276 : ENDIF
277 :
278 : ! Reset settings to standard
279 0 : HcoState%Buffer3D%Val => NULL()
280 0 : HcoState%Options%FillBuffer = .FALSE.
281 0 : HcoState%Options%ExtNr = 0
282 0 : HcoState%Options%AutoFillDiagn = .TRUE.
283 :
284 : ! Calculate production rates of NO, HNO3 and O3 based upon ship NO
285 : ! emissions and add these values to the respective emission
286 : ! arrays.
287 : ! Note: For O3, it is possible to get negative emissions (i.e.
288 : ! deposition), in which case these values will be added to the
289 : ! drydep array.
290 0 : CALL Evolve_Plume( ExtState, Inst%ShipNO, HcoState, Inst, RC )
291 0 : IF ( RC /= HCO_SUCCESS ) THEN
292 0 : CALL HCO_ERROR( 'ERROR 2', RC, THISLOC=LOC )
293 0 : RETURN
294 : ENDIF
295 :
296 : ! Leave w/ success
297 0 : Inst => NULL()
298 0 : CALL HCO_Leave( HcoState%Config%Err, RC )
299 :
300 : END SUBROUTINE HCOX_ParaNOx_Run
301 : !EOC
302 : !------------------------------------------------------------------------------
303 : ! Harmonized Emissions Component (HEMCO) !
304 : !------------------------------------------------------------------------------
305 : !BOP
306 : !
307 : ! !IROUTINE: Evolve_Plume
308 : !
309 : ! !DESCRIPTION: Subroutine EVOLVE\_PLUME performs plume dilution and chemistry
310 : ! of ship NO emissions for every grid box and writes the resulting NO, HNO3
311 : ! and O3 emission (production) rates into State\_Chm%NomixS.
312 : !\\
313 : !\\
314 : ! !INTERFACE:
315 : !
316 0 : SUBROUTINE Evolve_Plume( ExtState, ShipNoEmis, HcoState, Inst, RC )
317 : !
318 : ! !USES:
319 : !
320 : USE HCO_Types_Mod, ONLY : DiagnCont
321 : USE HCO_FluxArr_mod, ONLY : HCO_EmisAdd
322 : USE HCO_FluxArr_mod, ONLY : HCO_DepvAdd
323 : USE HCO_Clock_Mod, ONLY : HcoClock_First
324 : USE HCO_Calc_Mod, ONLY : HCO_CheckDepv
325 : USE HCO_GeoTools_Mod, ONLY : HCO_GetSUNCOS
326 : !
327 : ! !INPUT PARAMETERS:
328 : !
329 : TYPE(Ext_State), POINTER :: ExtState ! External data
330 : !
331 : ! !INPUT/OUTPUT PARAMETERS:
332 : !
333 : REAL(hp), INTENT(INOUT) :: ShipNoEmis(:,:,:) ! Emissions
334 : TYPE(HCO_State), POINTER :: HcoState ! HEMCO State obj
335 : TYPE(MyInst), POINTER :: Inst ! Local instance
336 : INTEGER, INTENT(INOUT) :: RC ! Success or failure
337 : !
338 : ! !REVISION HISTORY:
339 : ! 06 Aug 2013 - C. Keller - Initial Version
340 : ! See https://github.com/geoschem/hemco for complete history
341 : !EOP
342 : !------------------------------------------------------------------------------
343 : !BOC
344 : !
345 : ! !LOCAL VARIABLES:
346 : !
347 : INTEGER :: I, J, L
348 : LOGICAL :: ERR
349 : LOGICAL :: FILLED
350 : LOGICAL :: FIRST
351 : LOGICAL :: DefScaleEmis
352 : REAL(hp) :: iFlx, TMP
353 : CHARACTER(LEN=255) :: MSG, LOC
354 : CHARACTER(LEN=1) :: CHAR1
355 :
356 : ! Arrays
357 0 : REAL(hp), TARGET :: FLUXNO (HcoState%NX,HcoState%NY)
358 0 : REAL(hp), TARGET :: FLUXNO2 (HcoState%NX,HcoState%NY)
359 0 : REAL(hp), TARGET :: FLUXHNO3(HcoState%NX,HcoState%NY)
360 0 : REAL(hp), TARGET :: FLUXO3 (HcoState%NX,HcoState%NY)
361 : !%%% Comment out unused code
362 : !%%%! REAL(hp), TARGET :: DEPO3 (HcoState%NX,HcoState%NY)
363 : !%%%! REAL(hp), TARGET :: DEPHNO3 (HcoState%NX,HcoState%NY)
364 :
365 : ! Pointers
366 0 : REAL(hp), POINTER :: Arr2D(:,:)
367 :
368 : ! For diagnostics
369 0 : REAL(hp), TARGET :: DIAGN (HcoState%NX,HcoState%NY,5)
370 : LOGICAL, SAVE :: DODIAGN = .FALSE.
371 : CHARACTER(LEN=31) :: DiagnName
372 : TYPE(DiagnCont), POINTER :: TmpCnt
373 :
374 : ! Paranox update
375 : REAL(dp) :: SHIP_FNOx, SHIP_DNOx, SHIP_OPE, SHIP_MOE
376 : REAL(dp) :: FNO_NOx
377 : REAL(hp) :: iMass
378 : REAL(hp) :: ExpVal
379 : !%%% Comment out debug code
380 : !%%%! ! testing only
381 : !%%%! REAL*8 :: FRAC, TOTPRES, DELTPRES
382 : !%%%! INTEGER :: TOP
383 : !%%%! integer :: ix, jx
384 : !%%%! logical, parameter :: add2hemco = .true.
385 :
386 : !=================================================================
387 : ! EVOLVE_PLUME begins here!
388 : !=================================================================
389 0 : LOC = 'EVOLVE_PLUE (HCOX_PARANOX_MOD.F90)'
390 :
391 : ! Enter
392 0 : CALL HCO_ENTER(HcoState%Config%Err, LOC, RC)
393 0 : IF ( RC /= HCO_SUCCESS ) THEN
394 0 : CALL HCO_ERROR( 'ERROR 3', RC, THISLOC=LOC )
395 0 : RETURN
396 : ENDIF
397 :
398 : ! Leave here if none of the tracers defined
399 0 : IF ( Inst%IDTNO <= 0 .AND. Inst%IDTO3 <= 0 .AND. Inst%IDTHNO3 <= 0 ) THEN
400 0 : RC = HCO_SUCCESS
401 0 : RETURN
402 : ENDIF
403 :
404 : ! Nullify
405 0 : Arr2D => NULL()
406 0 : TmpCnt => NULL()
407 :
408 : ! ------------------------------------------------------------------
409 : ! First call: check for diagnostics to write and fill restart values
410 : ! ------------------------------------------------------------------
411 0 : FIRST = HcoClock_First( HcoState%Clock, .TRUE. )
412 :
413 0 : IF ( FIRST ) THEN
414 : ! See if we have to write out manual diagnostics
415 0 : IF ( .NOT. DoDiagn ) THEN
416 0 : DiagnName = 'PARANOX_NOXFRAC_REMAINING'
417 : CALL DiagnCont_Find ( HcoState%Diagn, -1, -1, -1, -1, -1, &
418 0 : DiagnName, 0, DoDiagn, TmpCnt )
419 0 : TmpCnt => NULL()
420 : ENDIF
421 0 : IF ( .NOT. DoDiagn ) THEN
422 0 : DiagnName = 'PARANOX_O3_PRODUCTION'
423 : CALL DiagnCont_Find ( HcoState%Diagn, -1, -1, -1, -1, -1, &
424 0 : DiagnName, 0, DoDiagn, TmpCnt )
425 0 : TmpCnt => NULL()
426 : ENDIF
427 0 : IF ( .NOT. DoDiagn ) THEN
428 0 : DiagnName = 'PARANOX_NO_PRODUCTION'
429 : CALL DiagnCont_Find ( HcoState%Diagn, -1, -1, -1, -1, -1, &
430 0 : DiagnName, 0, DoDiagn, TmpCnt )
431 0 : TmpCnt => NULL()
432 : ENDIF
433 0 : IF ( .NOT. DoDiagn ) THEN
434 0 : DiagnName = 'PARANOX_TOTAL_SHIPNOX'
435 : CALL DiagnCont_Find ( HcoState%Diagn, -1, -1, -1, -1, -1, &
436 0 : DiagnName, 0, DoDiagn, TmpCnt )
437 0 : TmpCnt => NULL()
438 : ENDIF
439 0 : IF ( .NOT. DoDiagn ) THEN
440 0 : DiagnName = 'PARANOX_OPE'
441 : CALL DiagnCont_Find ( HcoState%Diagn, -1, -1, -1, -1, -1, &
442 0 : DiagnName, 0, DoDiagn, TmpCnt )
443 0 : TmpCnt => NULL()
444 : ENDIF
445 : ENDIF
446 :
447 0 : IF ( DoDiagn ) DIAGN(:,:,:) = 0.0_hp
448 :
449 : ! ------------------------------------------------------------------
450 : ! Update SC5
451 : ! ------------------------------------------------------------------
452 : ! SC5 holds the SUNCOS values of 5 hours ago.
453 0 : CALL HCO_getSUNCOS( HcoState, Inst%SC5, -5, RC )
454 0 : IF ( RC /= HCO_SUCCESS ) THEN
455 0 : CALL HCO_ERROR( 'ERROR 4', RC, THISLOC=LOC )
456 0 : RETURN
457 : ENDIF
458 :
459 : ! Error check
460 0 : ERR = .FALSE.
461 :
462 : ! Initialize
463 0 : FLUXNO = 0.0_hp
464 0 : FLUXNO2 = 0.0_hp
465 0 : FLUXHNO3 = 0.0_hp
466 0 : FLUXO3 = 0.0_hp
467 :
468 : ! Deposition fluxes
469 0 : Inst%DEPO3 = 0.0_sp
470 0 : Inst%DEPHNO3 = 0.0_sp
471 :
472 : !%%% Comment out debug code
473 : !%%%! ! Debug
474 : !%%%! print*, '### In EVOLVE_PLUME:'
475 : !%%%! print*, '### JOH: ', SUM ( ExtState%JOH%Arr%Val ), &
476 : !%%%! MAXVAL( ExtState%JOH%Arr%Val )
477 : !%%%! print*, '### JNO2: ', SUM ( ExtState%JNO2%Arr%Val ), &
478 : !%%%! MAXVAL( ExtState%JNO2%Arr%Val )
479 : !%%%! print*, '### SC5 : ', SUM ( SC5 ), MAXVAL(SC5)
480 : !%%%! print*, '### EMIS: ', SUM ( SHIPNOEMIS(:,:,1) ), MAXVAL(SHIPNOEMIS(:,:,1))
481 :
482 : ! Loop over all grid boxes
483 : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
484 : ! Note: there seems to be a problem with the OMP loop in that
485 : ! the species concentrations (O3molec, NOmolec, NO2molec)
486 : ! differ slightly in a few grid boxes. Don't know exactly what
487 : ! is going on here, but uncomment for now! Needs more
488 : ! evaluation and testing.
489 : !
490 : ! Now use #if defined( 0 ) to block of this code (bmy, 6/6/14)
491 : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
492 : !!!$OMP PARALLEL DO &
493 : !!!$OMP DEFAULT( SHARED ) &
494 : !!!$OMP PRIVATE( I, J, L, MSG, iFlx, iMass, TMP ) &
495 : !!!$OMP PRIVATE( SHIP_FNOx, SHIP_DNOx, SHIP_OPE, SHIP_MOE, FNO_NOx ) &
496 : !!!$OMP SCHEDULE( DYNAMIC )
497 0 : DO J = 1, HcoState%NY
498 0 : DO I = 1, HcoState%NX
499 :
500 : ! Zero private variables for safety's sake
501 0 : FNO_NOx = 0.0_dp
502 0 : iFlx = 0.0_hp
503 0 : iMass = 0.0_hp
504 0 : SHIP_FNOx = 0.0_dp
505 0 : SHIP_DNOx = 0.0_dp
506 0 : SHIP_OPE = 0.0_dp
507 0 : SHIP_MOE = 0.0_dp
508 0 : TMP = 0.0_hp
509 :
510 : !---------------------------------------------------------------------
511 : ! Skip if no ship emissions in this grid box
512 : !---------------------------------------------------------------------
513 0 : IF ( .not. ( ShipNoEmis(I,J,1) > 0.0_hp ) ) CYCLE
514 :
515 : !---------------------------------------------------------------------
516 : ! Production Efficiency for ship emiss (gvinken,mpayer,2/7/12)
517 : ! Updated to include effects of wind speed (cdh, 3/25/2014)
518 : ! Updated for HEMCO (ckeller, 02/04/2015)
519 : !---------------------------------------------------------------------
520 : CALL PARANOX_LUT( ExtState, HcoState, Inst, I, &
521 : J, RC, SHIP_FNOx, SHIP_DNOx, &
522 0 : SHIP_OPE, SHIP_MOE )
523 0 : IF ( RC /= HCO_SUCCESS ) THEN
524 : ERR = .TRUE.; EXIT
525 : ENDIF
526 :
527 : !%%% Comment out debug code
528 : !%%%! ! for debugging only
529 : !%%%! if(I==3.and.J==35)then
530 : !%%%! write(*,*) 'PARANOX ship emissions @',I,J
531 : !%%%! write(*,*) 'Emis [kg/m2/s]: ', ShipNoEmis(I,J,1)
532 : !%%%! write(*,*) 'SHIP_FNOx: ', SHIP_FNOx
533 : !%%%! write(*,*) 'SHIP_DNOx: ', SHIP_DNOx
534 : !%%%! write(*,*) 'SHIP_OPE : ', SHIP_OPE
535 : !%%%! write(*,*) 'SHIP_MOE : ', SHIP_MOE
536 : !%%%! endif
537 :
538 : !---------------------------------------------------------------------
539 : ! Split the ship NOx emission into NO and NO2
540 : ! following the ambient ratio
541 : ! Now check for zero ambient air concentrations. In this
542 : ! case, arbitrarily emit everything as NO (ckeller, 04/10/15).
543 : !
544 : ! Bug fix: Make sure NO and NO2 are both positive
545 : ! See https://github.com/geoschem/geos-chem/issues/380
546 : ! -- Bob Yantosca (24 Jul 2020)
547 : !---------------------------------------------------------------------
548 0 : FNO_NOx = 1.0_hp
549 0 : IF ( ExtState%NO%Arr%Val (I,J,1) > 0.0_hp .and. &
550 : ExtState%NO2%Arr%Val(I,J,1) > 0.0_hp ) THEN
551 : FNO_NOx = ( ExtState%NO%Arr%Val(I,J,1) / Inst%MW_NO ) &
552 : / ( ( ExtState%NO%Arr%Val(I,J,1) / Inst%MW_NO ) &
553 0 : + ( ExtState%NO2%Arr%Val(I,J,1) / Inst%MW_NO2 ) )
554 : ENDIF
555 :
556 : !---------------------------------------------------------------------
557 : ! Calculate NO emissions
558 : !---------------------------------------------------------------------
559 0 : IF ( Inst%IDTNO > 0 ) THEN
560 :
561 : ! Of the total ship NOx, the fraction SHIP_FNOx
562 : ! survives after plume dilution and chemistry.
563 : ! FNO_NOx is the ratio of NO / NOx.
564 : ! Unit: kg/m2/s
565 0 : FLUXNO(I,J) = ShipNoEmis(I,J,1) * SHIP_FNOx * FNO_NOx
566 :
567 : ENDIF
568 :
569 : !---------------------------------------------------------------------
570 : ! Calculate NO2 emissions
571 : !---------------------------------------------------------------------
572 0 : IF ( Inst%IDTNO2 > 0 ) THEN
573 :
574 : ! NO2 emissions complement NO emissions, so that total NOx
575 : ! emissions are preserved.
576 0 : FLUXNO2(I,J) = ShipNoEmis(I,J,1) &
577 : * SHIP_FNOx &
578 : * ( 1.0d0 - FNO_NOx ) &
579 0 : * ( Inst%MW_NO2 / Inst%MW_NO )
580 : ENDIF
581 :
582 : !---------------------------------------------------------------------
583 : ! Calculate HNO3 emissions
584 : !---------------------------------------------------------------------
585 0 : IF ( Inst%IDTHNO3 > 0 ) THEN
586 :
587 : ! Of the total ship NOx, the fraction 1-SHIP_FNOx-SHIP_DNOx
588 : ! is converted to HNO3 during plume dilution and chemistry.
589 : ! Unit: kg/m2/s
590 0 : FLUXHNO3(I,J) = ShipNoEmis(I,J,1) &
591 : * ( 1d0 - SHIP_FNOx - SHIP_DNOx ) &
592 0 : * ( Inst%MW_HNO3 / Inst%MW_NO )
593 : ENDIF
594 :
595 : !--------------------------------------------------------------------
596 : ! NOy deposition (as HNO3) from the sub-grid plume.
597 : ! The calculated deposition flux is in kg/m2/s, which has to be
598 : ! converted to 1/s. The species mass is either the species mass in
599 : ! the first grid box or of the entire PBL column, depending on the
600 : ! HEMCO setting 'PBL_DRYDEP'.
601 : !
602 : ! As of 4/10/15, exchange loss rates in original units of kg/m2/s.
603 : ! (ckeller)
604 : !--------------------------------------------------------------------
605 0 : IF ( (Inst%IDTHNO3 > 0) .AND. (SHIP_DNOx > 0.0_dp) ) THEN
606 :
607 : ! Deposition flux in kg/m2/s.
608 0 : Inst%DEPHNO3(I,J) = ShipNoEmis(I,J,1) &
609 : * SHIP_DNOx &
610 0 : * ( Inst%MW_HNO3 / Inst%MW_NO )
611 :
612 : !%%% Comment out unused code
613 : !%%%! iFlx = ShipNoEmis(I,J,1) * SHIP_DNOx * ( MW_HNO3 / MW_NO )
614 : !%%%!
615 : !%%%! ! Get mass of species. This can either be the total PBL
616 : !%%%! ! column mass or the first layer only, depending on the
617 : !%%%! ! HEMCO setting.
618 : !%%%! iMass = ExtState%HNO3%Arr%Val(I,J,1) &
619 : !%%%! * ExtState%FRAC_OF_PBL%Arr%Val(I,J,1)
620 : !%%%! IF ( HcoState%Options%PBL_DRYDEP ) THEN
621 : !%%%! DO L = 1, HcoState%NZ
622 : !%%%! IF ( ExtState%FRAC_OF_PBL%Arr%Val(I,J,L) == 0.0_hp ) EXIT
623 : !%%%! iMass = iMass + ( ExtState%HNO3%Arr%Val(I,J,L) * &
624 : !%%%! ExtState%FRAC_OF_PBL%Arr%Val(I,J,L) )
625 : !%%%! ENDDO
626 : !%%%! ENDIF
627 : !%%%!
628 : !%%%! ! Calculate deposition velocity (1/s) from flux
629 : !%%%! ! Now avoid div-zero error (ckeller, 11/10/2014).
630 : !%%%! IF ( iMass > TINY(1.0_hp) ) THEN
631 : !%%%! TMP = ABS(iFlx) * HcoState%Grid%AREA_M2%Val(I,J)
632 : !%%%!
633 : !%%%! ! Check if it's safe to do division
634 : !%%%! IF ( (EXPONENT(TMP)-EXPONENT(iMass)) < MAXEXPONENT(TMP) ) THEN
635 : !%%%! DEPHNO3(I,J) = TMP / iMass
636 : !%%%! ENDIF
637 : !%%%!
638 : !%%%! ! Check deposition velocity
639 : !%%%! CALL HCO_CheckDepv( HcoState, DEPHNO3(I,J), RC )
640 : !%%%! ENDIF
641 :
642 : ENDIF
643 :
644 : !---------------------------------------------------------------------
645 : ! Calculate O3 emissions
646 : !---------------------------------------------------------------------
647 0 : IF ( Inst%IDTO3 > 0 ) THEN
648 :
649 : ! Of the total ship NOx, the fraction
650 : ! (1-FRACTION_NOX)*INT_OPE is converted to O3 during
651 : ! plume dilution and chemistry.
652 : ! Unit: kg/m2/s
653 0 : iFlx = ShipNoEmis(I,J,1) &
654 : * ( 1d0 - SHIP_FNOx ) &
655 : * SHIP_OPE &
656 0 : * ( Inst%MW_O3 / Inst%MW_NO )
657 :
658 : ! For positive fluxes, add to emission flux array
659 0 : IF ( iFlx >= 0.0_hp ) THEN
660 0 : FLUXO3(I,J) = iFlx
661 :
662 : ! For negative fluxes, calculate deposition velocity based
663 : ! on current surface O3 concentration and pass to deposition
664 : ! array. See comment on dry dep calculation of HNO3 above.
665 : ! As of 4/10/15, exchange loss rates in original units of kg/m2/s.
666 : ! (ckeller)
667 : ELSE
668 :
669 : ! Deposition flux in kg/m2/s.
670 : ! Make sure ozone deposition flux is positive (ckeller, 3/29/16).
671 : !DEPO3(I,J) = iFlx
672 0 : Inst%DEPO3(I,J) = ABS(iFlx)
673 :
674 : !%%% Comment out unused code
675 : !%%%! ! Get mass of species. This can either be the total PBL
676 : !%%%! ! column mass or the first layer only, depending on the
677 : !%%%! ! HEMCO setting.
678 : !%%%! iMass = ExtState%O3%Arr%Val(I,J,1) &
679 : !%%%! * ExtState%FRAC_OF_PBL%Arr%Val(I,J,1)
680 : !%%%! IF ( HcoState%Options%PBL_DRYDEP ) THEN
681 : !%%%! DO L = 1, HcoState%NZ
682 : !%%%! IF ( ExtState%FRAC_OF_PBL%Arr%Val(I,J,L) == 0.0_hp ) EXIT
683 : !%%%! iMass = iMass + ( ExtState%O3%Arr%Val(I,J,L) * &
684 : !%%%! ExtState%FRAC_OF_PBL%Arr%Val(I,J,L) )
685 : !%%%! ENDDO
686 : !%%%! ENDIF
687 : !%%%!
688 : !%%%! ! Calculate deposition velocity (1/s) from flux
689 : !%%%! ! Now avoid div-zero error (ckeller, 11/10/2014).
690 : !%%%! IF ( iMass > TINY(1.0_hp) ) THEN
691 : !%%%! TMP = ABS(iFlx) * HcoState%Grid%AREA_M2%Val(I,J)
692 : !%%%!
693 : !%%%! ! Check if it's safe to do division
694 : !%%%! IF ( (EXPONENT(TMP)-EXPONENT(iMass)) < MAXEXPONENT(TMP) ) THEN
695 : !%%%! DEPO3(I,J) = TMP / iMass
696 : !%%%! ENDIF
697 : !%%%!
698 : !%%%! ! Check deposition velocity
699 : !%%%! CALL HCO_CheckDepv( HcoState, DEPO3(I,J), RC )
700 : !%%%! ENDIF
701 :
702 : ENDIF
703 : ENDIF
704 :
705 : !---------------------------------------------------------------------
706 : ! Eventually write out into diagnostics array
707 : !---------------------------------------------------------------------
708 0 : IF ( DoDiagn ) THEN
709 0 : DIAGN(I,J,1) = SHIP_FNOx
710 0 : DIAGN(I,J,2) = SHIP_OPE
711 0 : DIAGN(I,J,3) = FLUXO3(I,J)
712 0 : DIAGN(I,J,4) = ShipNoEmis(I,J,1)
713 0 : DIAGN(I,J,5) = FLUXNO(I,J)
714 : ENDIF
715 :
716 : !---------------------------------------------------------------------
717 : ! Reset ship NO emissions to zero. Will be refilled on next
718 : ! emission step!
719 : !---------------------------------------------------------------------
720 0 : ShipNoEmis(I,J,1) = 0.0d0
721 :
722 : ENDDO !I
723 : ENDDO !J
724 : !!!$OMP END PARALLEL DO
725 :
726 : ! Error check
727 0 : IF ( ERR ) THEN
728 0 : RC = HCO_FAIL
729 0 : RETURN
730 : ENDIF
731 :
732 : !------------------------------------------------------------------------
733 : ! ! Debug
734 : ! print*, '### In EVOLVE_PLUME (B):'
735 : ! print*, '### DIAG 1: ', SUM ( DIAGN(:,:,1) ), &
736 : ! MAXVAL( DIAGN(:,:,1) )
737 : ! print*, '### DIAG 2: ', SUM ( DIAGN(:,:,2) ), &
738 : ! MAXVAL( DIAGN(:,:,2) )
739 : ! print*, '### DIAG 3: ', SUM ( DIAGN(:,:,3) ), &
740 : ! MAXVAL( DIAGN(:,:,3) )
741 : ! print*, '### DIAG 4: ', SUM ( DIAGN(:,:,4) ), &
742 : ! MAXVAL( DIAGN(:,:,4) )
743 : ! print*, '### DIAG 5: ', SUM ( DIAGN(:,:,5) ), &
744 : ! MAXVAL( DIAGN(:,:,5) )
745 : !------------------------------------------------------------------------
746 :
747 : !=======================================================================
748 : ! PASS TO HEMCO STATE AND UPDATE DIAGNOSTICS
749 : !=======================================================================
750 :
751 : ! Turn off emission scaling. We don't want the computed fluxes to be
752 : ! scaled any more. If a uniform scale factor is defined for NO, it
753 : ! has been applied to the ship NO emissions already (ckeller, 5/11/17).
754 0 : DefScaleEmis = HcoState%Options%ScaleEmis
755 0 : HcoState%Options%ScaleEmis = .FALSE.
756 :
757 : ! NO
758 0 : IF ( Inst%IDTNO > 0 ) THEN
759 :
760 : ! Add flux to emission array
761 : CALL HCO_EmisAdd( HcoState, FLUXNO, Inst%IDTNO, &
762 0 : RC, ExtNr=Inst%ExtNr )
763 0 : IF ( RC /= HCO_SUCCESS ) THEN
764 0 : CALL HCO_ERROR( 'HCO_EmisAdd error: FLUXNO', RC )
765 0 : RETURN
766 : ENDIF
767 : ENDIF
768 :
769 : ! NO2
770 0 : IF ( Inst%IDTNO2 > 0 ) THEN
771 :
772 : ! Add flux to emission array
773 : CALL HCO_EmisAdd( HcoState, FLUXNO2, Inst%IDTNO2, &
774 0 : RC, ExtNr=Inst%ExtNr )
775 : ENDIF
776 :
777 : ! HNO3
778 0 : IF ( Inst%IDTHNO3 > 0 ) THEN
779 :
780 : ! Add flux to emission array
781 : CALL HCO_EmisAdd( HcoState, FLUXHNO3, Inst%IDTHNO3, &
782 0 : RC, ExtNr=Inst%ExtNr )
783 : ENDIF
784 :
785 : ! O3
786 0 : IF ( Inst%IDTO3 > 0 ) THEN
787 :
788 : ! Add flux to emission array (kg/m2/s)
789 : CALL HCO_EmisAdd( HcoState, FLUXO3, Inst%IDTO3, &
790 0 : RC, ExtNr=Inst%ExtNr )
791 : ENDIF
792 :
793 :
794 : ! Eventually update manual diagnostics
795 0 : IF ( DoDiagn ) THEN
796 0 : DiagnName = 'PARANOX_NOXFRAC_REMAINING'
797 0 : Arr2D => DIAGN(:,:,1)
798 : CALL Diagn_Update( HcoState, ExtNr=Inst%ExtNr, &
799 0 : cName=TRIM(DiagnName), Array2D=Arr2D, RC=RC)
800 0 : IF ( RC /= HCO_SUCCESS ) THEN
801 0 : CALL HCO_ERROR( 'ERROR 5', RC, THISLOC=LOC )
802 0 : RETURN
803 : ENDIF
804 0 : Arr2D => NULL()
805 :
806 0 : DiagnName = 'PARANOX_OPE'
807 0 : Arr2D => DIAGN(:,:,2)
808 : CALL Diagn_Update( HcoState, ExtNr=Inst%ExtNr, &
809 0 : cName=TRIM(DiagnName), Array2D=Arr2D, RC=RC)
810 0 : IF ( RC /= HCO_SUCCESS ) THEN
811 0 : CALL HCO_ERROR( 'ERROR 6', RC, THISLOC=LOC )
812 0 : RETURN
813 : ENDIF
814 0 : Arr2D => NULL()
815 :
816 0 : DiagnName = 'PARANOX_O3_PRODUCTION'
817 0 : Arr2D => DIAGN(:,:,3)
818 : CALL Diagn_Update( HcoState, ExtNr=Inst%ExtNr, &
819 0 : cName=TRIM(DiagnName), Array2D=Arr2D, RC=RC)
820 0 : IF ( RC /= HCO_SUCCESS ) THEN
821 0 : CALL HCO_ERROR( 'ERROR 7', RC, THISLOC=LOC )
822 0 : RETURN
823 : ENDIF
824 0 : Arr2D => NULL()
825 :
826 0 : DiagnName = 'PARANOX_TOTAL_SHIPNOX'
827 0 : Arr2D => DIAGN(:,:,4)
828 : CALL Diagn_Update( HcoState, ExtNr=Inst%ExtNr, &
829 0 : cName=TRIM(DiagnName), Array2D=Arr2D, RC=RC)
830 0 : IF ( RC /= HCO_SUCCESS ) THEN
831 0 : CALL HCO_ERROR( 'ERROR 8', RC, THISLOC=LOC )
832 0 : RETURN
833 : ENDIF
834 0 : Arr2D => NULL()
835 :
836 0 : DiagnName = 'PARANOX_NO_PRODUCTION'
837 0 : Arr2D => DIAGN(:,:,5)
838 : CALL Diagn_Update( HcoState, ExtNr=Inst%ExtNr, &
839 0 : cName=TRIM(DiagnName), Array2D=Arr2D, RC=RC)
840 0 : IF ( RC /= HCO_SUCCESS ) THEN
841 0 : CALL HCO_ERROR( 'ERROR 9', RC, THISLOC=LOC )
842 0 : RETURN
843 : ENDIF
844 0 : Arr2D => NULL()
845 : ENDIF
846 :
847 : ! Reset option ScaleEmis to default value
848 0 : HcoState%Options%ScaleEmis = DefScaleEmis
849 :
850 : ! Return w/ success
851 0 : CALL HCO_LEAVE( HcoState%Config%Err,RC )
852 :
853 0 : END SUBROUTINE Evolve_Plume
854 : !EOC
855 : !------------------------------------------------------------------------------
856 : ! Harmonized Emissions Component (HEMCO) !
857 : !------------------------------------------------------------------------------
858 : !BOP
859 : !
860 : ! !IROUTINE: HCOX_ParaNOx_Init
861 : !
862 : ! !DESCRIPTION: Subroutine HcoX\_ParaNOx\_Init initializes the HEMCO
863 : ! PARANOX extension.
864 : !\\
865 : !\\
866 : ! !INTERFACE:
867 : !
868 0 : SUBROUTINE HCOX_ParaNOx_Init( HcoState, ExtName, ExtState, RC )
869 : !
870 : ! !USES:
871 : !
872 : USE HCO_Chartools_Mod, ONLY : HCO_CharParse
873 : USE HCO_State_MOD, ONLY : HCO_GetHcoID
874 : USE HCO_State_MOD, ONLY : HCO_GetExtHcoID
875 : USE HCO_ExtList_Mod, ONLY : GetExtNr
876 : USE HCO_ExtList_Mod, ONLY : GetExtOpt
877 : USE HCO_Restart_Mod, ONLY : HCO_RestartDefine
878 : ! USE ParaNOx_Util_Mod, ONLY : Read_ParaNOx_LUT
879 : !
880 : ! !INPUT PARAMETERS:
881 : !
882 : CHARACTER(LEN=*), INTENT(IN ) :: ExtName ! Extension name
883 : TYPE(Ext_State), POINTER :: ExtState ! Module options
884 : !
885 : ! !INPUT/OUTPUT PARAMETERS:
886 : !
887 : TYPE(HCO_State), POINTER :: HcoState ! HEMCO state object
888 : INTEGER, INTENT(INOUT) :: RC ! Success or failure?
889 : !
890 : ! !REVISION HISTORY:
891 : ! 06 Aug 2013 - C. Keller - Initial Version
892 : ! See https://github.com/geoschem/hemco for complete history
893 : !EOP
894 : !------------------------------------------------------------------------------
895 : !BOC
896 : !
897 : ! !LOCAL VARIABLES:
898 : !
899 : INTEGER :: ExtNr, I, NN, tmpID, nSpc
900 0 : INTEGER, ALLOCATABLE :: HcoIDs(:)
901 0 : CHARACTER(LEN=31), ALLOCATABLE :: SpcNames(:)
902 : CHARACTER(LEN=31) :: Dummy
903 : CHARACTER(LEN=31) :: DiagnName
904 : CHARACTER(LEN=255) :: MSG, LOC
905 : CHARACTER(LEN= 1) :: CHAR1
906 : TYPE(MyInst), POINTER :: Inst
907 :
908 : !========================================================================
909 : ! HCOX_PARANOX_INIT begins here!
910 : !========================================================================
911 0 : LOC = 'HCOX_PARANOX_INIT (HCOX_PARANOX_MOD.F90)'
912 :
913 : ! Assume success
914 0 : RC = HCO_SUCCESS
915 :
916 : ! Extension Nr.
917 0 : ExtNr = GetExtNr( HcoState%Config%ExtList, TRIM(ExtName) )
918 0 : IF ( ExtNr <= 0 ) RETURN
919 :
920 : ! Enter
921 0 : CALL HCO_ENTER( HcoState%Config%Err, LOC, RC )
922 0 : IF ( RC /= HCO_SUCCESS ) THEN
923 0 : CALL HCO_ERROR( 'ERROR 10', RC, THISLOC=LOC )
924 0 : RETURN
925 : ENDIF
926 :
927 : ! Create local instance
928 0 : Inst => NULL()
929 0 : CALL InstCreate( ExtNr, ExtState%ParaNOx, Inst, RC )
930 0 : IF ( RC /= HCO_SUCCESS ) THEN
931 : CALL HCO_ERROR( &
932 0 : 'Cannot create ParaNOx instance', RC )
933 0 : RETURN
934 : ENDIF
935 :
936 : !========================================================================
937 : ! Skip the following for GEOS-Chem dry-run or HEMCO-standalone dry-run
938 : !========================================================================
939 0 : IF ( .not. HcoState%Options%IsDryRun ) THEN
940 :
941 : !---------------------------------------------------------------------
942 : ! Initialize fields of Inst object for safety's sake (bmy, 10/17/18)
943 : !---------------------------------------------------------------------
944 0 : Inst%IDTNO = -1
945 0 : Inst%IDTNO2 = -1
946 0 : Inst%IDTO3 = -1
947 0 : Inst%IDTHNO3 = -1
948 0 : Inst%MW_O3 = 0.0d0
949 0 : Inst%MW_NO = 0.0d0
950 0 : Inst%MW_NO2 = 0.0d0
951 0 : Inst%MW_HNO3 = 0.0d0
952 0 : Inst%MW_AIR = 0.0d0
953 0 : Inst%Tlev = 0.0e0
954 0 : Inst%JNO2lev = 0.0e0
955 0 : Inst%O3lev = 0.0e0
956 0 : Inst%SEA0lev = 0.0e0
957 0 : Inst%SEA5lev = 0.0e0
958 0 : Inst%JRATIOlev = 0.0e0
959 0 : Inst%NOXlev = 0.0e0
960 0 : Inst%WSlev = 0.0e0
961 0 : Inst%ShipNO => NULL()
962 0 : Inst%SC5 => NULL()
963 0 : Inst%DEPO3 => NULL()
964 0 : Inst%DEPHNO3 => NULL()
965 0 : Inst%FRACNOX_LUT02 => NULL()
966 0 : Inst%FRACNOX_LUT06 => NULL()
967 0 : Inst%FRACNOX_LUT10 => NULL()
968 0 : Inst%FRACNOX_LUT14 => NULL()
969 0 : Inst%FRACNOX_LUT18 => NULL()
970 0 : Inst%OPE_LUT02 => NULL()
971 0 : Inst%OPE_LUT06 => NULL()
972 0 : Inst%OPE_LUT10 => NULL()
973 0 : Inst%OPE_LUT14 => NULL()
974 0 : Inst%OPE_LUT18 => NULL()
975 0 : Inst%MOE_LUT02 => NULL()
976 0 : Inst%MOE_LUT06 => NULL()
977 0 : Inst%MOE_LUT10 => NULL()
978 0 : Inst%MOE_LUT14 => NULL()
979 0 : Inst%MOE_LUT18 => NULL()
980 0 : Inst%DNOX_LUT02 => NULL()
981 0 : Inst%DNOX_LUT06 => NULL()
982 0 : Inst%DNOX_LUT10 => NULL()
983 0 : Inst%DNOX_LUT14 => NULL()
984 0 : Inst%DNOX_LUT18 => NULL()
985 :
986 : !------------------------------------------------------------------------
987 : ! Get species IDs
988 : !------------------------------------------------------------------------
989 :
990 : ! Get HEMCO species IDs
991 0 : CALL HCO_GetExtHcoID( HcoState, ExtNr, HcoIDs, SpcNames, nSpc, RC )
992 0 : IF ( RC /= HCO_SUCCESS ) THEN
993 0 : CALL HCO_ERROR( 'ERROR 11', RC, THISLOC=LOC )
994 0 : RETURN
995 : ENDIF
996 :
997 : ! Check for NO, NO2, O3, and HNO3
998 0 : DO I = 1, nSpc
999 0 : SELECT CASE ( TRIM(SpcNames(I)) )
1000 : CASE ( "NO" )
1001 0 : Inst%IDTNO = HcoIDs(I)
1002 : CASE ( "NO2" )
1003 0 : Inst%IDTNO2 = HcoIDs(I)
1004 : CASE ( "O3" )
1005 0 : Inst%IDTO3 = HcoIDs(I)
1006 : CASE ( "HNO3" )
1007 0 : Inst%IDTHNO3 = HcoIDs(I)
1008 : CASE DEFAULT
1009 : ! leave empty
1010 : END SELECT
1011 : ENDDO
1012 :
1013 : ! Get MW of all species. If species not set in the configuration
1014 : ! file (e.g. they are not being used by PARANOx), determine MW from
1015 : ! default values.
1016 :
1017 : ! O3
1018 0 : IF ( Inst%IDTO3 <= 0 ) THEN
1019 0 : tmpID = HCO_GetHcoID('O3', HcoState )
1020 0 : MSG = 'O3 not produced/removed in PARANOX'
1021 0 : CALL HCO_WARNING(HcoState%Config%Err, MSG, RC )
1022 : ELSE
1023 : tmpID = Inst%IDTO3
1024 : ENDIF
1025 0 : IF ( tmpID > 0 ) THEN
1026 0 : Inst%MW_O3 = HcoState%Spc(tmpID)%MW_g
1027 : ELSE
1028 0 : MSG = 'Use default O3 molecular weight of 48g/mol'
1029 0 : CALL HCO_WARNING(HcoState%Config%Err, MSG, RC )
1030 0 : Inst%MW_O3 = 48.0_dp
1031 : ENDIF
1032 :
1033 : ! NO
1034 0 : IF ( Inst%IDTNO <= 0 ) THEN
1035 0 : tmpID = HCO_GetHcoID('NO', HcoState )
1036 0 : MSG = 'NO not produced in PARANOX'
1037 0 : CALL HCO_WARNING(HcoState%Config%Err, MSG, RC )
1038 : ELSE
1039 : tmpID = Inst%IDTNO
1040 : ENDIF
1041 0 : IF ( tmpID > 0 ) THEN
1042 0 : Inst%MW_NO = HcoState%Spc(tmpID)%MW_g
1043 : ELSE
1044 0 : MSG = 'Use default NO molecular weight of 30g/mol'
1045 0 : CALL HCO_WARNING(HcoState%Config%Err, MSG, RC )
1046 0 : Inst%MW_NO = 30.0_dp
1047 : ENDIF
1048 :
1049 : ! NO2
1050 0 : IF ( Inst%IDTNO2 <= 0 ) THEN
1051 0 : tmpID = HCO_GetHcoID('NO2', HcoState )
1052 0 : MSG = 'NO2 not produced in PARANOX'
1053 0 : CALL HCO_WARNING(HcoState%Config%Err, MSG, RC )
1054 : ELSE
1055 : tmpID = Inst%IDTNO2
1056 : ENDIF
1057 0 : IF ( tmpID > 0 ) THEN
1058 0 : Inst%MW_NO2 = HcoState%Spc(tmpID)%MW_g
1059 : ELSE
1060 0 : MSG = 'Use default NO2 molecular weight of 46g/mol'
1061 0 : CALL HCO_WARNING(HcoState%Config%Err, MSG, RC )
1062 0 : Inst%MW_NO2 = 46.0_dp
1063 : ENDIF
1064 :
1065 : ! HNO3
1066 0 : IF ( Inst%IDTHNO3 <= 0 ) THEN
1067 0 : tmpID = HCO_GetHcoID('HNO3', HcoState )
1068 0 : MSG = 'HNO3 not produced/removed in PARANOX'
1069 0 : CALL HCO_WARNING(HcoState%Config%Err, MSG, RC )
1070 : ELSE
1071 : tmpID = Inst%IDTHNO3
1072 : ENDIF
1073 0 : IF ( tmpID > 0 ) THEN
1074 0 : Inst%MW_HNO3 = HcoState%Spc(tmpID)%MW_g
1075 : ELSE
1076 0 : MSG = 'Use default HNO3 molecular weight of 63g/mol'
1077 0 : CALL HCO_WARNING(HcoState%Config%Err, MSG, RC )
1078 0 : Inst%MW_HNO3 = 63.0_dp
1079 : ENDIF
1080 :
1081 : ! Verbose mode
1082 0 : IF ( HcoState%amIRoot ) THEN
1083 0 : MSG = 'Use ParaNOx ship emissions (extension module)'
1084 0 : CALL HCO_MSG(HcoState%Config%Err,MSG, SEP1='-' )
1085 0 : MSG = ' - Use the following species: (MW, emitted as HEMCO ID) '
1086 0 : CALL HCO_MSG(HcoState%Config%Err,MSG )
1087 0 : WRITE(MSG,"(a,F5.2,I5)") ' NO : ', Inst%MW_NO, Inst%IDTNO
1088 0 : CALL HCO_MSG(HcoState%Config%Err,MSG)
1089 0 : WRITE(MSG,"(a,F5.2,I5)") ' NO2 : ', Inst%MW_NO2, Inst%IDTNO2
1090 0 : CALL HCO_MSG(HcoState%Config%Err,MSG)
1091 0 : WRITE(MSG,"(a,F5.2,I5)") ' O3 : ', Inst%MW_O3, Inst%IDTO3
1092 0 : CALL HCO_MSG(HcoState%Config%Err,MSG)
1093 0 : WRITE(MSG,"(a,F5.2,I5)") ' HNO3: ', Inst%MW_HNO3, Inst%IDTHNO3
1094 0 : CALL HCO_MSG(HcoState%Config%Err,MSG)
1095 : ENDIF
1096 :
1097 : !--------------------------------
1098 : ! Allocate module arrays
1099 : !--------------------------------
1100 :
1101 : ! FNOX
1102 0 : ALLOCATE( Inst%FRACNOX_LUT02(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1103 0 : IF ( RC /= HCO_SUCCESS ) THEN
1104 0 : CALL HCO_ERROR ( 'FRACNOX_LUT02', RC )
1105 0 : RETURN
1106 : ENDIF
1107 0 : Inst%FRACNOX_LUT02 = 0.0_sp
1108 :
1109 0 : ALLOCATE( Inst%FRACNOX_LUT06(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1110 0 : IF ( RC /= HCO_SUCCESS ) THEN
1111 0 : CALL HCO_ERROR ( 'FRACNOX_LUT06', RC )
1112 0 : RETURN
1113 : ENDIF
1114 0 : Inst%FRACNOX_LUT06 = 0.0_sp
1115 :
1116 0 : ALLOCATE( Inst%FRACNOX_LUT10(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1117 0 : IF ( RC /= HCO_SUCCESS ) THEN
1118 0 : CALL HCO_ERROR ( 'FRACNOX_LUT10', RC )
1119 0 : RETURN
1120 : ENDIF
1121 0 : Inst%FRACNOX_LUT10 = 0.0_sp
1122 :
1123 0 : ALLOCATE( Inst%FRACNOX_LUT14(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1124 0 : IF ( RC /= HCO_SUCCESS ) THEN
1125 0 : CALL HCO_ERROR ( 'FRACNOX_LUT014', RC )
1126 0 : RETURN
1127 : ENDIF
1128 0 : Inst%FRACNOX_LUT14 = 0.0_sp
1129 :
1130 0 : ALLOCATE( Inst%FRACNOX_LUT18(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1131 0 : IF ( RC /= HCO_SUCCESS ) THEN
1132 0 : CALL HCO_ERROR ( 'FRACNOX_LUT18', RC )
1133 0 : RETURN
1134 : ENDIF
1135 0 : Inst%FRACNOX_LUT18 = 0.0_sp
1136 :
1137 : ! OPE
1138 0 : ALLOCATE( Inst%OPE_LUT02(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1139 0 : IF ( RC /= HCO_SUCCESS ) THEN
1140 0 : CALL HCO_ERROR ( 'OPE_LUT02', RC )
1141 0 : RETURN
1142 : ENDIF
1143 0 : Inst%OPE_LUT02 = 0.0_sp
1144 :
1145 0 : ALLOCATE( Inst%OPE_LUT06(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1146 0 : IF ( RC /= HCO_SUCCESS ) THEN
1147 0 : CALL HCO_ERROR ( 'OPE_LUT06', RC )
1148 0 : RETURN
1149 : ENDIF
1150 0 : Inst%OPE_LUT06 = 0.0_sp
1151 :
1152 0 : ALLOCATE( Inst%OPE_LUT10(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1153 0 : IF ( RC /= 0 ) THEN
1154 0 : CALL HCO_ERROR ( 'OPE_LUT10', RC )
1155 0 : RETURN
1156 : ENDIF
1157 0 : Inst%OPE_LUT10 = 0.0_sp
1158 :
1159 0 : ALLOCATE( Inst%OPE_LUT14(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1160 0 : IF ( RC /= 0 ) THEN
1161 0 : CALL HCO_ERROR ( 'OPE_LUT014', RC )
1162 0 : RETURN
1163 : ENDIF
1164 0 : Inst%OPE_LUT14 = 0.0_sp
1165 :
1166 0 : ALLOCATE( Inst%OPE_LUT18(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1167 0 : IF ( RC /= HCO_SUCCESS ) THEN
1168 0 : CALL HCO_ERROR ( 'OPE_LUT18', RC )
1169 0 : RETURN
1170 : ENDIF
1171 0 : Inst%OPE_LUT18 = 0.0_sp
1172 :
1173 : ! MOE
1174 0 : ALLOCATE( Inst%MOE_LUT02(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1175 0 : IF ( RC /= HCO_SUCCESS ) THEN
1176 0 : CALL HCO_ERROR ( 'MOE_LUT02', RC )
1177 0 : RETURN
1178 : ENDIF
1179 0 : Inst%MOE_LUT02 = 0.0_sp
1180 :
1181 0 : ALLOCATE( Inst%MOE_LUT06(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1182 0 : IF ( RC /= HCO_SUCCESS ) THEN
1183 0 : CALL HCO_ERROR ( 'MOE_LUT06', RC )
1184 0 : RETURN
1185 : ENDIF
1186 0 : Inst%MOE_LUT06 = 0.0_sp
1187 :
1188 0 : ALLOCATE( Inst%MOE_LUT10(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1189 0 : IF ( RC /= HCO_SUCCESS ) THEN
1190 0 : CALL HCO_ERROR ( 'MOE_LUT10', RC )
1191 0 : RETURN
1192 : ENDIF
1193 0 : Inst%MOE_LUT10 = 0.0_sp
1194 :
1195 0 : ALLOCATE( Inst%MOE_LUT14(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1196 0 : IF ( RC /= HCO_SUCCESS ) THEN
1197 0 : CALL HCO_ERROR ( 'MOE_LUT014', RC )
1198 0 : RETURN
1199 : ENDIF
1200 0 : Inst%MOE_LUT14 = 0.0_sp
1201 :
1202 0 : ALLOCATE( Inst%MOE_LUT18(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1203 0 : IF ( RC /= HCO_SUCCESS ) THEN
1204 0 : CALL HCO_ERROR ( 'MOE_LUT18', RC )
1205 0 : RETURN
1206 : ENDIF
1207 0 : Inst%MOE_LUT18 = 0.0_sp
1208 :
1209 : ! DNOx
1210 0 : ALLOCATE( Inst%DNOx_LUT02(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1211 0 : IF ( RC /= HCO_SUCCESS ) THEN
1212 0 : CALL HCO_ERROR ( 'DNOx_LUT02', RC )
1213 0 : RETURN
1214 : ENDIF
1215 0 : Inst%DNOx_LUT02 = 0.0_sp
1216 :
1217 0 : ALLOCATE( Inst%DNOx_LUT06(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1218 0 : IF ( RC /= HCO_SUCCESS ) THEN
1219 0 : CALL HCO_ERROR ( 'DNOx_LUT06', RC )
1220 0 : RETURN
1221 : ENDIF
1222 0 : Inst%DNOx_LUT06 = 0.0_sp
1223 :
1224 0 : ALLOCATE( Inst%DNOx_LUT10(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1225 0 : IF ( RC /= HCO_SUCCESS ) THEN
1226 0 : CALL HCO_ERROR ( 'DNOx_LUT10', RC )
1227 0 : RETURN
1228 : ENDIF
1229 0 : Inst%DNOx_LUT10 = 0.0_sp
1230 :
1231 0 : ALLOCATE( Inst%DNOx_LUT14(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1232 0 : IF ( RC /= HCO_SUCCESS ) THEN
1233 0 : CALL HCO_ERROR ( 'DNOx_LUT014', RC )
1234 0 : RETURN
1235 : ENDIF
1236 0 : Inst%DNOx_LUT14 = 0.0_sp
1237 :
1238 0 : ALLOCATE( Inst%DNOx_LUT18(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1239 0 : IF ( RC /= HCO_SUCCESS ) THEN
1240 0 : CALL HCO_ERROR ( 'DNOx_LUT18', RC )
1241 0 : RETURN
1242 : ENDIF
1243 0 : Inst%DNOx_LUT18 = 0.0_sp
1244 :
1245 : ALLOCATE(Inst%DEPO3 (HcoState%NX,HcoState%NY), &
1246 0 : Inst%DEPHNO3(HcoState%NX,HcoState%NY), STAT=RC )
1247 0 : IF ( RC /= HCO_SUCCESS ) THEN
1248 0 : CALL HCO_ERROR ( 'Deposition arrays', RC )
1249 0 : RETURN
1250 : ENDIF
1251 0 : Inst%DEPO3 = 0.0_sp
1252 0 : Inst%DEPHNO3 = 0.0_sp
1253 :
1254 : ! ! O3 loss and HNO3 deposition
1255 : ! ALLOCATE( Inst%SHIPO3LOSS(HcoState%NX,HcoState%NY), STAT=RC )
1256 : ! IF ( RC /= HCO_SUCCESS ) THEN
1257 : ! CALL HCO_ERROR ( 'SHIPO3LOSS', RC )
1258 : ! RETURN
1259 : ! ENDIF
1260 : ! Inst%SHIPO3LOSS = 0d0
1261 :
1262 : ! ALLOCATE( Inst%SHIPHNO3DEP(HcoState%NX,HcoState%NY), STAT=RC )
1263 : ! IF ( RC /= HCO_SUCCESS ) THEN
1264 : ! CALL HCO_ERROR ( 'SHIPHNO3DEP', RC ); RETURN
1265 : ! ENDIF
1266 : ! Inst%SHIPHNO3DEP = 0d0
1267 :
1268 : ENDIF
1269 :
1270 : !========================================================================
1271 : ! Initialize the PARANOX look-up tables
1272 : !========================================================================
1273 :
1274 : ! LUT data directory
1275 : CALL GetExtOpt( HcoState%Config, Inst%ExtNr, 'LUT source dir', &
1276 0 : OptValChar=Inst%LutDir, RC=RC)
1277 0 : IF ( RC /= HCO_SUCCESS ) THEN
1278 0 : CALL HCO_ERROR( 'ERROR 12', RC, THISLOC=LOC )
1279 0 : RETURN
1280 : ENDIF
1281 :
1282 : ! Call HEMCO parser to replace tokens such as $ROOT, $MET, or $RES.
1283 : ! There shouldn't be any date token in there ($YYYY, etc.), so just
1284 : ! provide some dummy variables here
1285 0 : CALL HCO_CharParse( HcoState%Config, Inst%LutDir, -999, -1, -1, -1, -1, RC )
1286 0 : IF ( RC /= HCO_SUCCESS ) THEN
1287 0 : CALL HCO_ERROR( 'ERROR 13', RC, THISLOC=LOC )
1288 0 : RETURN
1289 : ENDIF
1290 :
1291 : ! Data format: ncdf (default) or txt
1292 0 : Inst%IsNc = .TRUE.
1293 : CALL GetExtOpt( HcoState%Config, Inst%ExtNr, 'LUT data format', &
1294 0 : OptValChar=Dummy, RC=RC)
1295 0 : IF ( RC /= HCO_SUCCESS ) THEN
1296 0 : CALL HCO_ERROR( 'ERROR 14', RC, THISLOC=LOC )
1297 0 : RETURN
1298 : ENDIF
1299 0 : IF ( TRIM(Dummy) == 'txt' ) Inst%IsNc = .FALSE.
1300 :
1301 : ! Read PARANOX look-up tables from disk. This can be netCDF or txt
1302 : ! format, as determined above.
1303 : !
1304 : ! NOTE: For the GEOS-Chem dry-run or HEMCO-standalone dry-run,
1305 : ! these routines will print file paths to the dry-run log file,
1306 : ! but will not actually read any data.
1307 0 : IF ( Inst%IsNc ) THEN
1308 0 : CALL READ_PARANOX_LUT_NC( HcoState, Inst, RC )
1309 0 : IF ( RC /= HCO_SUCCESS ) THEN
1310 0 : CALL HCO_ERROR( 'ERROR 15', RC, THISLOC=LOC )
1311 0 : RETURN
1312 : ENDIF
1313 : ELSE
1314 0 : CALL READ_PARANOX_LUT_TXT( HcoState, Inst, RC )
1315 0 : IF ( RC /= HCO_SUCCESS ) THEN
1316 0 : CALL HCO_ERROR( 'ERROR 16', RC, THISLOC=LOC )
1317 0 : RETURN
1318 : ENDIF
1319 : ENDIF
1320 :
1321 : !========================================================================
1322 : ! Exit if this is a GEOS-Chem dry-run or HEMCO-standalone dry-run
1323 : !========================================================================
1324 0 : IF ( HcoState%Options%IsDryRun ) THEN
1325 0 : Inst => NULL()
1326 0 : CALL HCO_LEAVE( HcoState%Config%Err,RC )
1327 0 : RETURN
1328 : ENDIF
1329 :
1330 : !========================================================================
1331 : ! Continue initializing PARANOX for regular simulations
1332 : !========================================================================
1333 :
1334 : !------------------------------------------------------------------------
1335 : ! Set other module variables
1336 : !------------------------------------------------------------------------
1337 0 : ALLOCATE ( Inst%ShipNO(HcoState%NX,HcoState%NY,HcoState%NZ), STAT=RC )
1338 0 : IF ( RC /= HCO_SUCCESS ) THEN
1339 0 : CALL HCO_ERROR ( 'ShipNO', RC )
1340 0 : RETURN
1341 : ENDIF
1342 0 : Inst%ShipNO = 0.0_hp
1343 :
1344 : ! Allocate variables for SunCosMid from 5 hours ago.
1345 0 : ALLOCATE ( Inst%SC5(HcoState%NX,HcoState%NY), STAT=RC )
1346 0 : IF ( RC /= HCO_SUCCESS ) THEN
1347 0 : CALL HCO_ERROR ( 'SC5', RC )
1348 0 : RETURN
1349 : ENDIF
1350 0 : Inst%SC5 = 0.0_hp
1351 :
1352 : ! Prompt warning if chemistry time step is more than 60 mins
1353 0 : IF ( HcoState%TS_CHEM > 3600.0_hp ) THEN
1354 0 : IF ( HcoState%amIRoot ) THEN
1355 : MSG = ' Cannot properly store SUNCOS values ' // &
1356 0 : ' because chemistry time step is more than 60 mins!'
1357 0 : CALL HCO_WARNING(HcoState%Config%Err, MSG, RC )
1358 : ENDIF
1359 : ENDIF
1360 :
1361 : ! Molecular weight of AIR
1362 0 : Inst%MW_AIR = HcoState%Phys%AIRMW
1363 :
1364 : !------------------------------------------------------------------------
1365 : ! Define PARANOX diagnostics for the O3 and HNO3 deposition fluxes (in
1366 : ! kg/m2/s).
1367 : !------------------------------------------------------------------------
1368 : CALL Diagn_Create ( HcoState, &
1369 : cName = 'PARANOX_O3_DEPOSITION_FLUX', &
1370 : Trgt2D = Inst%DEPO3, &
1371 : SpaceDim = 2, &
1372 : OutUnit = 'kg/m2/s', &
1373 : COL = HcoState%Diagn%HcoDiagnIDManual, &
1374 0 : RC = RC )
1375 0 : IF ( RC /= HCO_SUCCESS ) THEN
1376 0 : CALL HCO_ERROR( 'ERROR 17', RC, THISLOC=LOC )
1377 0 : RETURN
1378 : ENDIF
1379 :
1380 : CALL Diagn_Create ( HcoState, &
1381 : cName = 'PARANOX_HNO3_DEPOSITION_FLUX', &
1382 : Trgt2D = Inst%DEPHNO3, &
1383 : SpaceDim = 2, &
1384 : OutUnit = 'kg/m2/s', &
1385 : COL = HcoState%Diagn%HcoDiagnIDManual, &
1386 0 : RC = RC )
1387 0 : IF ( RC /= HCO_SUCCESS ) THEN
1388 0 : CALL HCO_ERROR( 'ERROR 18', RC, THISLOC=LOC )
1389 0 : RETURN
1390 : ENDIF
1391 :
1392 : !------------------------------------------------------------------------
1393 : ! Set HEMCO extension variables
1394 : !------------------------------------------------------------------------
1395 :
1396 : ! Met. data required by module
1397 0 : ExtState%O3%DoUse = .TRUE.
1398 0 : ExtState%NO2%DoUse = .TRUE.
1399 0 : ExtState%NO%DoUse = .TRUE.
1400 0 : ExtState%AIR%DoUse = .TRUE.
1401 0 : ExtState%AIRVOL%DoUse = .TRUE.
1402 0 : ExtState%SUNCOS%DoUse = .TRUE.
1403 0 : ExtState%T2M%DoUse = .TRUE.
1404 0 : ExtState%U10M%DoUse = .TRUE.
1405 0 : ExtState%V10M%DoUse = .TRUE.
1406 0 : ExtState%FRAC_OF_PBL%DoUse = .TRUE.
1407 0 : IF ( Inst%IDTHNO3 > 0 ) THEN
1408 0 : ExtState%HNO3%DoUse = .TRUE.
1409 : ENDIF
1410 0 : ExtState%JNO2%DoUse = .TRUE.
1411 0 : ExtState%JOH%DoUse = .TRUE.
1412 :
1413 : !------------------------------------------------------------------------
1414 : ! Leave w/ success
1415 : !------------------------------------------------------------------------
1416 0 : IF ( ALLOCATED(HcoIDs ) ) DEALLOCATE(HcoIDs )
1417 0 : IF ( ALLOCATED(SpcNames) ) DEALLOCATE(SpcNames)
1418 0 : Inst => NULL()
1419 0 : CALL HCO_LEAVE( HcoState%Config%Err,RC )
1420 :
1421 0 : END SUBROUTINE HCOX_ParaNOx_Init
1422 : !EOC
1423 : !------------------------------------------------------------------------------
1424 : ! Harmonized Emissions Component (HEMCO) !
1425 : !------------------------------------------------------------------------------
1426 : !BOP
1427 : !
1428 : ! !IROUTINE: HCOX_ParaNOx_Final
1429 : !
1430 : ! !DESCRIPTION: Subroutine HcoX\_ParaNox\_Final finalizes the HEMCO
1431 : ! PARANOX extension.
1432 : !\\
1433 : !\\
1434 : ! !INTERFACE:
1435 : !
1436 0 : SUBROUTINE HCOX_ParaNOx_Final( HcoState, ExtState, RC )
1437 : !
1438 : ! !USES:
1439 : !
1440 : !
1441 : ! !INPUT PARAMETERS:
1442 : !
1443 : TYPE(HCO_State), POINTER :: HcoState ! HEMCO State obj
1444 : TYPE(Ext_State), POINTER :: ExtState ! Module options
1445 : !
1446 : ! !INPUT/OUTPUT PARAMETERS:
1447 : !
1448 : INTEGER, INTENT(INOUT) :: RC
1449 : !
1450 : ! !REVISION HISTORY:
1451 : ! 06 Aug 2013 - C. Keller - Initial Version
1452 : ! See https://github.com/geoschem/hemco for complete history
1453 : !EOP
1454 : !------------------------------------------------------------------------------
1455 : !BOC
1456 : !
1457 : ! LOCAL VARIABLES:
1458 : !
1459 :
1460 : !=================================================================
1461 : ! HCOX_PARANOX_FINAL begins here!
1462 : !=================================================================
1463 0 : CALL InstRemove( ExtState%ParaNOx )
1464 :
1465 0 : RC = HCO_SUCCESS
1466 :
1467 0 : END SUBROUTINE HCOX_ParaNOx_Final
1468 : !EOC
1469 : !------------------------------------------------------------------------------
1470 : ! GEOS-Chem Global Chemical Transport Model !
1471 : !------------------------------------------------------------------------------
1472 : !BOP
1473 : !
1474 : ! !IROUTINE: read_paranox_lut_nc
1475 : !
1476 : ! !DESCRIPTION: Subroutine READ\_PARANOX\_LUT\_NC reads look-up tables in
1477 : ! netCDF format for use in the PARANOX ship plume model (G.C.M. Vinken)
1478 : !\\
1479 : !\\
1480 : ! !INTERFACE:
1481 : !
1482 0 : SUBROUTINE READ_PARANOX_LUT_NC ( HcoState, Inst, RC )
1483 : !
1484 : ! !USES:
1485 : !
1486 : ! !INPUT ARGUMENTS:
1487 : !
1488 : TYPE(HCO_State), POINTER :: HcoState ! HEMCO State object
1489 : TYPE(MyInst), POINTER :: Inst
1490 : !
1491 : ! !INPUT/OUTPUT ARGUMENTS:
1492 : !
1493 : INTEGER, INTENT(INOUT) :: RC
1494 : !
1495 : ! !REVISION HISTORY:
1496 : ! 06 Feb 2012 - M. Payer - Initial version modified from code provided by
1497 : ! G.C.M. Vinken
1498 : ! See https://github.com/geoschem/hemco for complete history
1499 : !EOP
1500 : !------------------------------------------------------------------------------
1501 : !BOC
1502 : !
1503 : ! !LOCAL VARIABLES:
1504 : INTEGER :: IOS
1505 : CHARACTER(LEN=255) :: FILENAME
1506 : CHARACTER(LEN=255) :: MSG
1507 : INTEGER :: fID
1508 :
1509 : !=================================================================
1510 : ! READ_PARANOX_LUT_NC begins here
1511 : !=================================================================
1512 :
1513 : ! NetCDF reading of PARANOX LUT not supported in ESMF environment
1514 : #if defined(ESMF_)
1515 : MSG = 'In ESMF, cannot read PARANOX look-up-table in netCDF ' // &
1516 : 'format. Please set `LUT data format` to `txt` in the ' // &
1517 : 'HEMCO configuration file.'
1518 : CALL HCO_ERROR(MSG, RC, &
1519 : THISLOC = 'READ_PARANOX_LUT_NC (hcox_paranox_mod.F90)' )
1520 : RETURN
1521 : #else
1522 :
1523 : ! Clear FILENAME
1524 0 : FILENAME = ''
1525 :
1526 : ! FILENAME format string
1527 : 101 FORMAT( a, '/ship_plume_lut_', I2.2, 'ms.nc' )
1528 :
1529 : ! Wind speed levels correspond to the files we will read below
1530 0 : Inst%WSlev = (/ 2.0e0, 6.0e0, 10.0e0, 14.0e0, 18.0e0 /)
1531 :
1532 : ! Read 2m/s LUT
1533 0 : WRITE( FILENAME, 101 ) TRIM(Inst%LutDir), 2
1534 : CALL READ_LUT_NCFILE( HcoState, TRIM( FILENAME ), &
1535 : Inst%FRACNOX_LUT02, Inst%DNOx_LUT02, Inst%OPE_LUT02, Inst%MOE_LUT02, &
1536 : Inst%Tlev, Inst%JNO2lev, Inst%O3lev, Inst%SEA0lev, Inst%SEA5lev, &
1537 0 : Inst%JRATIOlev, Inst%NOXlev, RC=RC)
1538 :
1539 : ! Read 6 m/s LUT
1540 0 : WRITE( FILENAME, 101 ) TRIM(Inst%LutDir), 6
1541 : CALL READ_LUT_NCFILE( HcoState, TRIM( FILENAME ), &
1542 : Inst%FRACNOX_LUT06, Inst%DNOx_LUT06, Inst%OPE_LUT06, &
1543 0 : Inst%MOE_LUT06, RC=RC )
1544 :
1545 : ! Read 10 m/s LUT
1546 0 : WRITE( FILENAME, 101 ) TRIM(Inst%LutDir), 10
1547 : CALL READ_LUT_NCFILE( HcoState, TRIM( FILENAME ), &
1548 : Inst%FRACNOX_LUT10, Inst%DNOx_LUT10, Inst%OPE_LUT10, &
1549 0 : Inst%MOE_LUT10, RC=RC )
1550 :
1551 : ! Read 14 m/s LUT
1552 0 : WRITE( FILENAME, 101 ) TRIM(Inst%LutDir), 14
1553 : CALL READ_LUT_NCFILE( HcoState, TRIM( FILENAME ), &
1554 : Inst%FRACNOX_LUT14, Inst%DNOx_LUT14, Inst%OPE_LUT14, &
1555 0 : Inst%MOE_LUT14, RC=RC )
1556 :
1557 : ! Read 18 m/s LUT
1558 0 : WRITE( FILENAME, 101 ) TRIM(Inst%LutDir), 18
1559 : CALL READ_LUT_NCFILE( HcoState, TRIM( FILENAME ), &
1560 : Inst%FRACNOX_LUT18, Inst%DNOx_LUT18, Inst%OPE_LUT18, &
1561 0 : Inst%MOE_LUT18, RC=RC )
1562 :
1563 : ! ! To write into txt-file, uncomment the following lines
1564 : ! FILENAME = TRIM(LutDir)//'/ship_plume_lut_02ms.txt'
1565 : ! CALL WRITE_LUT_TXTFILE( HcoState, TRIM( FILENAME ), &
1566 : ! FRACNOX_LUT02, DNOx_LUT02, OPE_LUT02, MOE_LUT02, RC, &
1567 : ! Tlev, JNO2lev, O3lev, SEA0lev, SEA5lev, JRATIOlev, NOXlev )
1568 : ! IF ( RC /= HCO_SUCCESS ) RETURN
1569 : !
1570 : ! FILENAME = TRIM(LutDir)//'/ship_plume_lut_06ms.txt'
1571 : ! CALL WRITE_LUT_TXTFILE( HcoState, TRIM( FILENAME ), &
1572 : ! FRACNOX_LUT06, DNOx_LUT06, OPE_LUT06, MOE_LUT06, RC, &
1573 : ! Tlev, JNO2lev, O3lev, SEA0lev, SEA5lev, JRATIOlev, NOXlev )
1574 : ! IF ( RC /= HCO_SUCCESS ) RETURN
1575 : !
1576 : ! FILENAME = TRIM(LutDir)//'/ship_plume_lut_10ms.txt'
1577 : ! CALL WRITE_LUT_TXTFILE( HcoState, TRIM( FILENAME ), &
1578 : ! FRACNOX_LUT10, DNOx_LUT10, OPE_LUT10, MOE_LUT10, RC, &
1579 : ! Tlev, JNO2lev, O3lev, SEA0lev, SEA5lev, JRATIOlev, NOXlev )
1580 : ! IF ( RC /= HCO_SUCCESS ) RETURN
1581 : !
1582 : ! FILENAME = TRIM(LutDir)//'/ship_plume_lut_14ms.txt'
1583 : ! CALL WRITE_LUT_TXTFILE( HcoState, TRIM( FILENAME ), &
1584 : ! FRACNOX_LUT14, DNOx_LUT14, OPE_LUT14, MOE_LUT14, RC, &
1585 : ! Tlev, JNO2lev, O3lev, SEA0lev, SEA5lev, JRATIOlev, NOXlev )
1586 : ! IF ( RC /= HCO_SUCCESS ) RETURN
1587 : !
1588 : ! FILENAME = TRIM(LutDir)//'/ship_plume_lut_14ms.txt'
1589 : ! CALL WRITE_LUT_TXTFILE( HcoState, TRIM( FILENAME ), &
1590 : ! FRACNOX_LUT18, DNOx_LUT18, OPE_LUT18, MOE_LUT18, RC, &
1591 : ! Tlev, JNO2lev, O3lev, SEA0lev, SEA5lev, JRATIOlev, NOXlev )
1592 : ! IF ( RC /= HCO_SUCCESS ) RETURN
1593 :
1594 : ! Return w/ success
1595 0 : RC = HCO_SUCCESS
1596 : #endif
1597 :
1598 0 : END SUBROUTINE READ_PARANOX_LUT_NC
1599 : !EOC
1600 : #if !defined(ESMF_)
1601 : !------------------------------------------------------------------------------
1602 : ! GEOS-Chem Global Chemical Transport Model !
1603 : !------------------------------------------------------------------------------
1604 : !BOP
1605 : !
1606 : ! !IROUTINE: read_lut_ncfile
1607 : !
1608 : ! !DESCRIPTION: Subroutine READ\_LUT\_NCFILE reads look up tables for use in
1609 : ! the PARANOX ship plume model (C. Holmes)
1610 : !\\
1611 : !\\
1612 : ! !INTERFACE:
1613 : !
1614 0 : SUBROUTINE READ_LUT_NCFILE( HcoState, FILENAME, FNOX, &
1615 0 : DNOx, OPE, MOE, T, &
1616 0 : JNO2, O3, SEA0, SEA5, &
1617 0 : JRATIO, NOX, RC )
1618 : !
1619 : ! !USES:
1620 : !
1621 : ! Modules for netCDF read
1622 : USE HCO_m_netcdf_io_open
1623 : USE HCO_m_netcdf_io_get_dimlen
1624 : USE HCO_m_netcdf_io_read
1625 : USE HCO_m_netcdf_io_readattr
1626 : USE HCO_m_netcdf_io_close
1627 :
1628 : # include "netcdf.inc"
1629 : !
1630 : ! !INPUT PARAMETERS:
1631 : !
1632 : TYPE(HCO_State), POINTER :: HcoState ! HEMCO State object
1633 : CHARACTER(LEN=*),INTENT(IN) :: FILENAME
1634 : !
1635 : ! !OUTPUT PARAMETERS:
1636 : !
1637 : REAL*4, INTENT(OUT), DIMENSION(:,:,:,:,:,:,:) :: FNOX,OPE,MOE,DNOx
1638 : REAL*4, INTENT(OUT), OPTIONAL :: T(:), JNO2(:), O3(:)
1639 : REAL*4, INTENT(OUT), OPTIONAL :: SEA0(:), SEA5(:)
1640 : REAL*4, INTENT(OUT), OPTIONAL :: JRATIO(:), NOX(:)
1641 : INTEGER, INTENT(OUT), OPTIONAL :: RC
1642 : !
1643 : ! !REVISION HISTORY:
1644 : ! 06 Feb 2012 - M. Payer - Initial version modified from code provided by
1645 : ! G.C.M. Vinken
1646 : ! See https://github.com/geoschem/hemco for complete history
1647 : !EOP
1648 : !------------------------------------------------------------------------------
1649 : !BOC
1650 : !
1651 : ! !LOCAL VARIABLES:
1652 :
1653 : ! Scalars
1654 : LOGICAL :: FileExists
1655 : INTEGER :: AS, IOS
1656 : INTEGER :: fID, HMRC
1657 :
1658 : ! arrays
1659 : INTEGER :: st1d(1), ct1d(1)
1660 : INTEGER :: st7d(7), ct7d(7)
1661 :
1662 : CHARACTER(LEN=255) :: MSG, FileMsg
1663 :
1664 : !=================================================================
1665 : ! In dry-run mode, print file path to dryrun log and exit.
1666 : ! Otherwise, print file path to the HEMCO log file and continue.
1667 : !=================================================================
1668 :
1669 : ! Test if the file exists
1670 0 : INQUIRE( FILE=TRIM( FileName ), EXIST=FileExists )
1671 :
1672 : ! Create a display string based on whether or not the file is found
1673 0 : IF ( FileExists ) THEN
1674 0 : FileMsg = 'HEMCO (PARANOX): Opening'
1675 : ELSE
1676 0 : FileMsg = 'HEMCO (PARANOX): REQUIRED FILE NOT FOUND'
1677 : ENDIF
1678 :
1679 : ! Print file status to stdout and the HEMCO log
1680 0 : IF ( HcoState%amIRoot ) THEN
1681 0 : WRITE( 6, 300 ) TRIM( FileMsg ), TRIM( FileName )
1682 0 : WRITE( MSG, 300 ) TRIM( FileMsg ), TRIM( FileName )
1683 0 : CALL HCO_MSG( HcoState%Config%Err, MSG )
1684 : 300 FORMAT( a, ' ', a )
1685 : ENDIF
1686 :
1687 : ! For dry-run simulations, return to calling program.
1688 : ! For regular simulations, throw an error if we can't find the file.
1689 0 : IF ( HcoState%Options%IsDryRun ) THEN
1690 0 : RETURN
1691 : ELSE
1692 0 : IF ( .not. FileExists ) THEN
1693 0 : WRITE( MSG, 300 ) TRIM( FileMsg ), TRIM( FileName )
1694 0 : CALL HCO_ERROR(MSG, HMRC )
1695 0 : IF ( PRESENT( RC ) ) RC = HMRC
1696 0 : RETURN
1697 : ENDIF
1698 : ENDIF
1699 :
1700 : !=================================================================
1701 : ! READ_LUT_NCFILE begins here!
1702 : !=================================================================
1703 :
1704 : ! Open file for reading
1705 0 : CALL Ncop_Rd( fId, TRIM(FILENAME) )
1706 :
1707 : !-----------------------------------------------------------------
1708 : ! Read reference values used to construct the Look-up table
1709 : ! These are 1d values
1710 : !-----------------------------------------------------------------
1711 :
1712 : ! Temperature, K
1713 0 : IF ( PRESENT(T) ) THEN
1714 0 : st1d = (/ 1 /)
1715 0 : ct1d = (/ nT /)
1716 0 : CALL NcRd( T, fId, 'T', st1d, ct1d )
1717 : ENDIF
1718 :
1719 : ! J(NO2), 1/s
1720 0 : IF ( PRESENT(JNO2) ) THEN
1721 0 : st1d = (/ 1 /)
1722 0 : ct1d = (/ nJ /)
1723 0 : CALL NcRd( JNO2, fId, 'JNO2', st1d, ct1d )
1724 : ENDIF
1725 :
1726 : ! Ambient O3, ppb
1727 0 : IF ( PRESENT(O3) ) THEN
1728 0 : st1d = (/ 1 /)
1729 0 : ct1d = (/ nO3 /)
1730 0 : CALL NcRd( O3, fId, 'O3', st1d, ct1d )
1731 : ENDIF
1732 :
1733 : ! Solar elevation angle at emission time t=0, deg
1734 0 : IF ( PRESENT(SEA0) ) THEN
1735 0 : st1d = (/ 1 /)
1736 0 : ct1d = (/ nSEA /)
1737 0 : CALL NcRd( SEA0, fId, 'SEA0', st1d, ct1d )
1738 : ENDIF
1739 :
1740 : ! Solar elevation angle at time t=5h, deg
1741 0 : IF ( PRESENT(SEA5) ) THEN
1742 0 : st1d = (/ 1 /)
1743 0 : ct1d = (/ nSEA /)
1744 0 : CALL NcRd( SEA5, fId, 'SEA5', st1d, ct1d )
1745 : ENDIF
1746 :
1747 : ! J(OH) / J(NO2) ratio, s/s
1748 0 : IF ( PRESENT(JRATIO) ) THEN
1749 0 : st1d = (/ 1 /)
1750 0 : ct1d = (/ nJ /)
1751 0 : CALL NcRd( JRatio, fId, 'Jratio', st1d, ct1d )
1752 : ENDIF
1753 :
1754 : ! Ambient NOx, ppt
1755 0 : IF ( PRESENT(NOX) ) THEN
1756 0 : st1d = (/ 1 /)
1757 0 : ct1d = (/ nNOx /)
1758 0 : CALL NcRd( NOX, fId, 'NOx', st1d, ct1d )
1759 : ENDIF
1760 :
1761 : ! Define 7D variables used below
1762 0 : st7d = (/ 1, 1, 1, 1, 1, 1, 1 /)
1763 0 : ct7d = (/ nT,nJ,nO3,nSEA,nSEA,nJ,nNOx /)
1764 :
1765 : !-----------------------------------------------------------------
1766 : ! Read look up table for fraction of NOx remaining for ship
1767 : ! emissions after 5 h [unitless]
1768 : !-----------------------------------------------------------------
1769 :
1770 0 : CALL NcRd( FNOx, fId, 'FNOx', st7d, ct7d )
1771 :
1772 : ! testing only
1773 : ! PRINT*, "binary_fracnox: ", Fnox(1:4,1,1,2,3,4,4)
1774 :
1775 : !-----------------------------------------------------------------
1776 : ! Read look up table for 5-h integrated Ozone Production Efficiency
1777 : ! for ship emissions [molec O3 produced / molec NOx lost]
1778 : !-----------------------------------------------------------------
1779 :
1780 0 : CALL NcRd( OPE, fId, 'OPE', st7d, ct7d )
1781 :
1782 : ! testing only
1783 : ! PRINT*, "binary_intope: ", OPE(1:4,1,1,2,3,4,4)
1784 :
1785 : !-----------------------------------------------------------------
1786 : ! Read look up table for 5-h integrated Methane Oxidation Efficiency
1787 : ! for ship emissions [molec CH4 oxidized / molec NOx emitted]
1788 : !-----------------------------------------------------------------
1789 :
1790 0 : CALL NcRd( MOE, fId, 'MOE', st7d, ct7d )
1791 :
1792 : ! testing only
1793 : ! PRINT*, "binary_intmoe: ", MOE(1:4,1,1,2,3,4,4)
1794 :
1795 : !-----------------------------------------------------------------
1796 : ! Read look up table for 5-h integrated NOx deposition fraction
1797 : ! for ship emissions [molec NOx deposited / molec NOx emitted]
1798 : !-----------------------------------------------------------------
1799 :
1800 0 : CALL NcRd( DNOx, fId, 'DNOx', st7d, ct7d )
1801 :
1802 : ! testing only
1803 : ! PRINT*, "binary_depnox: ", DNOx(1:4,1,1,2,3,4,4)
1804 :
1805 : ! Close netCDF file
1806 0 : CALL NcCl( fId )
1807 :
1808 : END SUBROUTINE READ_LUT_NCFILE
1809 : !EOC
1810 : #endif
1811 : !------------------------------------------------------------------------------
1812 : ! GEOS-Chem Global Chemical Transport Model !
1813 : !------------------------------------------------------------------------------
1814 : !BOP
1815 : !
1816 : ! !IROUTINE: read_paranox_lut_txt
1817 : !
1818 : ! !DESCRIPTION: Subroutine READ\_PARANOX\_LUT\_TXT reads look-up tables in
1819 : ! txt format for use in the PARANOX ship plume model (G.C.M. Vinken)
1820 : !\\
1821 : !\\
1822 : ! !INTERFACE:
1823 : !
1824 0 : SUBROUTINE READ_PARANOX_LUT_TXT ( HcoState, Inst, RC )
1825 : !
1826 : ! !USES:
1827 : !
1828 : ! !INPUT ARGUMENTS:
1829 : !
1830 : TYPE(HCO_State), POINTER :: HcoState ! HEMCO State object
1831 : TYPE(MyInst), POINTER :: Inst
1832 : !
1833 : ! !INPUT/OUTPUT ARGUMENTS:
1834 : !
1835 : INTEGER, INTENT(INOUT) :: RC
1836 : !
1837 : ! !REVISION HISTORY:
1838 : ! 05 Feb 2015 - C. Keller - Initial version modified from code provided by
1839 : ! G.C.M. Vinken and C. Holmes
1840 : ! See https://github.com/geoschem/hemco for complete history
1841 : !EOP
1842 : !------------------------------------------------------------------------------
1843 : !BOC
1844 : !
1845 : ! !LOCAL VARIABLES:
1846 : INTEGER :: IOS
1847 : CHARACTER(LEN=255) :: FILENAME
1848 : CHARACTER(LEN=255) :: MSG, LOC
1849 : INTEGER :: fID
1850 :
1851 : !=================================================================
1852 : ! READ_PARANOX_LUT_TXT begins here
1853 : !=================================================================
1854 0 : LOC = 'READ_PARANOX_LUT_TXT (HCOX_PARANOX_MOD.F90)'
1855 :
1856 : ! Clear FILENAME
1857 0 : FILENAME = ''
1858 :
1859 : ! FILENAME format string
1860 : 101 FORMAT( a, '/ship_plume_lut_', I2.2, 'ms.txt' )
1861 :
1862 : ! Wind speed levels correspond to the files that we will read below
1863 0 : Inst%WSlev = (/ 2.0e0, 6.0e0, 10.0e0, 14.0e0, 18.0e0 /)
1864 :
1865 : ! Read 2m/s LUT
1866 0 : WRITE( FILENAME, 101 ) TRIM(Inst%LutDir), 2
1867 : CALL READ_LUT_TXTFILE( HcoState, TRIM( FILENAME ), &
1868 : Inst%FRACNOX_LUT02, Inst%DNOx_LUT02, Inst%OPE_LUT02, Inst%MOE_LUT02, RC, &
1869 0 : Inst%Tlev, Inst%JNO2lev, Inst%O3lev, Inst%SEA0lev, Inst%SEA5lev, Inst%JRATIOlev, Inst%NOXlev )
1870 0 : IF ( RC /= HCO_SUCCESS ) THEN
1871 0 : CALL HCO_ERROR( 'ERROR 19', RC, THISLOC=LOC )
1872 0 : RETURN
1873 : ENDIF
1874 :
1875 : ! Read 6 m/s LUT
1876 0 : WRITE( FILENAME, 101 ) TRIM(Inst%LutDir), 6
1877 : CALL READ_LUT_TXTFILE( HcoState, TRIM( FILENAME ), &
1878 0 : Inst%FRACNOX_LUT06, Inst%DNOx_LUT06, Inst%OPE_LUT06, Inst%MOE_LUT06, RC )
1879 0 : IF ( RC /= HCO_SUCCESS ) THEN
1880 0 : CALL HCO_ERROR( 'ERROR 20', RC, THISLOC=LOC )
1881 0 : RETURN
1882 : ENDIF
1883 :
1884 : ! Read 10 m/s LUT
1885 0 : WRITE( FILENAME, 101 ) TRIM(Inst%LutDir), 10
1886 : CALL READ_LUT_TXTFILE( HcoState, TRIM( FILENAME ), &
1887 0 : Inst%FRACNOX_LUT10, Inst%DNOx_LUT10, Inst%OPE_LUT10, Inst%MOE_LUT10, RC )
1888 0 : IF ( RC /= HCO_SUCCESS ) THEN
1889 0 : CALL HCO_ERROR( 'ERROR 21', RC, THISLOC=LOC )
1890 0 : RETURN
1891 : ENDIF
1892 :
1893 : ! Read 14 m/s LUT
1894 0 : WRITE( FILENAME, 101 ) TRIM(Inst%LutDir), 14
1895 : CALL READ_LUT_TXTFILE( HcoState, TRIM( FILENAME ), &
1896 0 : Inst%FRACNOX_LUT14, Inst%DNOx_LUT14, Inst%OPE_LUT14, Inst%MOE_LUT14, RC )
1897 0 : IF ( RC /= HCO_SUCCESS ) THEN
1898 0 : CALL HCO_ERROR( 'ERROR 22', RC, THISLOC=LOC )
1899 0 : RETURN
1900 : ENDIF
1901 :
1902 : ! Read 18 m/s LUT
1903 0 : WRITE( FILENAME, 101 ) TRIM(Inst%LutDir), 18
1904 : CALL READ_LUT_TXTFILE( HcoState, TRIM( FILENAME ), &
1905 0 : Inst%FRACNOX_LUT18, Inst%DNOx_LUT18, Inst%OPE_LUT18, Inst%MOE_LUT18, RC )
1906 0 : IF ( RC /= HCO_SUCCESS ) THEN
1907 0 : CALL HCO_ERROR( 'ERROR 23', RC, THISLOC=LOC )
1908 0 : RETURN
1909 : ENDIF
1910 :
1911 : ! Return w/ success
1912 0 : RC = HCO_SUCCESS
1913 :
1914 : END SUBROUTINE READ_PARANOX_LUT_TXT
1915 : !EOC
1916 : !------------------------------------------------------------------------------
1917 : ! GEOS-Chem Global Chemical Transport Model !
1918 : !------------------------------------------------------------------------------
1919 : !BOP
1920 : !
1921 : ! !IROUTINE: read_lut_txtfile
1922 : !
1923 : ! !DESCRIPTION: Subroutine READ\_LUT\_TXTFILE reads look up tables for use in
1924 : ! the PARANOX ship plume model (C. Holmes)
1925 : !\\
1926 : !\\
1927 : ! !INTERFACE:
1928 : !
1929 0 : SUBROUTINE READ_LUT_TXTFILE( HcoState, FILENAME, FNOX, DNOx, OPE, MOE, RC, &
1930 0 : T, JNO2, O3, SEA0, SEA5, JRATIO, NOX )
1931 : !
1932 : ! !USES:
1933 : !
1934 : USE HCO_inquireMod, ONLY : findFreeLUN
1935 : !
1936 : ! !INPUT PARAMETERS:
1937 : !
1938 : TYPE(HCO_State), POINTER :: HcoState ! HEMCO State object
1939 : CHARACTER(LEN=*),INTENT(IN) :: FILENAME
1940 : !
1941 : ! !INPUT/OUTPUT PARAMETERS:
1942 : !
1943 : INTEGER, INTENT(INOUT) :: RC
1944 : !
1945 : ! !OUTPUT PARAMETERS:
1946 : !
1947 : REAL*4, INTENT(OUT), TARGET, DIMENSION(:,:,:,:,:,:,:) :: FNOX,OPE,MOE,DNOx
1948 : REAL*4, INTENT(OUT), OPTIONAL :: T(:), JNO2(:), O3(:)
1949 : REAL*4, INTENT(OUT), OPTIONAL :: SEA0(:), SEA5(:)
1950 : REAL*4, INTENT(OUT), OPTIONAL :: JRATIO(:), NOX(:)
1951 : !
1952 : ! !REVISION HISTORY:
1953 : ! 05 Feb 2015 - C. Keller - Initial version modified from code provided by
1954 : ! G.C.M. Vinken and C. Holmes
1955 : ! See https://github.com/geoschem/hemco for complete history
1956 : !EOP
1957 : !------------------------------------------------------------------------------
1958 : !BOC
1959 : !
1960 : ! !LOCAL VARIABLES:
1961 : !
1962 : ! Scalars
1963 : LOGICAL :: FileExists
1964 : INTEGER :: fId, IOS
1965 :
1966 : ! INTEGER :: I, I1, I2, I3, I4, I5, I6, I7
1967 : ! REAL*4, POINTER :: TMPARR(:,:,:,:,:,:,:) => NULL()
1968 :
1969 : CHARACTER(LEN=255) :: MSG, FileMsg
1970 : !
1971 : ! !DEFINED PARAMETERS:
1972 : !
1973 : CHARACTER(LEN=255), PARAMETER :: FMAT = "(E40.32)"
1974 : CHARACTER(LEN=255), PARAMETER :: LOC = &
1975 : 'READ_LUT_TXTFILE (hcox_paranox_mod.F90)'
1976 :
1977 : !=================================================================
1978 : ! In dry-run mode, print file path to dryrun log and exit.
1979 : ! Otherwise, print file path to the HEMCO log file and continue.
1980 : !=================================================================
1981 :
1982 : ! Test if the file exists
1983 0 : INQUIRE ( FILE=TRIM( FileName ), EXIST=FileExists )
1984 :
1985 : ! Create a display string based on whether or not the file is found
1986 0 : IF ( FileExists ) THEN
1987 0 : FileMsg = 'HEMCO (PARANOX): Opening'
1988 : ELSE
1989 0 : FileMsg = 'HEMCO (PARANOX): REQUIRED FILE NOT FOUND'
1990 : ENDIF
1991 :
1992 : ! Print file status to stdout and the HEMCO log file
1993 0 : IF ( HcoState%amIRoot ) THEN
1994 0 : WRITE( 6, 300 ) TRIM( FileMsg ), TRIM( FileName )
1995 0 : WRITE( MSG, 300 ) TRIM( FileMsg ), TRIM( FileName )
1996 0 : CALL HCO_MSG(HcoState%Config%Err,MSG)
1997 : 300 FORMAT( a, ' ', a )
1998 : ENDIF
1999 :
2000 : ! For dry-run simulations, return to calling program.
2001 : ! For regular simulations, throw an error if we can't find the file.
2002 0 : IF ( HcoState%Options%IsDryRun ) THEN
2003 0 : RETURN
2004 : ELSE
2005 0 : IF ( .not. FileExists ) THEN
2006 0 : WRITE( MSG, 300 ) TRIM( FileMsg ), TRIM( FileName )
2007 0 : CALL HCO_ERROR(MSG, RC )
2008 0 : RETURN
2009 : ENDIF
2010 : ENDIF
2011 :
2012 : !=================================================================
2013 : ! READ_LUT_TXTFILE begins here
2014 : !=================================================================
2015 :
2016 : ! Find a free file LUN
2017 0 : fId = findFreeLUN()
2018 :
2019 : ! Open file for reading
2020 0 : OPEN ( fID, FILE=TRIM(FILENAME), FORM="FORMATTED", IOSTAT=IOS )
2021 0 : IF ( IOS /= 0 ) THEN
2022 0 : CALL HCO_ERROR( 'read_lut_txtfile:1', RC, THISLOC=LOC )
2023 0 : RETURN
2024 : ENDIF
2025 :
2026 : ! Read FNOx
2027 0 : READ( fId, FMT=FMAT, IOSTAT=IOS ) FNOx
2028 0 : IF ( IOS /= 0 ) THEN
2029 0 : CALL HCO_ERROR( 'read_lut_txtfile: FNOx', RC, THISLOC=LOC )
2030 0 : RETURN
2031 : ENDIF
2032 :
2033 : ! Read OPE
2034 0 : READ( fId, FMT=FMAT, IOSTAT=IOS ) OPE
2035 0 : IF ( IOS /= 0 ) THEN
2036 0 : CALL HCO_ERROR( 'read_lut_txtfile: OPE', RC, THISLOC=LOC )
2037 0 : RETURN
2038 : ENDIF
2039 :
2040 : ! Read MOE
2041 0 : READ( fId, FMT=FMAT, IOSTAT=IOS ) MOE
2042 0 : IF ( IOS /= 0 ) THEN
2043 0 : CALL HCO_ERROR( 'read_lut_txtfile: MOE', RC, THISLOC=LOC )
2044 0 : RETURN
2045 : ENDIF
2046 :
2047 : ! Read DNOx
2048 0 : READ( fId, FMT=FMAT, IOSTAT=IOS ) DNOx
2049 0 : IF ( IOS /= 0 ) THEN
2050 0 : CALL HCO_ERROR( 'read_lut_txtfile: DNOx', RC, THISLOC=LOC )
2051 0 : RETURN
2052 : ENDIF
2053 :
2054 : ! Read optional values
2055 0 : IF ( PRESENT(T) ) THEN
2056 0 : READ( fId, FMT=FMAT, IOSTAT=IOS ) T
2057 0 : IF ( IOS /= 0 ) THEN
2058 0 : CALL HCO_ERROR( 'read_lut_txtfile: T', RC, THISLOC=LOC )
2059 0 : RETURN
2060 : ENDIF
2061 : ENDIF
2062 :
2063 0 : IF ( PRESENT(JNO2) ) THEN
2064 0 : READ( fId, FMT=FMAT, IOSTAT=IOS ) JNO2
2065 0 : IF ( IOS /= 0 ) THEN
2066 0 : CALL HCO_ERROR( 'read_lut_txtfile: JNO2', RC, THISLOC=LOC )
2067 0 : RETURN
2068 : ENDIF
2069 : ENDIF
2070 :
2071 0 : IF ( PRESENT(O3) ) THEN
2072 0 : READ( fId, FMT=FMAT, IOSTAT=IOS ) O3
2073 0 : IF ( IOS /= 0 ) THEN
2074 0 : CALL HCO_ERROR( 'read_lut_txtfile: O3', RC, THISLOC=LOC )
2075 0 : RETURN
2076 : ENDIF
2077 : ENDIF
2078 :
2079 0 : IF ( PRESENT(SEA0) ) THEN
2080 0 : READ( fId, FMT=FMAT, IOSTAT=IOS ) SEA0
2081 0 : IF ( IOS /= 0 ) THEN
2082 0 : CALL HCO_ERROR( 'read_lut_txtfile: SEA0', RC, THISLOC=LOC )
2083 0 : RETURN
2084 : ENDIF
2085 : ENDIF
2086 :
2087 0 : IF ( PRESENT(SEA5) ) THEN
2088 0 : READ( fId, FMT=FMAT, IOSTAT=IOS ) SEA5
2089 0 : IF ( IOS /= 0 ) THEN
2090 0 : CALL HCO_ERROR( 'read_lut_txtfile: SEA5', RC, THISLOC=LOC )
2091 0 : RETURN
2092 : ENDIF
2093 : ENDIF
2094 :
2095 0 : IF ( PRESENT(JRATIO) ) THEN
2096 0 : READ( fId, FMT=FMAT, IOSTAT=IOS ) JRATIO
2097 0 : IF ( IOS /= 0 ) THEN
2098 0 : CALL HCO_ERROR( 'read_lut_txtfile: JRATIO', RC, THISLOC=LOC )
2099 0 : RETURN
2100 : ENDIF
2101 : ENDIF
2102 :
2103 0 : IF ( PRESENT(NOX) ) THEN
2104 0 : READ( fId, FMT=FMAT, IOSTAT=IOS ) NOX
2105 0 : IF ( IOS /= 0 ) THEN
2106 0 : CALL HCO_ERROR( 'read_lut_txtfile: NOX', RC, THISLOC=LOC )
2107 0 : RETURN
2108 : ENDIF
2109 : ENDIF
2110 :
2111 : ! ! Read
2112 : ! DO I = 1,4
2113 : !
2114 : ! ! Set pointer to output array
2115 : ! SELECT CASE ( I )
2116 : ! CASE ( 1 )
2117 : ! TMPARR => FNOX
2118 : ! CASE ( 2 )
2119 : ! TMPARR => OPE
2120 : ! CASE ( 3 )
2121 : ! TMPARR => MOE
2122 : ! CASE ( 4 )
2123 : ! TMPARR => DNOx
2124 : ! CASE DEFAULT
2125 : ! CALL HCO_ERROR( 'I > 4', RC, THISLOC=LOC )
2126 : ! RETURN
2127 : ! END SELECT
2128 : !
2129 : ! DO I1 = 1,nT
2130 : ! DO I2 = 1,nJ
2131 : ! DO I3 = 1,nO3
2132 : ! DO I4 = 1,nSEA
2133 : ! DO I5 = 1,nSEA
2134 : ! DO I6 = 1,nJ
2135 : ! DO I7 = 1,nNOx
2136 : ! READ( fId, FMT=FMAT, IOSTAT=IOS ) TMPARR(I1,I2,I3,I4,I5,I6,I7)
2137 : ! ENDDO
2138 : ! ENDDO
2139 : ! ENDDO
2140 : ! ENDDO
2141 : ! ENDDO
2142 : ! ENDDO
2143 : ! ENDDO
2144 : ! ENDDO
2145 :
2146 : ! Close file
2147 0 : CLOSE( fId )
2148 :
2149 : ! Return w/ success
2150 0 : RC = HCO_SUCCESS
2151 :
2152 : END SUBROUTINE READ_LUT_TXTFILE
2153 : !EOC
2154 : !------------------------------------------------------------------------------
2155 : ! GEOS-Chem Global Chemical Transport Model !
2156 : !------------------------------------------------------------------------------
2157 : !BOP
2158 : !
2159 : ! !IROUTINE: write_lut_txtfile
2160 : !
2161 : ! !DESCRIPTION: write\_lut\_txtfile
2162 : !\\
2163 : !\\
2164 : ! !INTERFACE:
2165 : !
2166 : SUBROUTINE WRITE_LUT_TXTFILE( HcoState, FILENAME, FNOX, DNOx, OPE, MOE, RC, &
2167 : T, JNO2, O3, SEA0, SEA5, JRATIO, NOX )
2168 : !
2169 : ! !USES:
2170 : !
2171 : USE HCO_inquireMod, ONLY : findFreeLUN
2172 : !
2173 : ! !INPUT PARAMETERS:
2174 : !
2175 : TYPE(HCO_State), INTENT(INOUT) :: HcoState ! HEMCO state obj
2176 : CHARACTER(LEN=*),INTENT(IN) :: FILENAME
2177 : !
2178 : ! !INPUT/OUTPUT PARAMETERS:
2179 : !
2180 : INTEGER, INTENT(INOUT) :: RC
2181 : !
2182 : ! !OUTPUT PARAMETERS:
2183 : !
2184 : REAL*4, INTENT(IN), TARGET, DIMENSION(:,:,:,:,:,:,:) :: FNOX,OPE,MOE,DNOx
2185 : REAL*4, INTENT(IN), OPTIONAL :: T(:), JNO2(:), O3(:)
2186 : REAL*4, INTENT(IN), OPTIONAL :: SEA0(:), SEA5(:)
2187 : REAL*4, INTENT(IN), OPTIONAL :: JRATIO(:), NOX(:)
2188 : !
2189 : ! !REVISION HISTORY:
2190 : ! 05 Feb 2015 - C. Keller - Initial version modified from code provided by
2191 : ! G.C.M. Vinken and C. Holmes
2192 : ! See https://github.com/geoschem/hemco for complete history
2193 : !EOP
2194 : !------------------------------------------------------------------------------
2195 : !BOC
2196 : !
2197 : ! !LOCAL VARIABLES:
2198 : INTEGER :: fId, IOS
2199 :
2200 : ! INTEGER :: I, I1, I2, I3, I4, I5, I6, I7
2201 : ! REAL*4, POINTER :: TMPARR(:,:,:,:,:,:,:) => NULL()
2202 :
2203 : CHARACTER(LEN=255) :: MSG
2204 : CHARACTER(LEN=255), PARAMETER :: FMAT = "(E40.32)"
2205 : CHARACTER(LEN=255), PARAMETER :: LOC = &
2206 : 'WRITE_LUT_TXTFILE (hcox_paranox_mod.F90)'
2207 :
2208 : !=================================================================
2209 : ! WRITE_LUT_TXTFILE begins here
2210 : !=================================================================
2211 :
2212 : ! Echo info
2213 : IF ( Hcostate%amIRoot ) THEN
2214 : WRITE( MSG, 100 ) TRIM( FILENAME )
2215 : CALL HCO_MSG(HcoState%Config%Err,MSG)
2216 : 100 FORMAT( 'WRITE_LUT_TXTFILE: Writing ', a )
2217 : ENDIF
2218 :
2219 : ! Find a free file LUN
2220 : fId = findFreeLUN()
2221 :
2222 : ! Open file for reading
2223 : OPEN ( fID, FILE=TRIM(FILENAME), ACTION="WRITE", FORM="FORMATTED", IOSTAT=IOS )
2224 : IF ( IOS /= 0 ) THEN
2225 : CALL HCO_ERROR( 'write_lut_txtfile:1', RC, THISLOC=LOC )
2226 : RETURN
2227 : ENDIF
2228 :
2229 : ! Read FNOx
2230 : WRITE( fId, FMT=FMAT, IOSTAT=IOS ) FNOx
2231 : IF ( IOS /= 0 ) THEN
2232 : CALL HCO_ERROR( 'write_lut_txtfile: FNOx', RC, THISLOC=LOC )
2233 : RETURN
2234 : ENDIF
2235 :
2236 : ! Read OPE
2237 : WRITE( fId, FMT=FMAT, IOSTAT=IOS ) OPE
2238 : IF ( IOS /= 0 ) THEN
2239 : CALL HCO_ERROR( 'read_lut_txtfile: OPE', RC, THISLOC=LOC )
2240 : RETURN
2241 : ENDIF
2242 :
2243 : ! Read MOE
2244 : WRITE( fId, FMT=FMAT, IOSTAT=IOS ) MOE
2245 : IF ( IOS /= 0 ) THEN
2246 : CALL HCO_ERROR( 'read_lut_txtfile: MOE', RC, THISLOC=LOC )
2247 : RETURN
2248 : ENDIF
2249 :
2250 : ! Read DNOx
2251 : WRITE( fId, FMT=FMAT, IOSTAT=IOS ) DNOx
2252 : IF ( IOS /= 0 ) THEN
2253 : CALL HCO_ERROR( 'read_lut_txtfile: DNOx', RC, THISLOC=LOC )
2254 : RETURN
2255 : ENDIF
2256 :
2257 : ! Read optional values
2258 : IF ( PRESENT(T) ) THEN
2259 : WRITE( fId, FMT=FMAT, IOSTAT=IOS ) T
2260 : IF ( IOS /= 0 ) THEN
2261 : CALL HCO_ERROR( 'read_lut_txtfile: T', RC, THISLOC=LOC )
2262 : RETURN
2263 : ENDIF
2264 : ENDIF
2265 :
2266 : IF ( PRESENT(JNO2) ) THEN
2267 : WRITE( fId, FMT=FMAT, IOSTAT=IOS ) JNO2
2268 : IF ( IOS /= 0 ) THEN
2269 : CALL HCO_ERROR( 'read_lut_txtfile: JNO2', RC, THISLOC=LOC )
2270 : RETURN
2271 : ENDIF
2272 : ENDIF
2273 :
2274 : IF ( PRESENT(O3) ) THEN
2275 : WRITE( fId, FMT=FMAT, IOSTAT=IOS ) O3
2276 : IF ( IOS /= 0 ) THEN
2277 : CALL HCO_ERROR( 'read_lut_txtfile: O3', RC, THISLOC=LOC )
2278 : RETURN
2279 : ENDIF
2280 : ENDIF
2281 :
2282 : IF ( PRESENT(SEA0) ) THEN
2283 : WRITE( fId, FMT=FMAT, IOSTAT=IOS ) SEA0
2284 : IF ( IOS /= 0 ) THEN
2285 : CALL HCO_ERROR( 'read_lut_txtfile: SEA0', RC, THISLOC=LOC )
2286 : RETURN
2287 : ENDIF
2288 : ENDIF
2289 :
2290 : IF ( PRESENT(SEA5) ) THEN
2291 : WRITE( fId, FMT=FMAT, IOSTAT=IOS ) SEA5
2292 : IF ( IOS /= 0 ) THEN
2293 : CALL HCO_ERROR( 'read_lut_txtfile: SEA5', RC, THISLOC=LOC )
2294 : RETURN
2295 : ENDIF
2296 : ENDIF
2297 :
2298 : IF ( PRESENT(JRATIO) ) THEN
2299 : WRITE( fId, FMT=FMAT, IOSTAT=IOS ) JRATIO
2300 : IF ( IOS /= 0 ) THEN
2301 : CALL HCO_ERROR( 'read_lut_txtfile: JRATIO', RC, THISLOC=LOC )
2302 : RETURN
2303 : ENDIF
2304 : ENDIF
2305 :
2306 : IF ( PRESENT(NOX) ) THEN
2307 : WRITE( fId, FMT=FMAT, IOSTAT=IOS ) NOX
2308 : IF ( IOS /= 0 ) THEN
2309 : CALL HCO_ERROR( 'read_lut_txtfile: NOX', RC, THISLOC=LOC )
2310 : RETURN
2311 : ENDIF
2312 : ENDIF
2313 :
2314 : ! Close file
2315 : CLOSE( fId )
2316 :
2317 : ! Return w/ success
2318 : RC = HCO_SUCCESS
2319 :
2320 : END SUBROUTINE WRITE_LUT_TXTFILE
2321 : !EOC
2322 : !------------------------------------------------------------------------------
2323 : ! GEOS-Chem Global Chemical Transport Model !
2324 : !------------------------------------------------------------------------------
2325 : !BOP
2326 : !
2327 : ! !IROUTINE: INTERPOL_LINWEIGHTS
2328 : !
2329 : ! !DESCRIPTION: Subroutine INTERPOL\_LINWEIGHTS finds the array elements and
2330 : ! weights for piecewise 1-D linear interpolation. The input array of NODES
2331 : ! must be in monotonic ascending order. (C. Holmes 3/27/2014)
2332 : !
2333 : ! If Y is an array containing values of a function evaluated at the points
2334 : ! given in NODES, then its interpolated value at the point VALUESIN will be
2335 : ! Y(VALUEIN) = Y(INDICES(1)) * WEIGHTS(1) +
2336 : ! Y(INDICES(2)) * WEIGHTS(2)
2337 : !
2338 : ! This subroutine finds indices of consecutive nodes that bracket VALUEIN and
2339 : ! weights such that
2340 : ! VALUEIN = NODES(INDICES(1)) * WEIGHTS(1) +
2341 : ! NODES(INDICES(1)+1) * (1-WEIGHTS(1))
2342 : !
2343 : ! For convenience, the returned values of INDICES and WEIGHTS are 2-element
2344 : ! arrays, where
2345 : ! INDICES(2) = INDICES(1)+1 and
2346 : ! WEIGHTS(2) = 1 - WEIGHTS(1)
2347 : !\\
2348 : !\\
2349 : ! !INTERFACE:
2350 : !
2351 0 : SUBROUTINE INTERPOL_LINWEIGHTS( NODES, VALUEIN, INDICES, WEIGHTS )
2352 : !
2353 : ! !USES:
2354 : !
2355 : !
2356 : ! !INPUT PARAMETERS:
2357 : !
2358 : REAL*4,INTENT(IN) :: NODES(:), VALUEIN
2359 : !
2360 : ! !OUTPUT PARAMETERS:
2361 : !
2362 : ! These arrays are always 2 elements each, but declaring
2363 : ! as deferred shape avoids array temporaries
2364 : INTEGER,INTENT(OUT) :: INDICES(:)
2365 : REAL*4, INTENT(OUT) :: WEIGHTS(:)
2366 : !
2367 : ! !REVISION HISTORY:
2368 : ! 03 Jun 2013 - C. Holmes - Initial version
2369 : ! See https://github.com/geoschem/hemco for complete history
2370 : !EOP
2371 : !------------------------------------------------------------------------------
2372 : !BOC
2373 : !
2374 : ! !LOCAL VARIABLES:
2375 : INTEGER :: I
2376 : REAL*8 :: VALUE
2377 :
2378 : !=================================================================
2379 : ! INTERPOL_LINWEIGHTS begins here!
2380 : !=================================================================
2381 :
2382 : ! If larger than largest in LUT, assign largest level values
2383 0 : VALUE = MIN( VALUEIN, MAXVAL( NODES ) )
2384 :
2385 : ! If smaller, assign smallest level value
2386 : !GanLuo+VALUE = MAX( VALUE, MINVAL( NODES ) )
2387 0 : VALUE = MAX( VALUE, MINVAL( NODES )*1.d0 )
2388 :
2389 : ! Initialize
2390 0 : INDICES = (/ 1, 1 /)
2391 :
2392 : ! Loop over interpolation nodes until we find the largest node value
2393 : ! that is less than the desired value
2394 0 : DO I=1, SIZE(NODES)
2395 0 : INDICES(1) = I
2396 0 : IF ( VALUE <= NODES(I+1) ) EXIT
2397 : END DO
2398 :
2399 : ! The next node
2400 0 : INDICES(2) = INDICES(1) + 1
2401 :
2402 : ! Weights for the corresponding node indices
2403 0 : WEIGHTS(1) = ( NODES(I+1) - VALUE ) / ( NODES(I+1) - NODES(I) )
2404 0 : WEIGHTS(2) = 1.0 - WEIGHTS(1)
2405 :
2406 0 : END SUBROUTINE INTERPOL_LINWEIGHTS
2407 : !EOC
2408 : !------------------------------------------------------------------------------
2409 : ! GEOS-Chem Global Chemical Transport Model !
2410 : !------------------------------------------------------------------------------
2411 : !BOP
2412 : !
2413 : ! !IROUTINE: paranox_lut
2414 : !
2415 : ! !DESCRIPTION: Subroutine PARANOX\_LUT returns fractional remainder of
2416 : ! ship NOx (FNOx), fraction of NOx that dry deposits as NOy species (DNOx),
2417 : ! ozone production efficiency (OPE), and methane oxidation
2418 : ! efficiency (MOE) after 5-hrs of plume aging. Values are taken taken from a
2419 : ! lookup table using piecewise linear interpolation. The look-up table is derived
2420 : ! from the PARANOx gaussian plume model (Vinken et al. 2011; Holmes et al. 2014)
2421 : ! (G.C.M. Vinken, KNMI, June 2010; C. Holmes June 2013)
2422 : !
2423 : ! The lookup table uses 8 input variables:
2424 : ! TEMP : model temperature, K
2425 : ! JNO2 : J(NO2) value, 1/s
2426 : ! O3 : concentration O3 in ambient air, ppb
2427 : ! SEA0 : solar elevation angle at emission time 5 hours ago, degree
2428 : ! SEA5 : solar elevation angle at this time, degree
2429 : ! JRatio : ratio J(OH)/J(NO2), unitless
2430 : ! NOx : concentration NOx in ambient air, ppt
2431 : ! WS : wind speed, m/s
2432 : !
2433 : ! In GEOS-Chem v9-01-03 through v9-02, the effects of wind speed on FNOx and OPE
2434 : ! were not included (wind speed set at 6 m/s). The JRatio also used J(O1D)
2435 : ! rather than J(OH); this has only a small effect on interpolated values.
2436 : ! To reproduce the behavior of these earlier versions, modify code below marked
2437 : ! with ******* and call READ\_PARANOX\_LUT\_v913 in emissions\_mod.F
2438 : !\\
2439 : !\\
2440 : ! !INTERFACE:
2441 : !
2442 0 : SUBROUTINE PARANOX_LUT( ExtState, HcoState, Inst, &
2443 : I, J, RC, FNOX, DNOx, OPE, MOE_OUT )
2444 : !
2445 : ! !USES:
2446 : !
2447 : USE HCO_STATE_MOD, ONLY : HCO_State
2448 : USE HCOX_STATE_MOD, ONLY : Ext_State
2449 : !
2450 : ! !INPUT PARAMETERS:
2451 : !
2452 : TYPE(Ext_State), POINTER :: ExtState
2453 : TYPE(HCO_State), POINTER :: HcoState
2454 : TYPE(MyInst), POINTER :: Inst
2455 : INTEGER, INTENT(IN) :: I, J ! Grid indices
2456 : !
2457 : ! !OUTPUT PARAMETERS:
2458 : !
2459 : REAL*8, INTENT(OUT) :: FNOX ! fraction of NOx remaining, mol/mol
2460 : REAL*8, INTENT(OUT) :: DNOX ! fraction of NOx deposited, mol/mol
2461 : REAL*8, INTENT(OUT) :: OPE ! net OPE, mol(net P(O3))/mol(P(HNO3))
2462 : REAL*8, INTENT(OUT), OPTIONAL :: MOE_OUT ! net MOE, mol(L(CH4))/mol(E(NOx))
2463 : !
2464 : ! !INPUT/OUTPUT PARAMETERS:
2465 : !
2466 : INTEGER, INTENT(INOUT) :: RC ! Return code
2467 : !
2468 : ! !REVISION HISTORY:
2469 : ! Jun 2010 - G.C.M. Vinken - Initial version
2470 : ! See https://github.com/geoschem/hemco for complete history
2471 : !EOP
2472 : !------------------------------------------------------------------------------
2473 : !BOC
2474 : !
2475 : ! !LOCAL VARIABLES:
2476 : !
2477 : INTEGER :: I1,I2,I3,I4,I5,I6,I7,I8
2478 : REAL(sp) :: FNOX_TMP, DNOX_TMP, OPE_TMP, MOE_TMP
2479 : REAL(sp) :: WEIGHT
2480 : REAL(sp) :: JNO2, JOH, TAIR
2481 : REAL(sp) :: AIR
2482 : REAL*8 :: MOE
2483 :
2484 : ! Interpolation variables, indices, and weights
2485 : REAL(sp), DIMENSION(8) :: VARS
2486 : INTEGER, DIMENSION(8,2) :: INDX
2487 : REAL(sp), DIMENSION(8,2) :: WTS
2488 :
2489 0 : REAL(sp), POINTER :: FRACNOX_LUT(:,:,:,:,:,:,:)
2490 0 : REAL(sp), POINTER :: DNOX_LUT (:,:,:,:,:,:,:)
2491 0 : REAL(sp), POINTER :: OPE_LUT (:,:,:,:,:,:,:)
2492 0 : REAL(sp), POINTER :: MOE_LUT (:,:,:,:,:,:,:)
2493 :
2494 : CHARACTER(LEN=255) :: MSG
2495 : CHARACTER(LEN=255) :: LOC = 'PARANOX_LUT'
2496 :
2497 : !=================================================================
2498 : ! PARANOX_LUT begins here!
2499 : !=================================================================
2500 :
2501 : ! Initialize for safety's sake
2502 0 : FRACNOX_LUT => NULL()
2503 0 : DNOX_LUT => NULL()
2504 0 : OPE_LUT => NULL()
2505 0 : MOE_LUT => NULL()
2506 0 : FNOX = 0.0_sp
2507 0 : DNOX = 0.0_sp
2508 0 : OPE = 0.0_sp
2509 0 : IF ( PRESENT( MOE_OUT ) ) THEN
2510 0 : MOE_OUT = 0.0_sp
2511 : ENDIF
2512 :
2513 : ! Air mass [kg]
2514 0 : AIR = ExtState%AIR%Arr%Val(I,J,1)
2515 :
2516 : ! Air temperature, K
2517 0 : Tair = ExtState%T2M%Arr%Val(I,J)
2518 :
2519 : ! ! for debugging only
2520 : ! if(I==3.and.J==35)then
2521 : ! write(*,*) 'Call PARANOX_LUT @ ',I,J
2522 : ! write(*,*) 'Tair: ', Tair
2523 : ! write(*,*) 'SUNCOSmid: ', SC5(I,J)
2524 : ! endif
2525 :
2526 : ! Check if sun is up
2527 0 : IF ( ExtState%SUNCOS%Arr%Val(I,J) > 0.0_hp ) THEN
2528 :
2529 : ! J(NO2), 1/s
2530 0 : JNO2 = ExtState%JNO2%Arr%Val(I,J)
2531 :
2532 : ! ! J(O1D), 1/s
2533 : ! JO1D = ExtState%JO1D%Arr%Val(I,J)
2534 : !
2535 : ! ! H2O, molec/cm3. Get from specific humidity, which is in kg/kg.
2536 : ! ! NOTE: SPHU is mass H2O / mass total air so use of dry air molecular
2537 : ! ! weight is slightly inaccurate. C (ewl, 9/11/15)
2538 : ! H2O = ExtState%SPHU%Arr%Val(I,J,1) * DENS &
2539 : ! * HcoState%Phys%AIRMW / MWH2O
2540 : !
2541 : ! ! Calculate J(OH), the effective rate for O3+hv -> OH+OH,
2542 : ! ! assuming steady state for O(1D).
2543 : ! ! Rate coefficients are cm3/molec/s; concentrations are molec/cm3
2544 : ! ! This should match the O3+hv (+H2O) -> OH+OH kinetics in calcrate.F
2545 : ! JOH = JO1D * &
2546 : ! 1.63e-10 * EXP( 60.e0/Tair) * H2O / &
2547 : ! ( 1.63e-10 * EXP( 60.e0/Tair) * H2O + &
2548 : ! 1.20e-10 * DENS * 0.5000e-6 + &
2549 : ! 2.15e-11 * EXP(110.e0/Tair) * DENS * 0.7808e0 + &
2550 : ! 3.30e-11 * EXP( 55.e0/Tair) * DENS * 0.2095e0 )
2551 :
2552 : ! J(OH) - effective rate for O3+hv(+H2O)-> OH+OH, 1/s
2553 0 : JOH = ExtState%JOH%Arr%Val(I,J)
2554 :
2555 : ELSE
2556 :
2557 : ! J-values are zero when sun is down
2558 : JNO2 = 0e0_sp
2559 : JOH = 0e0_sp
2560 :
2561 : ENDIF
2562 :
2563 : ! ! for debugging only
2564 : ! if(I==3.and.J==35)then
2565 : ! write(*,*) 'JNO2: ', JNO2
2566 : ! write(*,*) 'JOH : ', JOH
2567 : ! endif
2568 :
2569 : !========================================================================
2570 : ! Load all variables into a single array
2571 : !========================================================================
2572 :
2573 : ! Temperature, K
2574 0 : VARS(1) = Tair
2575 :
2576 : ! J(NO2), 1/s
2577 0 : VARS(2) = JNO2
2578 :
2579 : ! old
2580 : ! ! O3 concentration in ambient air, ppb
2581 : ! VARS(3) = ExtState%O3%Arr%Val(I,J,1) / AIR &
2582 : ! * HcoState%Phys%AIRMW / MW_O3 &
2583 : ! * 1.e9_sp
2584 : ! new
2585 : ! O3 concentration in ambient air, ppb
2586 : ! NOTE: ExtState%O3 units are now kg/kg dry air (ewl, 9/11/15)
2587 0 : VARS(3) = ExtState%O3%Arr%Val(I,J,1) &
2588 : * HcoState%Phys%AIRMW / Inst%MW_O3 &
2589 0 : * 1.e9_sp
2590 : ! end new (ewl)
2591 :
2592 : ! Solar elevation angle, degree
2593 : ! SEA0 = SEA when emitted from ship, 5-h before current model time
2594 : ! SEA5 = SEA at current model time, 5-h after emission from ship
2595 : ! Note: Since SEA = 90 - SZA, then cos(SZA) = sin(SEA) and
2596 : ! thus SEA = arcsin( cos( SZA ) )
2597 : !VARS(4) = ASIND( SC5(I,J) )
2598 : !VARS(5) = ASIND( ExtState%SUNCOS%Arr%Val(I,J) )
2599 0 : VARS(4) = ASIN( Inst%SC5(I,J) ) / HcoState%Phys%PI_180
2600 0 : VARS(5) = ASIN( ExtState%SUNCOS%Arr%Val(I,J) ) / HcoState%Phys%PI_180
2601 :
2602 : ! J(OH)/J(NO2), unitless
2603 : ! Note J(OH) is the loss rate (1/s) of O3 to OH, which accounts for
2604 : ! the temperature, pressure and water vapor dependence of these reactions
2605 0 : VARS(6) = 0.0_sp
2606 0 : IF ( JNO2 /= 0.0_sp ) THEN
2607 0 : IF ( (EXPONENT(JOH)-EXPONENT(JNO2)) < MAXEXPONENT(JOH) ) THEN
2608 0 : VARS(6) = JOH / JNO2
2609 : ENDIF
2610 : ENDIF
2611 :
2612 : ! old
2613 : ! ! NOx concetration in ambient air, ppt
2614 : ! VARS(7) = ( ( ExtState%NO%Arr%Val(I,J,1) / AIR &
2615 : ! * HcoState%Phys%AIRMW / MW_NO ) &
2616 : ! + ( ExtState%NO2%Arr%Val(I,J,1) / AIR &
2617 : ! * HcoState%Phys%AIRMW / MW_NO2 ) ) &
2618 : ! * 1.e12_sp
2619 : ! new
2620 : ! NOx concetration in ambient air, ppt
2621 : ! NOTE: ExtState vars NO and NO2 units are now kg/kg dry air (ewl, 9/11/15)
2622 0 : VARS(7) = ( ( ExtState%NO%Arr%Val(I,J,1) &
2623 : * HcoState%Phys%AIRMW / Inst%MW_NO ) &
2624 0 : + ( ExtState%NO2%Arr%Val(I,J,1) &
2625 : * HcoState%Phys%AIRMW / Inst%MW_NO2 ) ) &
2626 0 : * 1.e12_sp
2627 : ! end new (ewl)
2628 :
2629 : ! Wind speed, m/s
2630 0 : VARS(8) = SQRT( ExtState%U10M%Arr%Val(I,J)**2 &
2631 0 : + ExtState%V10M%Arr%Val(I,J)**2 )
2632 :
2633 : ! ! for debugging only
2634 : ! if(I==1.and.J==35)then
2635 : ! write(*,*) 'VARS(1) : ', VARS(1)
2636 : ! write(*,*) 'VARS(2) : ', VARS(2)
2637 : ! write(*,*) 'VARS(3) : ', VARS(3)
2638 : ! write(*,*) 'VARS(4) : ', VARS(4)
2639 : ! write(*,*) 'VARS(5) : ', VARS(5)
2640 : ! write(*,*) 'VARS(6) : ', VARS(6)
2641 : ! write(*,*) 'VARS(7) : ', VARS(7)
2642 : ! write(*,*) 'VARS(8) : ', VARS(8)
2643 : ! write(*,*) 'AIR : ', ExtState%AIR%Arr%Val(I,J,1)
2644 : ! write(*,*) 'AIRMW : ', HcoState%Phys%AIRMW
2645 : ! write(*,*) 'MWNO, NO2, O3: ', Inst%MW_NO, Inst%MW_NO2, Inst%MW_O3
2646 : ! write(*,*) 'O3conc : ', ExtState%O3%Arr%Val(I,J,1)
2647 : ! write(*,*) 'NO,NO2 conc : ', ExtState%NO%Arr%Val(I,J,1), &
2648 : ! ExtState%NO2%Arr%Val(I,J,1)
2649 : ! write(*,*) 'U, V : ', ExtState%U10M%Arr%Val(I,J), &
2650 : ! ExtState%V10M%Arr%Val(I,J)
2651 : ! endif
2652 :
2653 : !*****************************************************
2654 : ! Restoring the following lines reproduces the behavior of
2655 : ! GEOS-Chem v9-01-03 through v9-02
2656 : ! the LUT was indexed with the ratio J(O1D)/J(NO2) (cdh, 3/27/2014)
2657 : ! VARS(6) = SAFE_DIV( JO1D, JNO2, 0D0 )
2658 : ! JRATIOlev = (/ 5.e-4, 0.0015, 0.0025, 0.0055 /)
2659 : !*****************************************************
2660 :
2661 : !========================================================================
2662 : ! Find the indices of nodes and their corresponding weights for the
2663 : ! interpolation
2664 : !========================================================================
2665 :
2666 : ! Temperature:
2667 0 : CALL INTERPOL_LINWEIGHTS( Inst%Tlev, VARS(1), INDX(1,:), WTS(1,:) )
2668 :
2669 : ! J(NO2):
2670 0 : CALL INTERPOL_LINWEIGHTS( Inst%JNO2lev, VARS(2), INDX(2,:), WTS(2,:) )
2671 :
2672 : ! [O3]:
2673 0 : CALL INTERPOL_LINWEIGHTS( Inst%O3lev, VARS(3), INDX(3,:), WTS(3,:) )
2674 :
2675 : ! SEA0:
2676 0 : CALL INTERPOL_LINWEIGHTS( Inst%SEA0lev, VARS(4), INDX(4,:), WTS(4,:) )
2677 :
2678 : ! SEA5:
2679 0 : CALL INTERPOL_LINWEIGHTS( Inst%SEA5lev, VARS(5), INDX(5,:), WTS(5,:) )
2680 :
2681 : ! JRATIO:
2682 0 : CALL INTERPOL_LINWEIGHTS( Inst%JRATIOlev, VARS(6), INDX(6,:), WTS(6,:))
2683 :
2684 : ! [NOx]:
2685 0 : CALL INTERPOL_LINWEIGHTS( Inst%NOXlev, VARS(7), INDX(7,:), WTS(7,:) )
2686 :
2687 : ! Wind speed:
2688 0 : CALL INTERPOL_LINWEIGHTS( Inst%WSlev, VARS(8), INDX(8,:), WTS(8,:) )
2689 :
2690 : !========================================================================
2691 : ! Piecewise linear interpolation
2692 : !========================================================================
2693 :
2694 : ! Initialize
2695 0 : FNOX = 0.0d0
2696 0 : DNOx = 0.0d0
2697 0 : OPE = 0.0d0
2698 0 : MOE = 0.0d0
2699 :
2700 : ! Loop over wind speed
2701 0 : DO I8=1,2
2702 :
2703 : ! Point at the LUT for this wind speed
2704 : ! Last two digits in fortran variable names indicate wind speed in m/s
2705 0 : SELECT CASE ( NINT( Inst%WSlev(INDX(8,I8)) ) )
2706 : CASE ( 2 )
2707 0 : FRACNOX_LUT => Inst%FRACNOX_LUT02
2708 0 : DNOx_LUT => Inst%DNOx_LUT02
2709 0 : OPE_LUT => Inst%OPE_LUT02
2710 0 : MOE_LUT => Inst%MOE_LUT02
2711 : CASE ( 6 )
2712 0 : FRACNOX_LUT => Inst%FRACNOX_LUT06
2713 0 : DNOx_LUT => Inst%DNOx_LUT06
2714 0 : OPE_LUT => Inst%OPE_LUT06
2715 0 : MOE_LUT => Inst%MOE_LUT06
2716 : CASE ( 10 )
2717 0 : FRACNOX_LUT => Inst%FRACNOX_LUT10
2718 0 : DNOx_LUT => Inst%DNOx_LUT10
2719 0 : OPE_LUT => Inst%OPE_LUT10
2720 0 : MOE_LUT => Inst%MOE_LUT10
2721 : CASE ( 14 )
2722 0 : FRACNOX_LUT => Inst%FRACNOX_LUT14
2723 0 : DNOx_LUT => Inst%DNOx_LUT14
2724 0 : OPE_LUT => Inst%OPE_LUT14
2725 0 : MOE_LUT => Inst%MOE_LUT14
2726 : CASE ( 18 )
2727 0 : FRACNOX_LUT => Inst%FRACNOX_LUT18
2728 0 : DNOx_LUT => Inst%DNOx_LUT18
2729 0 : OPE_LUT => Inst%OPE_LUT18
2730 0 : MOE_LUT => Inst%MOE_LUT18
2731 : CASE DEFAULT
2732 0 : MSG = 'LUT error: Wind speed interpolation error!'
2733 0 : CALL HCO_ERROR(MSG, RC, THISLOC=LOC )
2734 0 : RETURN
2735 : END SELECT
2736 :
2737 : !*****************************************************
2738 : ! Restoring the following lines reproduces the behavior of
2739 : ! GEOS-Chem v9-01-03 through v9-02 in which wind speed
2740 : ! effects on FNOx and OPE are neglected (cdh, 3/25/2014)
2741 : ! FRACNOX_LUT => FRACNOX_LUT_v913
2742 : ! OPE_LUT => OPE_LUT_v913
2743 : !*****************************************************
2744 :
2745 : ! loop over all other variables
2746 0 : DO I7=1,2
2747 0 : DO I6=1,2
2748 0 : DO I5=1,2
2749 0 : DO I4=1,2
2750 0 : DO I3=1,2
2751 0 : DO I2=1,2
2752 0 : DO I1=1,2
2753 :
2754 : !------------------------------------------
2755 : ! Nodes and weights used in the interpolation
2756 : !------------------------------------------
2757 :
2758 : ! Fraction NOx from the LUT
2759 0 : FNOX_TMP = FRACNOX_LUT( INDX(1,I1), INDX(2,I2), INDX(3,I3), &
2760 0 : INDX(4,I4), INDX(5,I5), &
2761 0 : INDX(6,I6), INDX(7,I7) )
2762 :
2763 0 : DNOX_TMP = DNOx_LUT( INDX(1,I1), INDX(2,I2), INDX(3,I3), &
2764 : INDX(4,I4), INDX(5,I5), &
2765 0 : INDX(6,I6), INDX(7,I7) )
2766 :
2767 : ! OPE from the LUT
2768 0 : OPE_TMP = OPE_LUT( INDX(1,I1), INDX(2,I2), INDX(3,I3), &
2769 : INDX(4,I4), INDX(5,I5), &
2770 0 : INDX(6,I6), INDX(7,I7) )
2771 :
2772 : ! MOE from the LUT
2773 0 : MOE_TMP = MOE_LUT( INDX(1,I1), INDX(2,I2), INDX(3,I3), &
2774 : INDX(4,I4), INDX(5,I5), &
2775 0 : INDX(6,I6), INDX(7,I7) )
2776 :
2777 : ! Interpolation weight for this element
2778 : WEIGHT = WTS(1,I1) * WTS(2,I2) * WTS(3,I3) * WTS(4,I4) * &
2779 0 : WTS(5,I5) * WTS(6,I6) * WTS(7,I7) * WTS(8,I8)
2780 :
2781 : !-----------------------------------
2782 : ! Error Check
2783 : !-----------------------------------
2784 :
2785 : !IF ENCOUNTER -999 IN THE LUT PRINT ERROR!!
2786 0 : IF ( ( FNOX_TMP < 0. ) .or. ( FNOX_TMP > 1. ) ) THEN
2787 :
2788 0 : PRINT*, 'PARANOX_LUT: fracnox = ', FNOX_TMP
2789 0 : PRINT*, 'This occured at grid box ', I, J
2790 0 : PRINT*, 'Lon/Lat: ', HcoState%Grid%XMID%Val(I,J), HcoState%Grid%YMID%Val(I,J)
2791 0 : PRINT*, 'SZA 5 hours ago : ', VARS(4)
2792 0 : PRINT*, 'SZA at this time: ', VARS(5)
2793 0 : PRINT*, 'The two SZAs should not be more than 75 deg apart!'
2794 0 : PRINT*, 'If they are, your restart SZA might be wrong.'
2795 0 : PRINT*, 'You can try to coldstart PARANOx by commenting'
2796 0 : PRINT*, 'all PARANOX_SUNCOS entries in your HEMCO'
2797 0 : PRINT*, 'configuration file.'
2798 :
2799 : !print*, I1, I2, I3, I4, I5, I6, I7, I8
2800 : !print*, INDX(1,I1), INDX(2,I2), INDX(3,I3), INDX(4,I4), &
2801 : ! INDX(5,I5), INDX(6,I6), INDX(7,I7), INDX(8,I8)
2802 : !print*, VARS
2803 :
2804 0 : MSG = 'LUT error: Fracnox should be between 0 and 1!'
2805 0 : CALL HCO_ERROR(MSG, RC, THISLOC=LOC )
2806 0 : RETURN
2807 : ENDIF
2808 :
2809 : !-----------------------------------
2810 : ! Final interpolated values
2811 : !-----------------------------------
2812 :
2813 : ! Weighted sum of FNOx from the LUT
2814 0 : FNOx = FNOx + FNOX_TMP * WEIGHT
2815 :
2816 : ! Weighted sum of DNOx from the LUT
2817 0 : DNOx = DNOx + DNOX_TMP * WEIGHT
2818 :
2819 : ! Weighted sum of OPE from the LUT
2820 0 : OPE = OPE + OPE_TMP * WEIGHT
2821 :
2822 : ! Weighted sum of MOE from the LUT
2823 0 : MOE = MOE + MOE_TMP * WEIGHT
2824 :
2825 : END DO
2826 : END DO
2827 : END DO
2828 : END DO
2829 : END DO
2830 : END DO
2831 : END DO
2832 :
2833 : ! Free pointers
2834 0 : FRACNOX_LUT => NULL()
2835 0 : DNOx_LUT => NULL()
2836 0 : OPE_LUT => NULL()
2837 0 : MOE_LUT => NULL()
2838 : END DO
2839 :
2840 : ! Transfer MOE if optional output parameter is present
2841 0 : IF ( PRESENT( MOE_OUT ) ) MOE_OUT = MOE
2842 :
2843 : ! Nullify pointers
2844 0 : NULLIFY( FRACNOX_LUT )
2845 0 : NULLIFY( DNOx_LUT )
2846 0 : NULLIFY( OPE_LUT )
2847 0 : NULLIFY( MOE_LUT )
2848 :
2849 : ! Return w/ success
2850 0 : RC = HCO_SUCCESS
2851 :
2852 0 : END SUBROUTINE PARANOX_LUT
2853 : !EOC
2854 : !------------------------------------------------------------------------------
2855 : ! Harmonized Emissions Component (HEMCO) !
2856 : !------------------------------------------------------------------------------
2857 : !BOP
2858 : !
2859 : ! !IROUTINE: InstGet
2860 : !
2861 : ! !DESCRIPTION: Subroutine InstGet returns a poiner to the desired instance.
2862 : !\\
2863 : !\\
2864 : ! !INTERFACE:
2865 : !
2866 0 : SUBROUTINE InstGet ( Instance, Inst, RC, PrevInst )
2867 : !
2868 : ! !INPUT PARAMETERS:
2869 : !
2870 : INTEGER :: Instance
2871 : TYPE(MyInst), POINTER :: Inst
2872 : INTEGER :: RC
2873 : TYPE(MyInst), POINTER, OPTIONAL :: PrevInst
2874 : !
2875 : ! !REVISION HISTORY:
2876 : ! 18 Feb 2016 - C. Keller - Initial version
2877 : ! See https://github.com/geoschem/hemco for complete history
2878 : !EOP
2879 : !------------------------------------------------------------------------------
2880 : !BOC
2881 : TYPE(MyInst), POINTER :: PrvInst
2882 :
2883 : !=================================================================
2884 : ! InstGet begins here!
2885 : !=================================================================
2886 :
2887 : ! Get instance. Also archive previous instance.
2888 0 : PrvInst => NULL()
2889 0 : Inst => AllInst
2890 0 : DO WHILE ( ASSOCIATED(Inst) )
2891 0 : IF ( Inst%Instance == Instance ) EXIT
2892 0 : PrvInst => Inst
2893 0 : Inst => Inst%NextInst
2894 : END DO
2895 0 : IF ( .NOT. ASSOCIATED( Inst ) ) THEN
2896 0 : RC = HCO_FAIL
2897 0 : RETURN
2898 : ENDIF
2899 :
2900 : ! Pass output arguments
2901 0 : IF ( PRESENT(PrevInst) ) PrevInst => PrvInst
2902 :
2903 : ! Cleanup & Return
2904 0 : PrvInst => NULL()
2905 0 : RC = HCO_SUCCESS
2906 :
2907 : END SUBROUTINE InstGet
2908 : !EOC
2909 : !------------------------------------------------------------------------------
2910 : ! Harmonized Emissions Component (HEMCO) !
2911 : !------------------------------------------------------------------------------
2912 : !BOP
2913 : !
2914 : ! !IROUTINE: InstCreate
2915 : !
2916 : ! !DESCRIPTION: Subroutine InstCreate creates a new instance.
2917 : !\\
2918 : !\\
2919 : ! !INTERFACE:
2920 : !
2921 0 : SUBROUTINE InstCreate ( ExtNr, Instance, Inst, RC )
2922 : !
2923 : ! !INPUT PARAMETERS:
2924 : !
2925 : INTEGER, INTENT(IN) :: ExtNr
2926 : !
2927 : ! !OUTPUT PARAMETERS:
2928 : !
2929 : INTEGER, INTENT( OUT) :: Instance
2930 : TYPE(MyInst), POINTER :: Inst
2931 : !
2932 : ! !INPUT/OUTPUT PARAMETERS:
2933 : !
2934 : INTEGER, INTENT(INOUT) :: RC
2935 : !
2936 : ! !REVISION HISTORY:
2937 : ! 18 Feb 2016 - C. Keller - Initial version
2938 : ! See https://github.com/geoschem/hemco for complete history
2939 : !EOP
2940 : !------------------------------------------------------------------------------
2941 : !BOC
2942 : TYPE(MyInst), POINTER :: TmpInst
2943 : INTEGER :: nnInst
2944 :
2945 : !=================================================================
2946 : ! InstCreate begins here!
2947 : !=================================================================
2948 :
2949 : ! ----------------------------------------------------------------
2950 : ! Generic instance initialization
2951 : ! ----------------------------------------------------------------
2952 :
2953 : ! Initialize
2954 0 : Inst => NULL()
2955 :
2956 : ! Get number of already existing instances
2957 0 : TmpInst => AllInst
2958 0 : nnInst = 0
2959 0 : DO WHILE ( ASSOCIATED(TmpInst) )
2960 0 : nnInst = nnInst + 1
2961 0 : TmpInst => TmpInst%NextInst
2962 : END DO
2963 :
2964 : ! Create new instance
2965 0 : ALLOCATE(Inst)
2966 0 : Inst%Instance = nnInst + 1
2967 0 : Inst%ExtNr = ExtNr
2968 :
2969 : ! Attach to instance list
2970 0 : Inst%NextInst => AllInst
2971 0 : AllInst => Inst
2972 :
2973 : ! Update output instance
2974 0 : Instance = Inst%Instance
2975 :
2976 : ! ----------------------------------------------------------------
2977 : ! Type specific initialization statements follow below
2978 : ! ----------------------------------------------------------------
2979 :
2980 : ! Return w/ success
2981 0 : RC = HCO_SUCCESS
2982 :
2983 0 : END SUBROUTINE InstCreate
2984 : !EOC
2985 : !------------------------------------------------------------------------------
2986 : ! Harmonized Emissions Component (HEMCO) !
2987 : !------------------------------------------------------------------------------
2988 : !BOP
2989 : !
2990 : ! !IROUTINE: InstRemove
2991 : !
2992 : ! !DESCRIPTION: Subroutine InstRemove creates a new instance.
2993 : !\\
2994 : !\\
2995 : ! !INTERFACE:
2996 : !
2997 0 : SUBROUTINE InstRemove ( Instance )
2998 : !
2999 : ! !INPUT PARAMETERS:
3000 : !
3001 : INTEGER :: Instance
3002 : !
3003 : ! !REVISION HISTORY:
3004 : ! 18 Feb 2016 - C. Keller - Initial version
3005 : ! See https://github.com/geoschem/hemco for complete history
3006 : !EOP
3007 : !------------------------------------------------------------------------------
3008 : !BOC
3009 : INTEGER :: RC
3010 : TYPE(MyInst), POINTER :: PrevInst
3011 : TYPE(MyInst), POINTER :: Inst
3012 :
3013 : !=================================================================
3014 : ! InstRemove begins here!
3015 : !=================================================================
3016 :
3017 : ! Init
3018 0 : PrevInst => NULL()
3019 0 : Inst => NULL()
3020 :
3021 : ! Get instance. Also archive previous instance.
3022 0 : CALL InstGet ( Instance, Inst, RC, PrevInst=PrevInst )
3023 :
3024 : ! Instance-specific deallocation
3025 0 : IF ( ASSOCIATED(Inst) ) THEN
3026 :
3027 : !---------------------------------------------------------------------
3028 : ! Deallocate fields of Inst before popping off from the list
3029 : ! in order to avoid memory leaks (Bob Yantosca (17 Aug 2022)
3030 : !---------------------------------------------------------------------
3031 0 : IF ( ASSOCIATED( Inst%ShipNO ) ) THEN
3032 0 : DEALLOCATE ( Inst%ShipNO )
3033 : ENDIF
3034 0 : Inst%ShipNO => NULL()
3035 :
3036 0 : IF ( ASSOCIATED( Inst%SC5 ) ) THEN
3037 0 : DEALLOCATE ( Inst%SC5 )
3038 : ENDIF
3039 0 : Inst%SC5 => NULL()
3040 :
3041 0 : IF ( ASSOCIATED( Inst%FRACNOX_LUT02 ) ) THEN
3042 0 : DEALLOCATE( Inst%FRACNOX_LUT02 )
3043 : ENDIF
3044 0 : Inst%FRACNOX_LUT02 => NULL()
3045 :
3046 0 : IF ( ASSOCIATED( Inst%FRACNOX_LUT06 ) ) THEN
3047 0 : DEALLOCATE( Inst%FRACNOX_LUT06 )
3048 : ENDIF
3049 0 : Inst%FRACNOX_LUT06 => NULL()
3050 :
3051 0 : IF ( ASSOCIATED( Inst%FRACNOX_LUT10 ) ) THEN
3052 0 : DEALLOCATE( Inst%FRACNOX_LUT10 )
3053 : ENDIF
3054 0 : Inst%FRACNOX_LUT10 => NULL()
3055 :
3056 0 : IF ( ASSOCIATED( Inst%FRACNOX_LUT14 ) ) THEN
3057 0 : DEALLOCATE( Inst%FRACNOX_LUT14 )
3058 : ENDIF
3059 0 : Inst%FRACNOX_LUT14 => NULL()
3060 :
3061 0 : IF ( ASSOCIATED( Inst%FRACNOX_LUT18 ) ) THEN
3062 0 : DEALLOCATE( Inst%FRACNOX_LUT18 )
3063 : ENDIF
3064 0 : Inst%FRACNOX_LUT18 => NULL()
3065 :
3066 0 : IF ( ASSOCIATED( Inst%OPE_LUT02 ) ) THEN
3067 0 : DEALLOCATE( Inst%OPE_LUT02 )
3068 : ENDIF
3069 0 : Inst%OPE_LUT02 => NULL()
3070 :
3071 0 : IF ( ASSOCIATED( Inst%OPE_LUT06 ) ) THEN
3072 0 : DEALLOCATE( Inst%OPE_LUT06 )
3073 : ENDIF
3074 0 : Inst%OPE_LUT06 => NULL()
3075 :
3076 0 : IF ( ASSOCIATED( Inst%OPE_LUT10 ) ) THEN
3077 0 : DEALLOCATE( Inst%OPE_LUT10 )
3078 : ENDIF
3079 0 : Inst%OPE_LUT10 => NULL()
3080 :
3081 0 : IF ( ASSOCIATED( Inst%OPE_LUT14 ) ) THEN
3082 0 : DEALLOCATE( Inst%OPE_LUT14 )
3083 : ENDIF
3084 0 : Inst%OPE_LUT14 => NULL()
3085 :
3086 0 : IF ( ASSOCIATED( Inst%OPE_LUT18 ) ) THEN
3087 0 : DEALLOCATE( Inst%OPE_LUT18 )
3088 : ENDIF
3089 0 : Inst%OPE_LUT18 => NULL()
3090 :
3091 0 : IF ( ASSOCIATED( Inst%MOE_LUT02 ) ) THEN
3092 0 : DEALLOCATE( Inst%MOE_LUT02 )
3093 : ENDIF
3094 0 : Inst%MOE_LUT02 => NULL()
3095 :
3096 0 : IF ( ASSOCIATED( Inst%MOE_LUT06 ) ) THEN
3097 0 : DEALLOCATE( Inst%MOE_LUT06 )
3098 : ENDIF
3099 0 : Inst%MOE_LUT06 => NULL()
3100 :
3101 0 : IF ( ASSOCIATED( Inst%MOE_LUT10 ) ) THEN
3102 0 : DEALLOCATE( Inst%MOE_LUT10 )
3103 : ENDIF
3104 0 : Inst%MOE_LUT10 => NULL()
3105 :
3106 0 : IF ( ASSOCIATED( Inst%MOE_LUT14 ) ) THEN
3107 0 : DEALLOCATE( Inst%MOE_LUT14 )
3108 : ENDIF
3109 0 : Inst%MOE_LUT14 => NULL()
3110 :
3111 0 : IF ( ASSOCIATED( Inst%MOE_LUT18 ) ) THEN
3112 0 : DEALLOCATE( Inst%MOE_LUT18 )
3113 : ENDIF
3114 0 : Inst%MOE_LUT18 => NULL()
3115 :
3116 0 : IF ( ASSOCIATED( Inst%DNOx_LUT02 ) ) THEN
3117 0 : DEALLOCATE( Inst%DNOx_LUT02 )
3118 : ENDIF
3119 0 : Inst%DNOx_LUT02 => NULL()
3120 :
3121 0 : IF ( ASSOCIATED( Inst%DNOx_LUT06 ) ) THEN
3122 0 : DEALLOCATE( Inst%DNOx_LUT06 )
3123 : ENDIF
3124 0 : Inst%DNOx_LUT06 => NULL()
3125 :
3126 0 : IF ( ASSOCIATED( Inst%DNOx_LUT10 ) ) THEN
3127 0 : DEALLOCATE( Inst%DNOx_LUT10 )
3128 : ENDIF
3129 0 : Inst%DNOx_LUT10 => NULL()
3130 :
3131 0 : IF ( ASSOCIATED( Inst%DNOx_LUT14 ) ) THEN
3132 0 : DEALLOCATE( Inst%DNOx_LUT14 )
3133 : ENDIF
3134 0 : Inst%DNOx_LUT14 => NULL()
3135 :
3136 0 : IF ( ASSOCIATED( Inst%DNOx_LUT18 ) ) THEN
3137 0 : DEALLOCATE( Inst%DNOx_LUT18 )
3138 : ENDIF
3139 0 : Inst%DNOx_LUT18 => NULL()
3140 :
3141 0 : IF ( ASSOCIATED( Inst%DEPO3 ) ) THEN
3142 0 : DEALLOCATE( Inst%DEPO3 )
3143 : ENDIF
3144 0 : Inst%DEPO3 => NULL()
3145 :
3146 0 : IF ( ASSOCIATED( Inst%DEPHNO3 ) ) THEN
3147 0 : DEALLOCATE( Inst%DEPHNO3 )
3148 : ENDIF
3149 0 : Inst%DEPHNO3 => NULL()
3150 :
3151 : !---------------------------------------------------------------------
3152 : ! Pop off instance from list
3153 : !---------------------------------------------------------------------
3154 0 : IF ( ASSOCIATED(PrevInst) ) THEN
3155 0 : PrevInst%NextInst => Inst%NextInst
3156 : ELSE
3157 0 : AllInst => Inst%NextInst
3158 : ENDIF
3159 0 : DEALLOCATE(Inst)
3160 : ENDIF
3161 :
3162 : ! Free pointers before exiting
3163 0 : PrevInst => NULL()
3164 0 : Inst => NULL()
3165 :
3166 0 : END SUBROUTINE InstRemove
3167 : !EOC
3168 0 : END MODULE HCOX_ParaNOx_mod
|
