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 :
1084 : ! Write the name of the extension regardless of the verbose setting
1085 0 : msg = 'Using HEMCO extension: ParaNOx (ship emission plumes)'
1086 0 : IF ( HCO_IsVerb( HcoState%Config%Err ) ) THEN
1087 0 : CALL HCO_Msg( HcoState%Config%Err, msg, sep1='-' ) ! with separator
1088 : ELSE
1089 0 : CALL HCO_Msg( msg, verb=.TRUE. ) ! w/o separator
1090 : ENDIF
1091 :
1092 : ! Write the rest of the information only when verbose is set
1093 0 : MSG = ' - Use the following species: (MW, emitted as HEMCO ID) '
1094 0 : CALL HCO_MSG(HcoState%Config%Err,MSG )
1095 0 : WRITE(MSG,"(a,F5.2,I5)") ' NO : ', Inst%MW_NO, Inst%IDTNO
1096 0 : CALL HCO_MSG(HcoState%Config%Err,MSG)
1097 0 : WRITE(MSG,"(a,F5.2,I5)") ' NO2 : ', Inst%MW_NO2, Inst%IDTNO2
1098 0 : CALL HCO_MSG(HcoState%Config%Err,MSG)
1099 0 : WRITE(MSG,"(a,F5.2,I5)") ' O3 : ', Inst%MW_O3, Inst%IDTO3
1100 0 : CALL HCO_MSG(HcoState%Config%Err,MSG)
1101 0 : WRITE(MSG,"(a,F5.2,I5)") ' HNO3: ', Inst%MW_HNO3, Inst%IDTHNO3
1102 0 : CALL HCO_MSG(HcoState%Config%Err,MSG)
1103 : ENDIF
1104 :
1105 : !--------------------------------
1106 : ! Allocate module arrays
1107 : !--------------------------------
1108 :
1109 : ! FNOX
1110 0 : ALLOCATE( Inst%FRACNOX_LUT02(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1111 0 : IF ( RC /= HCO_SUCCESS ) THEN
1112 0 : CALL HCO_ERROR ( 'FRACNOX_LUT02', RC )
1113 0 : RETURN
1114 : ENDIF
1115 0 : Inst%FRACNOX_LUT02 = 0.0_sp
1116 :
1117 0 : ALLOCATE( Inst%FRACNOX_LUT06(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1118 0 : IF ( RC /= HCO_SUCCESS ) THEN
1119 0 : CALL HCO_ERROR ( 'FRACNOX_LUT06', RC )
1120 0 : RETURN
1121 : ENDIF
1122 0 : Inst%FRACNOX_LUT06 = 0.0_sp
1123 :
1124 0 : ALLOCATE( Inst%FRACNOX_LUT10(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1125 0 : IF ( RC /= HCO_SUCCESS ) THEN
1126 0 : CALL HCO_ERROR ( 'FRACNOX_LUT10', RC )
1127 0 : RETURN
1128 : ENDIF
1129 0 : Inst%FRACNOX_LUT10 = 0.0_sp
1130 :
1131 0 : ALLOCATE( Inst%FRACNOX_LUT14(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1132 0 : IF ( RC /= HCO_SUCCESS ) THEN
1133 0 : CALL HCO_ERROR ( 'FRACNOX_LUT014', RC )
1134 0 : RETURN
1135 : ENDIF
1136 0 : Inst%FRACNOX_LUT14 = 0.0_sp
1137 :
1138 0 : ALLOCATE( Inst%FRACNOX_LUT18(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1139 0 : IF ( RC /= HCO_SUCCESS ) THEN
1140 0 : CALL HCO_ERROR ( 'FRACNOX_LUT18', RC )
1141 0 : RETURN
1142 : ENDIF
1143 0 : Inst%FRACNOX_LUT18 = 0.0_sp
1144 :
1145 : ! OPE
1146 0 : ALLOCATE( Inst%OPE_LUT02(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1147 0 : IF ( RC /= HCO_SUCCESS ) THEN
1148 0 : CALL HCO_ERROR ( 'OPE_LUT02', RC )
1149 0 : RETURN
1150 : ENDIF
1151 0 : Inst%OPE_LUT02 = 0.0_sp
1152 :
1153 0 : ALLOCATE( Inst%OPE_LUT06(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1154 0 : IF ( RC /= HCO_SUCCESS ) THEN
1155 0 : CALL HCO_ERROR ( 'OPE_LUT06', RC )
1156 0 : RETURN
1157 : ENDIF
1158 0 : Inst%OPE_LUT06 = 0.0_sp
1159 :
1160 0 : ALLOCATE( Inst%OPE_LUT10(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1161 0 : IF ( RC /= 0 ) THEN
1162 0 : CALL HCO_ERROR ( 'OPE_LUT10', RC )
1163 0 : RETURN
1164 : ENDIF
1165 0 : Inst%OPE_LUT10 = 0.0_sp
1166 :
1167 0 : ALLOCATE( Inst%OPE_LUT14(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1168 0 : IF ( RC /= 0 ) THEN
1169 0 : CALL HCO_ERROR ( 'OPE_LUT014', RC )
1170 0 : RETURN
1171 : ENDIF
1172 0 : Inst%OPE_LUT14 = 0.0_sp
1173 :
1174 0 : ALLOCATE( Inst%OPE_LUT18(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1175 0 : IF ( RC /= HCO_SUCCESS ) THEN
1176 0 : CALL HCO_ERROR ( 'OPE_LUT18', RC )
1177 0 : RETURN
1178 : ENDIF
1179 0 : Inst%OPE_LUT18 = 0.0_sp
1180 :
1181 : ! MOE
1182 0 : ALLOCATE( Inst%MOE_LUT02(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1183 0 : IF ( RC /= HCO_SUCCESS ) THEN
1184 0 : CALL HCO_ERROR ( 'MOE_LUT02', RC )
1185 0 : RETURN
1186 : ENDIF
1187 0 : Inst%MOE_LUT02 = 0.0_sp
1188 :
1189 0 : ALLOCATE( Inst%MOE_LUT06(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1190 0 : IF ( RC /= HCO_SUCCESS ) THEN
1191 0 : CALL HCO_ERROR ( 'MOE_LUT06', RC )
1192 0 : RETURN
1193 : ENDIF
1194 0 : Inst%MOE_LUT06 = 0.0_sp
1195 :
1196 0 : ALLOCATE( Inst%MOE_LUT10(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1197 0 : IF ( RC /= HCO_SUCCESS ) THEN
1198 0 : CALL HCO_ERROR ( 'MOE_LUT10', RC )
1199 0 : RETURN
1200 : ENDIF
1201 0 : Inst%MOE_LUT10 = 0.0_sp
1202 :
1203 0 : ALLOCATE( Inst%MOE_LUT14(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1204 0 : IF ( RC /= HCO_SUCCESS ) THEN
1205 0 : CALL HCO_ERROR ( 'MOE_LUT014', RC )
1206 0 : RETURN
1207 : ENDIF
1208 0 : Inst%MOE_LUT14 = 0.0_sp
1209 :
1210 0 : ALLOCATE( Inst%MOE_LUT18(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1211 0 : IF ( RC /= HCO_SUCCESS ) THEN
1212 0 : CALL HCO_ERROR ( 'MOE_LUT18', RC )
1213 0 : RETURN
1214 : ENDIF
1215 0 : Inst%MOE_LUT18 = 0.0_sp
1216 :
1217 : ! DNOx
1218 0 : ALLOCATE( Inst%DNOx_LUT02(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1219 0 : IF ( RC /= HCO_SUCCESS ) THEN
1220 0 : CALL HCO_ERROR ( 'DNOx_LUT02', RC )
1221 0 : RETURN
1222 : ENDIF
1223 0 : Inst%DNOx_LUT02 = 0.0_sp
1224 :
1225 0 : ALLOCATE( Inst%DNOx_LUT06(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1226 0 : IF ( RC /= HCO_SUCCESS ) THEN
1227 0 : CALL HCO_ERROR ( 'DNOx_LUT06', RC )
1228 0 : RETURN
1229 : ENDIF
1230 0 : Inst%DNOx_LUT06 = 0.0_sp
1231 :
1232 0 : ALLOCATE( Inst%DNOx_LUT10(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1233 0 : IF ( RC /= HCO_SUCCESS ) THEN
1234 0 : CALL HCO_ERROR ( 'DNOx_LUT10', RC )
1235 0 : RETURN
1236 : ENDIF
1237 0 : Inst%DNOx_LUT10 = 0.0_sp
1238 :
1239 0 : ALLOCATE( Inst%DNOx_LUT14(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1240 0 : IF ( RC /= HCO_SUCCESS ) THEN
1241 0 : CALL HCO_ERROR ( 'DNOx_LUT014', RC )
1242 0 : RETURN
1243 : ENDIF
1244 0 : Inst%DNOx_LUT14 = 0.0_sp
1245 :
1246 0 : ALLOCATE( Inst%DNOx_LUT18(nT,nJ,nO3,nSEA,nSEA,nJ,nNOx), STAT=RC )
1247 0 : IF ( RC /= HCO_SUCCESS ) THEN
1248 0 : CALL HCO_ERROR ( 'DNOx_LUT18', RC )
1249 0 : RETURN
1250 : ENDIF
1251 0 : Inst%DNOx_LUT18 = 0.0_sp
1252 :
1253 : ALLOCATE(Inst%DEPO3 (HcoState%NX,HcoState%NY), &
1254 0 : Inst%DEPHNO3(HcoState%NX,HcoState%NY), STAT=RC )
1255 0 : IF ( RC /= HCO_SUCCESS ) THEN
1256 0 : CALL HCO_ERROR ( 'Deposition arrays', RC )
1257 0 : RETURN
1258 : ENDIF
1259 0 : Inst%DEPO3 = 0.0_sp
1260 0 : Inst%DEPHNO3 = 0.0_sp
1261 :
1262 : ! ! O3 loss and HNO3 deposition
1263 : ! ALLOCATE( Inst%SHIPO3LOSS(HcoState%NX,HcoState%NY), STAT=RC )
1264 : ! IF ( RC /= HCO_SUCCESS ) THEN
1265 : ! CALL HCO_ERROR ( 'SHIPO3LOSS', RC )
1266 : ! RETURN
1267 : ! ENDIF
1268 : ! Inst%SHIPO3LOSS = 0d0
1269 :
1270 : ! ALLOCATE( Inst%SHIPHNO3DEP(HcoState%NX,HcoState%NY), STAT=RC )
1271 : ! IF ( RC /= HCO_SUCCESS ) THEN
1272 : ! CALL HCO_ERROR ( 'SHIPHNO3DEP', RC ); RETURN
1273 : ! ENDIF
1274 : ! Inst%SHIPHNO3DEP = 0d0
1275 :
1276 : ENDIF
1277 :
1278 : !========================================================================
1279 : ! Initialize the PARANOX look-up tables
1280 : !========================================================================
1281 :
1282 : ! LUT data directory
1283 : CALL GetExtOpt( HcoState%Config, Inst%ExtNr, 'LUT source dir', &
1284 0 : OptValChar=Inst%LutDir, RC=RC)
1285 0 : IF ( RC /= HCO_SUCCESS ) THEN
1286 : CALL HCO_ERROR( &
1287 0 : 'PARANOX: Could not read "LUT source dir"!', RC, THISLOC=LOC )
1288 0 : RETURN
1289 : ENDIF
1290 :
1291 : ! Call HEMCO parser to replace tokens such as $ROOT, $MET, or $RES.
1292 : ! There shouldn't be any date token in there ($YYYY, etc.), so just
1293 : ! provide some dummy variables here
1294 0 : CALL HCO_CharParse( HcoState%Config, Inst%LutDir, -999, -1, -1, -1, -1, RC )
1295 0 : IF ( RC /= HCO_SUCCESS ) THEN
1296 : CALL HCO_ERROR( &
1297 0 : 'PARANOX: Error encountered in "HCO_CharParse"', RC, THISLOC=LOC )
1298 0 : RETURN
1299 : ENDIF
1300 :
1301 : ! Data format: ncdf (default) or txt
1302 0 : Inst%IsNc = .TRUE.
1303 : CALL GetExtOpt( HcoState%Config, Inst%ExtNr, 'LUT data format', &
1304 0 : OptValChar=Dummy, RC=RC)
1305 0 : IF ( RC /= HCO_SUCCESS ) THEN
1306 : CALL HCO_ERROR( &
1307 0 : 'PARANOX: Could not read "LUT data format"', RC, THISLOC=LOC )
1308 0 : RETURN
1309 : ENDIF
1310 0 : IF ( TRIM(Dummy) == 'txt' ) Inst%IsNc = .FALSE.
1311 :
1312 : ! Read PARANOX look-up tables from disk. This can be netCDF or txt
1313 : ! format, as determined above.
1314 : !
1315 : ! NOTE: For the GEOS-Chem dry-run or HEMCO-standalone dry-run,
1316 : ! these routines will print file paths to the dry-run log file,
1317 : ! but will not actually read any data.
1318 0 : IF ( Inst%IsNc ) THEN
1319 0 : CALL READ_PARANOX_LUT_NC( HcoState, Inst, RC )
1320 0 : IF ( RC /= HCO_SUCCESS ) THEN
1321 : CALL HCO_ERROR( &
1322 0 : 'PARANOX: Error in "READ_PARANOX_LUT_NC"!', RC, THISLOC=LOC )
1323 0 : RETURN
1324 : ENDIF
1325 : ELSE
1326 0 : CALL READ_PARANOX_LUT_TXT( HcoState, Inst, RC )
1327 0 : IF ( RC /= HCO_SUCCESS ) THEN
1328 : CALL HCO_ERROR( &
1329 0 : 'PARANOX: Error in "READ_PARANOX_LUT_NC"!', RC, THISLOC=LOC )
1330 0 : RETURN
1331 : ENDIF
1332 : ENDIF
1333 :
1334 : !========================================================================
1335 : ! Exit if this is a GEOS-Chem dry-run or HEMCO-standalone dry-run
1336 : !========================================================================
1337 0 : IF ( HcoState%Options%IsDryRun ) THEN
1338 0 : Inst => NULL()
1339 0 : CALL HCO_LEAVE( HcoState%Config%Err,RC )
1340 0 : RETURN
1341 : ENDIF
1342 :
1343 : !========================================================================
1344 : ! Continue initializing PARANOX for regular simulations
1345 : !========================================================================
1346 :
1347 : !------------------------------------------------------------------------
1348 : ! Set other module variables
1349 : !------------------------------------------------------------------------
1350 0 : ALLOCATE ( Inst%ShipNO(HcoState%NX,HcoState%NY,HcoState%NZ), STAT=RC )
1351 0 : IF ( RC /= HCO_SUCCESS ) THEN
1352 0 : CALL HCO_ERROR ( 'ShipNO', RC )
1353 0 : RETURN
1354 : ENDIF
1355 0 : Inst%ShipNO = 0.0_hp
1356 :
1357 : ! Allocate variables for SunCosMid from 5 hours ago.
1358 0 : ALLOCATE ( Inst%SC5(HcoState%NX,HcoState%NY), STAT=RC )
1359 0 : IF ( RC /= HCO_SUCCESS ) THEN
1360 0 : CALL HCO_ERROR ( 'SC5', RC )
1361 0 : RETURN
1362 : ENDIF
1363 0 : Inst%SC5 = 0.0_hp
1364 :
1365 : ! Prompt warning if chemistry time step is more than 60 mins
1366 0 : IF ( HcoState%TS_CHEM > 3600.0_hp ) THEN
1367 0 : IF ( HcoState%amIRoot ) THEN
1368 : MSG = ' Cannot properly store SUNCOS values ' // &
1369 0 : ' because chemistry time step is more than 60 mins!'
1370 0 : CALL HCO_WARNING(HcoState%Config%Err, MSG, RC )
1371 : ENDIF
1372 : ENDIF
1373 :
1374 : ! Molecular weight of AIR
1375 0 : Inst%MW_AIR = HcoState%Phys%AIRMW
1376 :
1377 : !------------------------------------------------------------------------
1378 : ! Define PARANOX diagnostics for the O3 and HNO3 deposition fluxes (in
1379 : ! kg/m2/s).
1380 : !------------------------------------------------------------------------
1381 : CALL Diagn_Create ( HcoState, &
1382 : cName = 'PARANOX_O3_DEPOSITION_FLUX', &
1383 : Trgt2D = Inst%DEPO3, &
1384 : SpaceDim = 2, &
1385 : OutUnit = 'kg/m2/s', &
1386 : COL = HcoState%Diagn%HcoDiagnIDManual, &
1387 0 : RC = RC )
1388 0 : IF ( RC /= HCO_SUCCESS ) THEN
1389 0 : CALL HCO_ERROR( 'ERROR 17', RC, THISLOC=LOC )
1390 0 : RETURN
1391 : ENDIF
1392 :
1393 : CALL Diagn_Create ( HcoState, &
1394 : cName = 'PARANOX_HNO3_DEPOSITION_FLUX', &
1395 : Trgt2D = Inst%DEPHNO3, &
1396 : SpaceDim = 2, &
1397 : OutUnit = 'kg/m2/s', &
1398 : COL = HcoState%Diagn%HcoDiagnIDManual, &
1399 0 : RC = RC )
1400 0 : IF ( RC /= HCO_SUCCESS ) THEN
1401 0 : CALL HCO_ERROR( 'ERROR 18', RC, THISLOC=LOC )
1402 0 : RETURN
1403 : ENDIF
1404 :
1405 : !------------------------------------------------------------------------
1406 : ! Set HEMCO extension variables
1407 : !------------------------------------------------------------------------
1408 :
1409 : ! Met. data required by module
1410 0 : ExtState%O3%DoUse = .TRUE.
1411 0 : ExtState%NO2%DoUse = .TRUE.
1412 0 : ExtState%NO%DoUse = .TRUE.
1413 0 : ExtState%AIR%DoUse = .TRUE.
1414 0 : ExtState%AIRVOL%DoUse = .TRUE.
1415 0 : ExtState%SUNCOS%DoUse = .TRUE.
1416 0 : ExtState%T2M%DoUse = .TRUE.
1417 0 : ExtState%U10M%DoUse = .TRUE.
1418 0 : ExtState%V10M%DoUse = .TRUE.
1419 0 : ExtState%FRAC_OF_PBL%DoUse = .TRUE.
1420 0 : IF ( Inst%IDTHNO3 > 0 ) THEN
1421 0 : ExtState%HNO3%DoUse = .TRUE.
1422 : ENDIF
1423 0 : ExtState%JNO2%DoUse = .TRUE.
1424 0 : ExtState%JOH%DoUse = .TRUE.
1425 :
1426 : !------------------------------------------------------------------------
1427 : ! Leave w/ success
1428 : !------------------------------------------------------------------------
1429 0 : IF ( ALLOCATED(HcoIDs ) ) DEALLOCATE(HcoIDs )
1430 0 : IF ( ALLOCATED(SpcNames) ) DEALLOCATE(SpcNames)
1431 0 : Inst => NULL()
1432 0 : CALL HCO_LEAVE( HcoState%Config%Err,RC )
1433 :
1434 0 : END SUBROUTINE HCOX_ParaNOx_Init
1435 : !EOC
1436 : !------------------------------------------------------------------------------
1437 : ! Harmonized Emissions Component (HEMCO) !
1438 : !------------------------------------------------------------------------------
1439 : !BOP
1440 : !
1441 : ! !IROUTINE: HCOX_ParaNOx_Final
1442 : !
1443 : ! !DESCRIPTION: Subroutine HcoX\_ParaNox\_Final finalizes the HEMCO
1444 : ! PARANOX extension.
1445 : !\\
1446 : !\\
1447 : ! !INTERFACE:
1448 : !
1449 0 : SUBROUTINE HCOX_ParaNOx_Final( HcoState, ExtState, RC )
1450 : !
1451 : ! !USES:
1452 : !
1453 : !
1454 : ! !INPUT PARAMETERS:
1455 : !
1456 : TYPE(HCO_State), POINTER :: HcoState ! HEMCO State obj
1457 : TYPE(Ext_State), POINTER :: ExtState ! Module options
1458 : !
1459 : ! !INPUT/OUTPUT PARAMETERS:
1460 : !
1461 : INTEGER, INTENT(INOUT) :: RC
1462 : !
1463 : ! !REVISION HISTORY:
1464 : ! 06 Aug 2013 - C. Keller - Initial Version
1465 : ! See https://github.com/geoschem/hemco for complete history
1466 : !EOP
1467 : !------------------------------------------------------------------------------
1468 : !BOC
1469 : !
1470 : ! LOCAL VARIABLES:
1471 : !
1472 :
1473 : !=================================================================
1474 : ! HCOX_PARANOX_FINAL begins here!
1475 : !=================================================================
1476 0 : CALL InstRemove( ExtState%ParaNOx )
1477 :
1478 0 : RC = HCO_SUCCESS
1479 :
1480 0 : END SUBROUTINE HCOX_ParaNOx_Final
1481 : !EOC
1482 : !------------------------------------------------------------------------------
1483 : ! GEOS-Chem Global Chemical Transport Model !
1484 : !------------------------------------------------------------------------------
1485 : !BOP
1486 : !
1487 : ! !IROUTINE: read_paranox_lut_nc
1488 : !
1489 : ! !DESCRIPTION: Subroutine READ\_PARANOX\_LUT\_NC reads look-up tables in
1490 : ! netCDF format for use in the PARANOX ship plume model (G.C.M. Vinken)
1491 : !\\
1492 : !\\
1493 : ! !INTERFACE:
1494 : !
1495 0 : SUBROUTINE READ_PARANOX_LUT_NC ( HcoState, Inst, RC )
1496 : !
1497 : ! !USES:
1498 : !
1499 : ! !INPUT ARGUMENTS:
1500 : !
1501 : TYPE(HCO_State), POINTER :: HcoState ! HEMCO State object
1502 : TYPE(MyInst), POINTER :: Inst
1503 : !
1504 : ! !INPUT/OUTPUT ARGUMENTS:
1505 : !
1506 : INTEGER, INTENT(INOUT) :: RC
1507 : !
1508 : ! !REVISION HISTORY:
1509 : ! 06 Feb 2012 - M. Payer - Initial version modified from code provided by
1510 : ! G.C.M. Vinken
1511 : ! See https://github.com/geoschem/hemco for complete history
1512 : !EOP
1513 : !------------------------------------------------------------------------------
1514 : !BOC
1515 : !
1516 : ! !LOCAL VARIABLES:
1517 : INTEGER :: IOS
1518 : CHARACTER(LEN=255) :: FILENAME
1519 : CHARACTER(LEN=255) :: MSG
1520 : INTEGER :: fID
1521 :
1522 : !=================================================================
1523 : ! READ_PARANOX_LUT_NC begins here
1524 : !=================================================================
1525 :
1526 : ! NetCDF reading of PARANOX LUT not supported in ESMF environment
1527 : #if defined(ESMF_)
1528 : MSG = 'In ESMF, cannot read PARANOX look-up-table in netCDF ' // &
1529 : 'format. Please set `LUT data format` to `txt` in the ' // &
1530 : 'HEMCO configuration file.'
1531 : CALL HCO_ERROR(MSG, RC, &
1532 : THISLOC = 'READ_PARANOX_LUT_NC (hcox_paranox_mod.F90)' )
1533 : RETURN
1534 : #else
1535 :
1536 : ! Clear FILENAME
1537 0 : FILENAME = ''
1538 :
1539 : ! FILENAME format string
1540 : 101 FORMAT( a, '/ship_plume_lut_', I2.2, 'ms.nc' )
1541 :
1542 : ! Wind speed levels correspond to the files we will read below
1543 0 : Inst%WSlev = (/ 2.0e0, 6.0e0, 10.0e0, 14.0e0, 18.0e0 /)
1544 :
1545 : ! Read 2m/s LUT
1546 0 : WRITE( FILENAME, 101 ) TRIM(Inst%LutDir), 2
1547 : CALL READ_LUT_NCFILE( HcoState, TRIM( FILENAME ), &
1548 : Inst%FRACNOX_LUT02, Inst%DNOx_LUT02, Inst%OPE_LUT02, Inst%MOE_LUT02, &
1549 : Inst%Tlev, Inst%JNO2lev, Inst%O3lev, Inst%SEA0lev, Inst%SEA5lev, &
1550 0 : Inst%JRATIOlev, Inst%NOXlev, RC=RC)
1551 :
1552 : ! Read 6 m/s LUT
1553 0 : WRITE( FILENAME, 101 ) TRIM(Inst%LutDir), 6
1554 : CALL READ_LUT_NCFILE( HcoState, TRIM( FILENAME ), &
1555 : Inst%FRACNOX_LUT06, Inst%DNOx_LUT06, Inst%OPE_LUT06, &
1556 0 : Inst%MOE_LUT06, RC=RC )
1557 :
1558 : ! Read 10 m/s LUT
1559 0 : WRITE( FILENAME, 101 ) TRIM(Inst%LutDir), 10
1560 : CALL READ_LUT_NCFILE( HcoState, TRIM( FILENAME ), &
1561 : Inst%FRACNOX_LUT10, Inst%DNOx_LUT10, Inst%OPE_LUT10, &
1562 0 : Inst%MOE_LUT10, RC=RC )
1563 :
1564 : ! Read 14 m/s LUT
1565 0 : WRITE( FILENAME, 101 ) TRIM(Inst%LutDir), 14
1566 : CALL READ_LUT_NCFILE( HcoState, TRIM( FILENAME ), &
1567 : Inst%FRACNOX_LUT14, Inst%DNOx_LUT14, Inst%OPE_LUT14, &
1568 0 : Inst%MOE_LUT14, RC=RC )
1569 :
1570 : ! Read 18 m/s LUT
1571 0 : WRITE( FILENAME, 101 ) TRIM(Inst%LutDir), 18
1572 : CALL READ_LUT_NCFILE( HcoState, TRIM( FILENAME ), &
1573 : Inst%FRACNOX_LUT18, Inst%DNOx_LUT18, Inst%OPE_LUT18, &
1574 0 : Inst%MOE_LUT18, RC=RC )
1575 :
1576 : ! ! To write into txt-file, uncomment the following lines
1577 : ! FILENAME = TRIM(LutDir)//'/ship_plume_lut_02ms.txt'
1578 : ! CALL WRITE_LUT_TXTFILE( HcoState, TRIM( FILENAME ), &
1579 : ! FRACNOX_LUT02, DNOx_LUT02, OPE_LUT02, MOE_LUT02, RC, &
1580 : ! Tlev, JNO2lev, O3lev, SEA0lev, SEA5lev, JRATIOlev, NOXlev )
1581 : ! IF ( RC /= HCO_SUCCESS ) RETURN
1582 : !
1583 : ! FILENAME = TRIM(LutDir)//'/ship_plume_lut_06ms.txt'
1584 : ! CALL WRITE_LUT_TXTFILE( HcoState, TRIM( FILENAME ), &
1585 : ! FRACNOX_LUT06, DNOx_LUT06, OPE_LUT06, MOE_LUT06, RC, &
1586 : ! Tlev, JNO2lev, O3lev, SEA0lev, SEA5lev, JRATIOlev, NOXlev )
1587 : ! IF ( RC /= HCO_SUCCESS ) RETURN
1588 : !
1589 : ! FILENAME = TRIM(LutDir)//'/ship_plume_lut_10ms.txt'
1590 : ! CALL WRITE_LUT_TXTFILE( HcoState, TRIM( FILENAME ), &
1591 : ! FRACNOX_LUT10, DNOx_LUT10, OPE_LUT10, MOE_LUT10, RC, &
1592 : ! Tlev, JNO2lev, O3lev, SEA0lev, SEA5lev, JRATIOlev, NOXlev )
1593 : ! IF ( RC /= HCO_SUCCESS ) RETURN
1594 : !
1595 : ! FILENAME = TRIM(LutDir)//'/ship_plume_lut_14ms.txt'
1596 : ! CALL WRITE_LUT_TXTFILE( HcoState, TRIM( FILENAME ), &
1597 : ! FRACNOX_LUT14, DNOx_LUT14, OPE_LUT14, MOE_LUT14, RC, &
1598 : ! Tlev, JNO2lev, O3lev, SEA0lev, SEA5lev, JRATIOlev, NOXlev )
1599 : ! IF ( RC /= HCO_SUCCESS ) RETURN
1600 : !
1601 : ! FILENAME = TRIM(LutDir)//'/ship_plume_lut_14ms.txt'
1602 : ! CALL WRITE_LUT_TXTFILE( HcoState, TRIM( FILENAME ), &
1603 : ! FRACNOX_LUT18, DNOx_LUT18, OPE_LUT18, MOE_LUT18, RC, &
1604 : ! Tlev, JNO2lev, O3lev, SEA0lev, SEA5lev, JRATIOlev, NOXlev )
1605 : ! IF ( RC /= HCO_SUCCESS ) RETURN
1606 :
1607 : ! Return w/ success
1608 0 : RC = HCO_SUCCESS
1609 : #endif
1610 :
1611 0 : END SUBROUTINE READ_PARANOX_LUT_NC
1612 : !EOC
1613 : #if !defined(ESMF_)
1614 : !------------------------------------------------------------------------------
1615 : ! GEOS-Chem Global Chemical Transport Model !
1616 : !------------------------------------------------------------------------------
1617 : !BOP
1618 : !
1619 : ! !IROUTINE: read_lut_ncfile
1620 : !
1621 : ! !DESCRIPTION: Subroutine READ\_LUT\_NCFILE reads look up tables for use in
1622 : ! the PARANOX ship plume model (C. Holmes)
1623 : !\\
1624 : !\\
1625 : ! !INTERFACE:
1626 : !
1627 0 : SUBROUTINE READ_LUT_NCFILE( HcoState, FILENAME, FNOX, &
1628 0 : DNOx, OPE, MOE, T, &
1629 0 : JNO2, O3, SEA0, SEA5, &
1630 0 : JRATIO, NOX, RC )
1631 : !
1632 : ! !USES:
1633 : !
1634 : ! Modules for netCDF read
1635 : USE HCO_m_netcdf_io_open
1636 : USE HCO_m_netcdf_io_get_dimlen
1637 : USE HCO_m_netcdf_io_read
1638 : USE HCO_m_netcdf_io_readattr
1639 : USE HCO_m_netcdf_io_close
1640 :
1641 : # include "netcdf.inc"
1642 : !
1643 : ! !INPUT PARAMETERS:
1644 : !
1645 : TYPE(HCO_State), POINTER :: HcoState ! HEMCO State object
1646 : CHARACTER(LEN=*),INTENT(IN) :: FILENAME
1647 : !
1648 : ! !OUTPUT PARAMETERS:
1649 : !
1650 : REAL*4, INTENT(OUT), DIMENSION(:,:,:,:,:,:,:) :: FNOX,OPE,MOE,DNOx
1651 : REAL*4, INTENT(OUT), OPTIONAL :: T(:), JNO2(:), O3(:)
1652 : REAL*4, INTENT(OUT), OPTIONAL :: SEA0(:), SEA5(:)
1653 : REAL*4, INTENT(OUT), OPTIONAL :: JRATIO(:), NOX(:)
1654 : INTEGER, INTENT(OUT), OPTIONAL :: RC
1655 : !
1656 : ! !REVISION HISTORY:
1657 : ! 06 Feb 2012 - M. Payer - Initial version modified from code provided by
1658 : ! G.C.M. Vinken
1659 : ! See https://github.com/geoschem/hemco for complete history
1660 : !EOP
1661 : !------------------------------------------------------------------------------
1662 : !BOC
1663 : !
1664 : ! !LOCAL VARIABLES:
1665 :
1666 : ! Scalars
1667 : LOGICAL :: FileExists
1668 : INTEGER :: AS, IOS
1669 : INTEGER :: fID, HMRC
1670 :
1671 : ! arrays
1672 : INTEGER :: st1d(1), ct1d(1)
1673 : INTEGER :: st7d(7), ct7d(7)
1674 :
1675 : CHARACTER(LEN=255) :: MSG, FileMsg
1676 :
1677 : !=================================================================
1678 : ! In dry-run mode, print file path to dryrun log and exit.
1679 : ! Otherwise, print file path to the HEMCO log file and continue.
1680 : !=================================================================
1681 :
1682 : ! Test if the file exists
1683 0 : INQUIRE( FILE=TRIM( FileName ), EXIST=FileExists )
1684 :
1685 : ! Create a display string based on whether or not the file is found
1686 0 : IF ( FileExists ) THEN
1687 0 : FileMsg = 'HEMCO (PARANOX): Opening'
1688 : ELSE
1689 0 : FileMsg = 'HEMCO (PARANOX): REQUIRED FILE NOT FOUND'
1690 : ENDIF
1691 :
1692 : ! Print file status to stdout and the HEMCO log
1693 0 : IF ( HcoState%amIRoot ) THEN
1694 0 : WRITE( 6, 300 ) TRIM( FileMsg ), TRIM( FileName )
1695 0 : WRITE( MSG, 300 ) TRIM( FileMsg ), TRIM( FileName )
1696 0 : CALL HCO_MSG( HcoState%Config%Err, MSG )
1697 : 300 FORMAT( a, ' ', a )
1698 : ENDIF
1699 :
1700 : ! For dry-run simulations, return to calling program.
1701 : ! For regular simulations, throw an error if we can't find the file.
1702 0 : IF ( HcoState%Options%IsDryRun ) THEN
1703 0 : RETURN
1704 : ELSE
1705 0 : IF ( .not. FileExists ) THEN
1706 0 : WRITE( MSG, 300 ) TRIM( FileMsg ), TRIM( FileName )
1707 0 : CALL HCO_ERROR(MSG, HMRC )
1708 0 : IF ( PRESENT( RC ) ) RC = HMRC
1709 0 : RETURN
1710 : ENDIF
1711 : ENDIF
1712 :
1713 : !=================================================================
1714 : ! READ_LUT_NCFILE begins here!
1715 : !=================================================================
1716 :
1717 : ! Open file for reading
1718 0 : CALL Ncop_Rd( fId, TRIM(FILENAME) )
1719 :
1720 : !-----------------------------------------------------------------
1721 : ! Read reference values used to construct the Look-up table
1722 : ! These are 1d values
1723 : !-----------------------------------------------------------------
1724 :
1725 : ! Temperature, K
1726 0 : IF ( PRESENT(T) ) THEN
1727 0 : st1d = (/ 1 /)
1728 0 : ct1d = (/ nT /)
1729 0 : CALL NcRd( T, fId, 'T', st1d, ct1d )
1730 : ENDIF
1731 :
1732 : ! J(NO2), 1/s
1733 0 : IF ( PRESENT(JNO2) ) THEN
1734 0 : st1d = (/ 1 /)
1735 0 : ct1d = (/ nJ /)
1736 0 : CALL NcRd( JNO2, fId, 'JNO2', st1d, ct1d )
1737 : ENDIF
1738 :
1739 : ! Ambient O3, ppb
1740 0 : IF ( PRESENT(O3) ) THEN
1741 0 : st1d = (/ 1 /)
1742 0 : ct1d = (/ nO3 /)
1743 0 : CALL NcRd( O3, fId, 'O3', st1d, ct1d )
1744 : ENDIF
1745 :
1746 : ! Solar elevation angle at emission time t=0, deg
1747 0 : IF ( PRESENT(SEA0) ) THEN
1748 0 : st1d = (/ 1 /)
1749 0 : ct1d = (/ nSEA /)
1750 0 : CALL NcRd( SEA0, fId, 'SEA0', st1d, ct1d )
1751 : ENDIF
1752 :
1753 : ! Solar elevation angle at time t=5h, deg
1754 0 : IF ( PRESENT(SEA5) ) THEN
1755 0 : st1d = (/ 1 /)
1756 0 : ct1d = (/ nSEA /)
1757 0 : CALL NcRd( SEA5, fId, 'SEA5', st1d, ct1d )
1758 : ENDIF
1759 :
1760 : ! J(OH) / J(NO2) ratio, s/s
1761 0 : IF ( PRESENT(JRATIO) ) THEN
1762 0 : st1d = (/ 1 /)
1763 0 : ct1d = (/ nJ /)
1764 0 : CALL NcRd( JRatio, fId, 'Jratio', st1d, ct1d )
1765 : ENDIF
1766 :
1767 : ! Ambient NOx, ppt
1768 0 : IF ( PRESENT(NOX) ) THEN
1769 0 : st1d = (/ 1 /)
1770 0 : ct1d = (/ nNOx /)
1771 0 : CALL NcRd( NOX, fId, 'NOx', st1d, ct1d )
1772 : ENDIF
1773 :
1774 : ! Define 7D variables used below
1775 0 : st7d = (/ 1, 1, 1, 1, 1, 1, 1 /)
1776 0 : ct7d = (/ nT,nJ,nO3,nSEA,nSEA,nJ,nNOx /)
1777 :
1778 : !-----------------------------------------------------------------
1779 : ! Read look up table for fraction of NOx remaining for ship
1780 : ! emissions after 5 h [unitless]
1781 : !-----------------------------------------------------------------
1782 :
1783 0 : CALL NcRd( FNOx, fId, 'FNOx', st7d, ct7d )
1784 :
1785 : ! testing only
1786 : ! PRINT*, "binary_fracnox: ", Fnox(1:4,1,1,2,3,4,4)
1787 :
1788 : !-----------------------------------------------------------------
1789 : ! Read look up table for 5-h integrated Ozone Production Efficiency
1790 : ! for ship emissions [molec O3 produced / molec NOx lost]
1791 : !-----------------------------------------------------------------
1792 :
1793 0 : CALL NcRd( OPE, fId, 'OPE', st7d, ct7d )
1794 :
1795 : ! testing only
1796 : ! PRINT*, "binary_intope: ", OPE(1:4,1,1,2,3,4,4)
1797 :
1798 : !-----------------------------------------------------------------
1799 : ! Read look up table for 5-h integrated Methane Oxidation Efficiency
1800 : ! for ship emissions [molec CH4 oxidized / molec NOx emitted]
1801 : !-----------------------------------------------------------------
1802 :
1803 0 : CALL NcRd( MOE, fId, 'MOE', st7d, ct7d )
1804 :
1805 : ! testing only
1806 : ! PRINT*, "binary_intmoe: ", MOE(1:4,1,1,2,3,4,4)
1807 :
1808 : !-----------------------------------------------------------------
1809 : ! Read look up table for 5-h integrated NOx deposition fraction
1810 : ! for ship emissions [molec NOx deposited / molec NOx emitted]
1811 : !-----------------------------------------------------------------
1812 :
1813 0 : CALL NcRd( DNOx, fId, 'DNOx', st7d, ct7d )
1814 :
1815 : ! testing only
1816 : ! PRINT*, "binary_depnox: ", DNOx(1:4,1,1,2,3,4,4)
1817 :
1818 : ! Close netCDF file
1819 0 : CALL NcCl( fId )
1820 :
1821 : END SUBROUTINE READ_LUT_NCFILE
1822 : !EOC
1823 : #endif
1824 : !------------------------------------------------------------------------------
1825 : ! GEOS-Chem Global Chemical Transport Model !
1826 : !------------------------------------------------------------------------------
1827 : !BOP
1828 : !
1829 : ! !IROUTINE: read_paranox_lut_txt
1830 : !
1831 : ! !DESCRIPTION: Subroutine READ\_PARANOX\_LUT\_TXT reads look-up tables in
1832 : ! txt format for use in the PARANOX ship plume model (G.C.M. Vinken)
1833 : !\\
1834 : !\\
1835 : ! !INTERFACE:
1836 : !
1837 0 : SUBROUTINE READ_PARANOX_LUT_TXT ( HcoState, Inst, RC )
1838 : !
1839 : ! !USES:
1840 : !
1841 : ! !INPUT ARGUMENTS:
1842 : !
1843 : TYPE(HCO_State), POINTER :: HcoState ! HEMCO State object
1844 : TYPE(MyInst), POINTER :: Inst
1845 : !
1846 : ! !INPUT/OUTPUT ARGUMENTS:
1847 : !
1848 : INTEGER, INTENT(INOUT) :: RC
1849 : !
1850 : ! !REVISION HISTORY:
1851 : ! 05 Feb 2015 - C. Keller - Initial version modified from code provided by
1852 : ! G.C.M. Vinken and C. Holmes
1853 : ! See https://github.com/geoschem/hemco for complete history
1854 : !EOP
1855 : !------------------------------------------------------------------------------
1856 : !BOC
1857 : !
1858 : ! !LOCAL VARIABLES:
1859 : INTEGER :: IOS
1860 : CHARACTER(LEN=255) :: FILENAME
1861 : CHARACTER(LEN=255) :: MSG, LOC
1862 : INTEGER :: fID
1863 :
1864 : !=================================================================
1865 : ! READ_PARANOX_LUT_TXT begins here
1866 : !=================================================================
1867 0 : LOC = 'READ_PARANOX_LUT_TXT (HCOX_PARANOX_MOD.F90)'
1868 :
1869 : ! Clear FILENAME
1870 0 : FILENAME = ''
1871 :
1872 : ! FILENAME format string
1873 : 101 FORMAT( a, '/ship_plume_lut_', I2.2, 'ms.txt' )
1874 :
1875 : ! Wind speed levels correspond to the files that we will read below
1876 0 : Inst%WSlev = (/ 2.0e0, 6.0e0, 10.0e0, 14.0e0, 18.0e0 /)
1877 :
1878 : ! Read 2m/s LUT
1879 0 : WRITE( FILENAME, 101 ) TRIM(Inst%LutDir), 2
1880 : CALL READ_LUT_TXTFILE( HcoState, TRIM( FILENAME ), &
1881 : Inst%FRACNOX_LUT02, Inst%DNOx_LUT02, Inst%OPE_LUT02, Inst%MOE_LUT02, RC, &
1882 0 : Inst%Tlev, Inst%JNO2lev, Inst%O3lev, Inst%SEA0lev, Inst%SEA5lev, Inst%JRATIOlev, Inst%NOXlev )
1883 0 : IF ( RC /= HCO_SUCCESS ) THEN
1884 0 : CALL HCO_ERROR( 'ERROR 19', RC, THISLOC=LOC )
1885 0 : RETURN
1886 : ENDIF
1887 :
1888 : ! Read 6 m/s LUT
1889 0 : WRITE( FILENAME, 101 ) TRIM(Inst%LutDir), 6
1890 : CALL READ_LUT_TXTFILE( HcoState, TRIM( FILENAME ), &
1891 0 : Inst%FRACNOX_LUT06, Inst%DNOx_LUT06, Inst%OPE_LUT06, Inst%MOE_LUT06, RC )
1892 0 : IF ( RC /= HCO_SUCCESS ) THEN
1893 0 : CALL HCO_ERROR( 'ERROR 20', RC, THISLOC=LOC )
1894 0 : RETURN
1895 : ENDIF
1896 :
1897 : ! Read 10 m/s LUT
1898 0 : WRITE( FILENAME, 101 ) TRIM(Inst%LutDir), 10
1899 : CALL READ_LUT_TXTFILE( HcoState, TRIM( FILENAME ), &
1900 0 : Inst%FRACNOX_LUT10, Inst%DNOx_LUT10, Inst%OPE_LUT10, Inst%MOE_LUT10, RC )
1901 0 : IF ( RC /= HCO_SUCCESS ) THEN
1902 0 : CALL HCO_ERROR( 'ERROR 21', RC, THISLOC=LOC )
1903 0 : RETURN
1904 : ENDIF
1905 :
1906 : ! Read 14 m/s LUT
1907 0 : WRITE( FILENAME, 101 ) TRIM(Inst%LutDir), 14
1908 : CALL READ_LUT_TXTFILE( HcoState, TRIM( FILENAME ), &
1909 0 : Inst%FRACNOX_LUT14, Inst%DNOx_LUT14, Inst%OPE_LUT14, Inst%MOE_LUT14, RC )
1910 0 : IF ( RC /= HCO_SUCCESS ) THEN
1911 0 : CALL HCO_ERROR( 'ERROR 22', RC, THISLOC=LOC )
1912 0 : RETURN
1913 : ENDIF
1914 :
1915 : ! Read 18 m/s LUT
1916 0 : WRITE( FILENAME, 101 ) TRIM(Inst%LutDir), 18
1917 : CALL READ_LUT_TXTFILE( HcoState, TRIM( FILENAME ), &
1918 0 : Inst%FRACNOX_LUT18, Inst%DNOx_LUT18, Inst%OPE_LUT18, Inst%MOE_LUT18, RC )
1919 0 : IF ( RC /= HCO_SUCCESS ) THEN
1920 0 : CALL HCO_ERROR( 'ERROR 23', RC, THISLOC=LOC )
1921 0 : RETURN
1922 : ENDIF
1923 :
1924 : ! Return w/ success
1925 0 : RC = HCO_SUCCESS
1926 :
1927 : END SUBROUTINE READ_PARANOX_LUT_TXT
1928 : !EOC
1929 : !------------------------------------------------------------------------------
1930 : ! GEOS-Chem Global Chemical Transport Model !
1931 : !------------------------------------------------------------------------------
1932 : !BOP
1933 : !
1934 : ! !IROUTINE: read_lut_txtfile
1935 : !
1936 : ! !DESCRIPTION: Subroutine READ\_LUT\_TXTFILE reads look up tables for use in
1937 : ! the PARANOX ship plume model (C. Holmes)
1938 : !\\
1939 : !\\
1940 : ! !INTERFACE:
1941 : !
1942 0 : SUBROUTINE READ_LUT_TXTFILE( HcoState, FILENAME, FNOX, DNOx, OPE, MOE, RC, &
1943 0 : T, JNO2, O3, SEA0, SEA5, JRATIO, NOX )
1944 : !
1945 : ! !USES:
1946 : !
1947 : USE HCO_inquireMod, ONLY : findFreeLUN
1948 : !
1949 : ! !INPUT PARAMETERS:
1950 : !
1951 : TYPE(HCO_State), POINTER :: HcoState ! HEMCO State object
1952 : CHARACTER(LEN=*),INTENT(IN) :: FILENAME
1953 : !
1954 : ! !INPUT/OUTPUT PARAMETERS:
1955 : !
1956 : INTEGER, INTENT(INOUT) :: RC
1957 : !
1958 : ! !OUTPUT PARAMETERS:
1959 : !
1960 : REAL*4, INTENT(OUT), TARGET, DIMENSION(:,:,:,:,:,:,:) :: FNOX,OPE,MOE,DNOx
1961 : REAL*4, INTENT(OUT), OPTIONAL :: T(:), JNO2(:), O3(:)
1962 : REAL*4, INTENT(OUT), OPTIONAL :: SEA0(:), SEA5(:)
1963 : REAL*4, INTENT(OUT), OPTIONAL :: JRATIO(:), NOX(:)
1964 : !
1965 : ! !REVISION HISTORY:
1966 : ! 05 Feb 2015 - C. Keller - Initial version modified from code provided by
1967 : ! G.C.M. Vinken and C. Holmes
1968 : ! See https://github.com/geoschem/hemco for complete history
1969 : !EOP
1970 : !------------------------------------------------------------------------------
1971 : !BOC
1972 : !
1973 : ! !LOCAL VARIABLES:
1974 : !
1975 : ! Scalars
1976 : LOGICAL :: FileExists
1977 : INTEGER :: fId, IOS
1978 :
1979 : ! INTEGER :: I, I1, I2, I3, I4, I5, I6, I7
1980 : ! REAL*4, POINTER :: TMPARR(:,:,:,:,:,:,:) => NULL()
1981 :
1982 : CHARACTER(LEN=255) :: MSG, FileMsg
1983 : !
1984 : ! !DEFINED PARAMETERS:
1985 : !
1986 : CHARACTER(LEN=255), PARAMETER :: FMAT = "(E40.32)"
1987 : CHARACTER(LEN=255), PARAMETER :: LOC = &
1988 : 'READ_LUT_TXTFILE (hcox_paranox_mod.F90)'
1989 :
1990 : !=================================================================
1991 : ! In dry-run mode, print file path to dryrun log and exit.
1992 : ! Otherwise, print file path to the HEMCO log file and continue.
1993 : !=================================================================
1994 :
1995 : ! Test if the file exists
1996 0 : INQUIRE ( FILE=TRIM( FileName ), EXIST=FileExists )
1997 :
1998 : ! Create a display string based on whether or not the file is found
1999 0 : IF ( FileExists ) THEN
2000 0 : FileMsg = 'HEMCO (PARANOX): Opening'
2001 : ELSE
2002 0 : FileMsg = 'HEMCO (PARANOX): REQUIRED FILE NOT FOUND'
2003 : ENDIF
2004 :
2005 : ! Print file status to stdout and the HEMCO log file
2006 0 : IF ( HcoState%amIRoot ) THEN
2007 0 : WRITE( 6, 300 ) TRIM( FileMsg ), TRIM( FileName )
2008 0 : WRITE( MSG, 300 ) TRIM( FileMsg ), TRIM( FileName )
2009 0 : CALL HCO_MSG(HcoState%Config%Err,MSG)
2010 : 300 FORMAT( a, ' ', a )
2011 : ENDIF
2012 :
2013 : ! For dry-run simulations, return to calling program.
2014 : ! For regular simulations, throw an error if we can't find the file.
2015 0 : IF ( HcoState%Options%IsDryRun ) THEN
2016 0 : RETURN
2017 : ELSE
2018 0 : IF ( .not. FileExists ) THEN
2019 0 : WRITE( MSG, 300 ) TRIM( FileMsg ), TRIM( FileName )
2020 0 : CALL HCO_ERROR(MSG, RC )
2021 0 : RETURN
2022 : ENDIF
2023 : ENDIF
2024 :
2025 : !=================================================================
2026 : ! READ_LUT_TXTFILE begins here
2027 : !=================================================================
2028 :
2029 : ! Find a free file LUN
2030 0 : fId = findFreeLUN()
2031 :
2032 : ! Open file for reading
2033 0 : OPEN ( fID, FILE=TRIM(FILENAME), FORM="FORMATTED", IOSTAT=IOS )
2034 0 : IF ( IOS /= 0 ) THEN
2035 0 : CALL HCO_ERROR( 'read_lut_txtfile:1', RC, THISLOC=LOC )
2036 0 : RETURN
2037 : ENDIF
2038 :
2039 : ! Read FNOx
2040 0 : READ( fId, FMT=FMAT, IOSTAT=IOS ) FNOx
2041 0 : IF ( IOS /= 0 ) THEN
2042 0 : CALL HCO_ERROR( 'read_lut_txtfile: FNOx', RC, THISLOC=LOC )
2043 0 : RETURN
2044 : ENDIF
2045 :
2046 : ! Read OPE
2047 0 : READ( fId, FMT=FMAT, IOSTAT=IOS ) OPE
2048 0 : IF ( IOS /= 0 ) THEN
2049 0 : CALL HCO_ERROR( 'read_lut_txtfile: OPE', RC, THISLOC=LOC )
2050 0 : RETURN
2051 : ENDIF
2052 :
2053 : ! Read MOE
2054 0 : READ( fId, FMT=FMAT, IOSTAT=IOS ) MOE
2055 0 : IF ( IOS /= 0 ) THEN
2056 0 : CALL HCO_ERROR( 'read_lut_txtfile: MOE', RC, THISLOC=LOC )
2057 0 : RETURN
2058 : ENDIF
2059 :
2060 : ! Read DNOx
2061 0 : READ( fId, FMT=FMAT, IOSTAT=IOS ) DNOx
2062 0 : IF ( IOS /= 0 ) THEN
2063 0 : CALL HCO_ERROR( 'read_lut_txtfile: DNOx', RC, THISLOC=LOC )
2064 0 : RETURN
2065 : ENDIF
2066 :
2067 : ! Read optional values
2068 0 : IF ( PRESENT(T) ) THEN
2069 0 : READ( fId, FMT=FMAT, IOSTAT=IOS ) T
2070 0 : IF ( IOS /= 0 ) THEN
2071 0 : CALL HCO_ERROR( 'read_lut_txtfile: T', RC, THISLOC=LOC )
2072 0 : RETURN
2073 : ENDIF
2074 : ENDIF
2075 :
2076 0 : IF ( PRESENT(JNO2) ) THEN
2077 0 : READ( fId, FMT=FMAT, IOSTAT=IOS ) JNO2
2078 0 : IF ( IOS /= 0 ) THEN
2079 0 : CALL HCO_ERROR( 'read_lut_txtfile: JNO2', RC, THISLOC=LOC )
2080 0 : RETURN
2081 : ENDIF
2082 : ENDIF
2083 :
2084 0 : IF ( PRESENT(O3) ) THEN
2085 0 : READ( fId, FMT=FMAT, IOSTAT=IOS ) O3
2086 0 : IF ( IOS /= 0 ) THEN
2087 0 : CALL HCO_ERROR( 'read_lut_txtfile: O3', RC, THISLOC=LOC )
2088 0 : RETURN
2089 : ENDIF
2090 : ENDIF
2091 :
2092 0 : IF ( PRESENT(SEA0) ) THEN
2093 0 : READ( fId, FMT=FMAT, IOSTAT=IOS ) SEA0
2094 0 : IF ( IOS /= 0 ) THEN
2095 0 : CALL HCO_ERROR( 'read_lut_txtfile: SEA0', RC, THISLOC=LOC )
2096 0 : RETURN
2097 : ENDIF
2098 : ENDIF
2099 :
2100 0 : IF ( PRESENT(SEA5) ) THEN
2101 0 : READ( fId, FMT=FMAT, IOSTAT=IOS ) SEA5
2102 0 : IF ( IOS /= 0 ) THEN
2103 0 : CALL HCO_ERROR( 'read_lut_txtfile: SEA5', RC, THISLOC=LOC )
2104 0 : RETURN
2105 : ENDIF
2106 : ENDIF
2107 :
2108 0 : IF ( PRESENT(JRATIO) ) THEN
2109 0 : READ( fId, FMT=FMAT, IOSTAT=IOS ) JRATIO
2110 0 : IF ( IOS /= 0 ) THEN
2111 0 : CALL HCO_ERROR( 'read_lut_txtfile: JRATIO', RC, THISLOC=LOC )
2112 0 : RETURN
2113 : ENDIF
2114 : ENDIF
2115 :
2116 0 : IF ( PRESENT(NOX) ) THEN
2117 0 : READ( fId, FMT=FMAT, IOSTAT=IOS ) NOX
2118 0 : IF ( IOS /= 0 ) THEN
2119 0 : CALL HCO_ERROR( 'read_lut_txtfile: NOX', RC, THISLOC=LOC )
2120 0 : RETURN
2121 : ENDIF
2122 : ENDIF
2123 :
2124 : ! ! Read
2125 : ! DO I = 1,4
2126 : !
2127 : ! ! Set pointer to output array
2128 : ! SELECT CASE ( I )
2129 : ! CASE ( 1 )
2130 : ! TMPARR => FNOX
2131 : ! CASE ( 2 )
2132 : ! TMPARR => OPE
2133 : ! CASE ( 3 )
2134 : ! TMPARR => MOE
2135 : ! CASE ( 4 )
2136 : ! TMPARR => DNOx
2137 : ! CASE DEFAULT
2138 : ! CALL HCO_ERROR( 'I > 4', RC, THISLOC=LOC )
2139 : ! RETURN
2140 : ! END SELECT
2141 : !
2142 : ! DO I1 = 1,nT
2143 : ! DO I2 = 1,nJ
2144 : ! DO I3 = 1,nO3
2145 : ! DO I4 = 1,nSEA
2146 : ! DO I5 = 1,nSEA
2147 : ! DO I6 = 1,nJ
2148 : ! DO I7 = 1,nNOx
2149 : ! READ( fId, FMT=FMAT, IOSTAT=IOS ) TMPARR(I1,I2,I3,I4,I5,I6,I7)
2150 : ! ENDDO
2151 : ! ENDDO
2152 : ! ENDDO
2153 : ! ENDDO
2154 : ! ENDDO
2155 : ! ENDDO
2156 : ! ENDDO
2157 : ! ENDDO
2158 :
2159 : ! Close file
2160 0 : CLOSE( fId )
2161 :
2162 : ! Return w/ success
2163 0 : RC = HCO_SUCCESS
2164 :
2165 : END SUBROUTINE READ_LUT_TXTFILE
2166 : !EOC
2167 : !------------------------------------------------------------------------------
2168 : ! GEOS-Chem Global Chemical Transport Model !
2169 : !------------------------------------------------------------------------------
2170 : !BOP
2171 : !
2172 : ! !IROUTINE: write_lut_txtfile
2173 : !
2174 : ! !DESCRIPTION: write\_lut\_txtfile
2175 : !\\
2176 : !\\
2177 : ! !INTERFACE:
2178 : !
2179 : SUBROUTINE WRITE_LUT_TXTFILE( HcoState, FILENAME, FNOX, DNOx, OPE, MOE, RC, &
2180 : T, JNO2, O3, SEA0, SEA5, JRATIO, NOX )
2181 : !
2182 : ! !USES:
2183 : !
2184 : USE HCO_inquireMod, ONLY : findFreeLUN
2185 : !
2186 : ! !INPUT PARAMETERS:
2187 : !
2188 : TYPE(HCO_State), INTENT(INOUT) :: HcoState ! HEMCO state obj
2189 : CHARACTER(LEN=*),INTENT(IN) :: FILENAME
2190 : !
2191 : ! !INPUT/OUTPUT PARAMETERS:
2192 : !
2193 : INTEGER, INTENT(INOUT) :: RC
2194 : !
2195 : ! !OUTPUT PARAMETERS:
2196 : !
2197 : REAL*4, INTENT(IN), TARGET, DIMENSION(:,:,:,:,:,:,:) :: FNOX,OPE,MOE,DNOx
2198 : REAL*4, INTENT(IN), OPTIONAL :: T(:), JNO2(:), O3(:)
2199 : REAL*4, INTENT(IN), OPTIONAL :: SEA0(:), SEA5(:)
2200 : REAL*4, INTENT(IN), OPTIONAL :: JRATIO(:), NOX(:)
2201 : !
2202 : ! !REVISION HISTORY:
2203 : ! 05 Feb 2015 - C. Keller - Initial version modified from code provided by
2204 : ! G.C.M. Vinken and C. Holmes
2205 : ! See https://github.com/geoschem/hemco for complete history
2206 : !EOP
2207 : !------------------------------------------------------------------------------
2208 : !BOC
2209 : !
2210 : ! !LOCAL VARIABLES:
2211 : INTEGER :: fId, IOS
2212 :
2213 : ! INTEGER :: I, I1, I2, I3, I4, I5, I6, I7
2214 : ! REAL*4, POINTER :: TMPARR(:,:,:,:,:,:,:) => NULL()
2215 :
2216 : CHARACTER(LEN=255) :: MSG
2217 : CHARACTER(LEN=255), PARAMETER :: FMAT = "(E40.32)"
2218 : CHARACTER(LEN=255), PARAMETER :: LOC = &
2219 : 'WRITE_LUT_TXTFILE (hcox_paranox_mod.F90)'
2220 :
2221 : !=================================================================
2222 : ! WRITE_LUT_TXTFILE begins here
2223 : !=================================================================
2224 :
2225 : ! Echo info
2226 : IF ( Hcostate%amIRoot ) THEN
2227 : WRITE( MSG, 100 ) TRIM( FILENAME )
2228 : CALL HCO_MSG(HcoState%Config%Err,MSG)
2229 : 100 FORMAT( 'WRITE_LUT_TXTFILE: Writing ', a )
2230 : ENDIF
2231 :
2232 : ! Find a free file LUN
2233 : fId = findFreeLUN()
2234 :
2235 : ! Open file for reading
2236 : OPEN ( fID, FILE=TRIM(FILENAME), ACTION="WRITE", FORM="FORMATTED", IOSTAT=IOS )
2237 : IF ( IOS /= 0 ) THEN
2238 : CALL HCO_ERROR( 'write_lut_txtfile:1', RC, THISLOC=LOC )
2239 : RETURN
2240 : ENDIF
2241 :
2242 : ! Read FNOx
2243 : WRITE( fId, FMT=FMAT, IOSTAT=IOS ) FNOx
2244 : IF ( IOS /= 0 ) THEN
2245 : CALL HCO_ERROR( 'write_lut_txtfile: FNOx', RC, THISLOC=LOC )
2246 : RETURN
2247 : ENDIF
2248 :
2249 : ! Read OPE
2250 : WRITE( fId, FMT=FMAT, IOSTAT=IOS ) OPE
2251 : IF ( IOS /= 0 ) THEN
2252 : CALL HCO_ERROR( 'read_lut_txtfile: OPE', RC, THISLOC=LOC )
2253 : RETURN
2254 : ENDIF
2255 :
2256 : ! Read MOE
2257 : WRITE( fId, FMT=FMAT, IOSTAT=IOS ) MOE
2258 : IF ( IOS /= 0 ) THEN
2259 : CALL HCO_ERROR( 'read_lut_txtfile: MOE', RC, THISLOC=LOC )
2260 : RETURN
2261 : ENDIF
2262 :
2263 : ! Read DNOx
2264 : WRITE( fId, FMT=FMAT, IOSTAT=IOS ) DNOx
2265 : IF ( IOS /= 0 ) THEN
2266 : CALL HCO_ERROR( 'read_lut_txtfile: DNOx', RC, THISLOC=LOC )
2267 : RETURN
2268 : ENDIF
2269 :
2270 : ! Read optional values
2271 : IF ( PRESENT(T) ) THEN
2272 : WRITE( fId, FMT=FMAT, IOSTAT=IOS ) T
2273 : IF ( IOS /= 0 ) THEN
2274 : CALL HCO_ERROR( 'read_lut_txtfile: T', RC, THISLOC=LOC )
2275 : RETURN
2276 : ENDIF
2277 : ENDIF
2278 :
2279 : IF ( PRESENT(JNO2) ) THEN
2280 : WRITE( fId, FMT=FMAT, IOSTAT=IOS ) JNO2
2281 : IF ( IOS /= 0 ) THEN
2282 : CALL HCO_ERROR( 'read_lut_txtfile: JNO2', RC, THISLOC=LOC )
2283 : RETURN
2284 : ENDIF
2285 : ENDIF
2286 :
2287 : IF ( PRESENT(O3) ) THEN
2288 : WRITE( fId, FMT=FMAT, IOSTAT=IOS ) O3
2289 : IF ( IOS /= 0 ) THEN
2290 : CALL HCO_ERROR( 'read_lut_txtfile: O3', RC, THISLOC=LOC )
2291 : RETURN
2292 : ENDIF
2293 : ENDIF
2294 :
2295 : IF ( PRESENT(SEA0) ) THEN
2296 : WRITE( fId, FMT=FMAT, IOSTAT=IOS ) SEA0
2297 : IF ( IOS /= 0 ) THEN
2298 : CALL HCO_ERROR( 'read_lut_txtfile: SEA0', RC, THISLOC=LOC )
2299 : RETURN
2300 : ENDIF
2301 : ENDIF
2302 :
2303 : IF ( PRESENT(SEA5) ) THEN
2304 : WRITE( fId, FMT=FMAT, IOSTAT=IOS ) SEA5
2305 : IF ( IOS /= 0 ) THEN
2306 : CALL HCO_ERROR( 'read_lut_txtfile: SEA5', RC, THISLOC=LOC )
2307 : RETURN
2308 : ENDIF
2309 : ENDIF
2310 :
2311 : IF ( PRESENT(JRATIO) ) THEN
2312 : WRITE( fId, FMT=FMAT, IOSTAT=IOS ) JRATIO
2313 : IF ( IOS /= 0 ) THEN
2314 : CALL HCO_ERROR( 'read_lut_txtfile: JRATIO', RC, THISLOC=LOC )
2315 : RETURN
2316 : ENDIF
2317 : ENDIF
2318 :
2319 : IF ( PRESENT(NOX) ) THEN
2320 : WRITE( fId, FMT=FMAT, IOSTAT=IOS ) NOX
2321 : IF ( IOS /= 0 ) THEN
2322 : CALL HCO_ERROR( 'read_lut_txtfile: NOX', RC, THISLOC=LOC )
2323 : RETURN
2324 : ENDIF
2325 : ENDIF
2326 :
2327 : ! Close file
2328 : CLOSE( fId )
2329 :
2330 : ! Return w/ success
2331 : RC = HCO_SUCCESS
2332 :
2333 : END SUBROUTINE WRITE_LUT_TXTFILE
2334 : !EOC
2335 : !------------------------------------------------------------------------------
2336 : ! GEOS-Chem Global Chemical Transport Model !
2337 : !------------------------------------------------------------------------------
2338 : !BOP
2339 : !
2340 : ! !IROUTINE: INTERPOL_LINWEIGHTS
2341 : !
2342 : ! !DESCRIPTION: Subroutine INTERPOL\_LINWEIGHTS finds the array elements and
2343 : ! weights for piecewise 1-D linear interpolation. The input array of NODES
2344 : ! must be in monotonic ascending order. (C. Holmes 3/27/2014)
2345 : !
2346 : ! If Y is an array containing values of a function evaluated at the points
2347 : ! given in NODES, then its interpolated value at the point VALUESIN will be
2348 : ! Y(VALUEIN) = Y(INDICES(1)) * WEIGHTS(1) +
2349 : ! Y(INDICES(2)) * WEIGHTS(2)
2350 : !
2351 : ! This subroutine finds indices of consecutive nodes that bracket VALUEIN and
2352 : ! weights such that
2353 : ! VALUEIN = NODES(INDICES(1)) * WEIGHTS(1) +
2354 : ! NODES(INDICES(1)+1) * (1-WEIGHTS(1))
2355 : !
2356 : ! For convenience, the returned values of INDICES and WEIGHTS are 2-element
2357 : ! arrays, where
2358 : ! INDICES(2) = INDICES(1)+1 and
2359 : ! WEIGHTS(2) = 1 - WEIGHTS(1)
2360 : !\\
2361 : !\\
2362 : ! !INTERFACE:
2363 : !
2364 0 : SUBROUTINE INTERPOL_LINWEIGHTS( NODES, VALUEIN, INDICES, WEIGHTS )
2365 : !
2366 : ! !USES:
2367 : !
2368 : !
2369 : ! !INPUT PARAMETERS:
2370 : !
2371 : REAL*4,INTENT(IN) :: NODES(:), VALUEIN
2372 : !
2373 : ! !OUTPUT PARAMETERS:
2374 : !
2375 : ! These arrays are always 2 elements each, but declaring
2376 : ! as deferred shape avoids array temporaries
2377 : INTEGER,INTENT(OUT) :: INDICES(:)
2378 : REAL*4, INTENT(OUT) :: WEIGHTS(:)
2379 : !
2380 : ! !REVISION HISTORY:
2381 : ! 03 Jun 2013 - C. Holmes - Initial version
2382 : ! See https://github.com/geoschem/hemco for complete history
2383 : !EOP
2384 : !------------------------------------------------------------------------------
2385 : !BOC
2386 : !
2387 : ! !LOCAL VARIABLES:
2388 : INTEGER :: I
2389 : REAL*8 :: VALUE
2390 :
2391 : !=================================================================
2392 : ! INTERPOL_LINWEIGHTS begins here!
2393 : !=================================================================
2394 :
2395 : ! If larger than largest in LUT, assign largest level values
2396 0 : VALUE = MIN( VALUEIN, MAXVAL( NODES ) )
2397 :
2398 : ! If smaller, assign smallest level value
2399 : !GanLuo+VALUE = MAX( VALUE, MINVAL( NODES ) )
2400 0 : VALUE = MAX( VALUE, MINVAL( NODES )*1.d0 )
2401 :
2402 : ! Initialize
2403 0 : INDICES = (/ 1, 1 /)
2404 :
2405 : ! Loop over interpolation nodes until we find the largest node value
2406 : ! that is less than the desired value
2407 0 : DO I=1, SIZE(NODES)
2408 0 : INDICES(1) = I
2409 0 : IF ( VALUE <= NODES(I+1) ) EXIT
2410 : END DO
2411 :
2412 : ! The next node
2413 0 : INDICES(2) = INDICES(1) + 1
2414 :
2415 : ! Weights for the corresponding node indices
2416 0 : WEIGHTS(1) = ( NODES(I+1) - VALUE ) / ( NODES(I+1) - NODES(I) )
2417 0 : WEIGHTS(2) = 1.0 - WEIGHTS(1)
2418 :
2419 0 : END SUBROUTINE INTERPOL_LINWEIGHTS
2420 : !EOC
2421 : !------------------------------------------------------------------------------
2422 : ! GEOS-Chem Global Chemical Transport Model !
2423 : !------------------------------------------------------------------------------
2424 : !BOP
2425 : !
2426 : ! !IROUTINE: paranox_lut
2427 : !
2428 : ! !DESCRIPTION: Subroutine PARANOX\_LUT returns fractional remainder of
2429 : ! ship NOx (FNOx), fraction of NOx that dry deposits as NOy species (DNOx),
2430 : ! ozone production efficiency (OPE), and methane oxidation
2431 : ! efficiency (MOE) after 5-hrs of plume aging. Values are taken taken from a
2432 : ! lookup table using piecewise linear interpolation. The look-up table is derived
2433 : ! from the PARANOx gaussian plume model (Vinken et al. 2011; Holmes et al. 2014)
2434 : ! (G.C.M. Vinken, KNMI, June 2010; C. Holmes June 2013)
2435 : !
2436 : ! The lookup table uses 8 input variables:
2437 : ! TEMP : model temperature, K
2438 : ! JNO2 : J(NO2) value, 1/s
2439 : ! O3 : concentration O3 in ambient air, ppb
2440 : ! SEA0 : solar elevation angle at emission time 5 hours ago, degree
2441 : ! SEA5 : solar elevation angle at this time, degree
2442 : ! JRatio : ratio J(OH)/J(NO2), unitless
2443 : ! NOx : concentration NOx in ambient air, ppt
2444 : ! WS : wind speed, m/s
2445 : !
2446 : ! In GEOS-Chem v9-01-03 through v9-02, the effects of wind speed on FNOx and OPE
2447 : ! were not included (wind speed set at 6 m/s). The JRatio also used J(O1D)
2448 : ! rather than J(OH); this has only a small effect on interpolated values.
2449 : ! To reproduce the behavior of these earlier versions, modify code below marked
2450 : ! with ******* and call READ\_PARANOX\_LUT\_v913 in emissions\_mod.F
2451 : !\\
2452 : !\\
2453 : ! !INTERFACE:
2454 : !
2455 0 : SUBROUTINE PARANOX_LUT( ExtState, HcoState, Inst, &
2456 : I, J, RC, FNOX, DNOx, OPE, MOE_OUT )
2457 : !
2458 : ! !USES:
2459 : !
2460 : USE HCO_STATE_MOD, ONLY : HCO_State
2461 : USE HCOX_STATE_MOD, ONLY : Ext_State
2462 : !
2463 : ! !INPUT PARAMETERS:
2464 : !
2465 : TYPE(Ext_State), POINTER :: ExtState
2466 : TYPE(HCO_State), POINTER :: HcoState
2467 : TYPE(MyInst), POINTER :: Inst
2468 : INTEGER, INTENT(IN) :: I, J ! Grid indices
2469 : !
2470 : ! !OUTPUT PARAMETERS:
2471 : !
2472 : REAL*8, INTENT(OUT) :: FNOX ! fraction of NOx remaining, mol/mol
2473 : REAL*8, INTENT(OUT) :: DNOX ! fraction of NOx deposited, mol/mol
2474 : REAL*8, INTENT(OUT) :: OPE ! net OPE, mol(net P(O3))/mol(P(HNO3))
2475 : REAL*8, INTENT(OUT), OPTIONAL :: MOE_OUT ! net MOE, mol(L(CH4))/mol(E(NOx))
2476 : !
2477 : ! !INPUT/OUTPUT PARAMETERS:
2478 : !
2479 : INTEGER, INTENT(INOUT) :: RC ! Return code
2480 : !
2481 : ! !REVISION HISTORY:
2482 : ! Jun 2010 - G.C.M. Vinken - Initial version
2483 : ! See https://github.com/geoschem/hemco for complete history
2484 : !EOP
2485 : !------------------------------------------------------------------------------
2486 : !BOC
2487 : !
2488 : ! !LOCAL VARIABLES:
2489 : !
2490 : INTEGER :: I1,I2,I3,I4,I5,I6,I7,I8
2491 : REAL(sp) :: FNOX_TMP, DNOX_TMP, OPE_TMP, MOE_TMP
2492 : REAL(sp) :: WEIGHT
2493 : REAL(sp) :: JNO2, JOH, TAIR
2494 : REAL(sp) :: AIR
2495 : REAL*8 :: MOE
2496 :
2497 : ! Interpolation variables, indices, and weights
2498 : REAL(sp), DIMENSION(8) :: VARS
2499 : INTEGER, DIMENSION(8,2) :: INDX
2500 : REAL(sp), DIMENSION(8,2) :: WTS
2501 :
2502 0 : REAL(sp), POINTER :: FRACNOX_LUT(:,:,:,:,:,:,:)
2503 0 : REAL(sp), POINTER :: DNOX_LUT (:,:,:,:,:,:,:)
2504 0 : REAL(sp), POINTER :: OPE_LUT (:,:,:,:,:,:,:)
2505 0 : REAL(sp), POINTER :: MOE_LUT (:,:,:,:,:,:,:)
2506 :
2507 : CHARACTER(LEN=255) :: MSG
2508 : CHARACTER(LEN=255) :: LOC = 'PARANOX_LUT'
2509 :
2510 : !=================================================================
2511 : ! PARANOX_LUT begins here!
2512 : !=================================================================
2513 :
2514 : ! Initialize for safety's sake
2515 0 : FRACNOX_LUT => NULL()
2516 0 : DNOX_LUT => NULL()
2517 0 : OPE_LUT => NULL()
2518 0 : MOE_LUT => NULL()
2519 0 : FNOX = 0.0_sp
2520 0 : DNOX = 0.0_sp
2521 0 : OPE = 0.0_sp
2522 0 : IF ( PRESENT( MOE_OUT ) ) THEN
2523 0 : MOE_OUT = 0.0_sp
2524 : ENDIF
2525 :
2526 : ! Air mass [kg]
2527 0 : AIR = ExtState%AIR%Arr%Val(I,J,1)
2528 :
2529 : ! Air temperature, K
2530 0 : Tair = ExtState%T2M%Arr%Val(I,J)
2531 :
2532 : ! ! for debugging only
2533 : ! if(I==3.and.J==35)then
2534 : ! write(*,*) 'Call PARANOX_LUT @ ',I,J
2535 : ! write(*,*) 'Tair: ', Tair
2536 : ! write(*,*) 'SUNCOSmid: ', SC5(I,J)
2537 : ! endif
2538 :
2539 : ! Check if sun is up
2540 0 : IF ( ExtState%SUNCOS%Arr%Val(I,J) > 0.0_hp ) THEN
2541 :
2542 : ! J(NO2), 1/s
2543 0 : JNO2 = ExtState%JNO2%Arr%Val(I,J)
2544 :
2545 : ! ! J(O1D), 1/s
2546 : ! JO1D = ExtState%JO1D%Arr%Val(I,J)
2547 : !
2548 : ! ! H2O, molec/cm3. Get from specific humidity, which is in kg/kg.
2549 : ! ! NOTE: SPHU is mass H2O / mass total air so use of dry air molecular
2550 : ! ! weight is slightly inaccurate. C (ewl, 9/11/15)
2551 : ! H2O = ExtState%SPHU%Arr%Val(I,J,1) * DENS &
2552 : ! * HcoState%Phys%AIRMW / MWH2O
2553 : !
2554 : ! ! Calculate J(OH), the effective rate for O3+hv -> OH+OH,
2555 : ! ! assuming steady state for O(1D).
2556 : ! ! Rate coefficients are cm3/molec/s; concentrations are molec/cm3
2557 : ! ! This should match the O3+hv (+H2O) -> OH+OH kinetics in calcrate.F
2558 : ! JOH = JO1D * &
2559 : ! 1.63e-10 * EXP( 60.e0/Tair) * H2O / &
2560 : ! ( 1.63e-10 * EXP( 60.e0/Tair) * H2O + &
2561 : ! 1.20e-10 * DENS * 0.5000e-6 + &
2562 : ! 2.15e-11 * EXP(110.e0/Tair) * DENS * 0.7808e0 + &
2563 : ! 3.30e-11 * EXP( 55.e0/Tair) * DENS * 0.2095e0 )
2564 :
2565 : ! J(OH) - effective rate for O3+hv(+H2O)-> OH+OH, 1/s
2566 0 : JOH = ExtState%JOH%Arr%Val(I,J)
2567 :
2568 : ELSE
2569 :
2570 : ! J-values are zero when sun is down
2571 : JNO2 = 0e0_sp
2572 : JOH = 0e0_sp
2573 :
2574 : ENDIF
2575 :
2576 : ! ! for debugging only
2577 : ! if(I==3.and.J==35)then
2578 : ! write(*,*) 'JNO2: ', JNO2
2579 : ! write(*,*) 'JOH : ', JOH
2580 : ! endif
2581 :
2582 : !========================================================================
2583 : ! Load all variables into a single array
2584 : !========================================================================
2585 :
2586 : ! Temperature, K
2587 0 : VARS(1) = Tair
2588 :
2589 : ! J(NO2), 1/s
2590 0 : VARS(2) = JNO2
2591 :
2592 : ! old
2593 : ! ! O3 concentration in ambient air, ppb
2594 : ! VARS(3) = ExtState%O3%Arr%Val(I,J,1) / AIR &
2595 : ! * HcoState%Phys%AIRMW / MW_O3 &
2596 : ! * 1.e9_sp
2597 : ! new
2598 : ! O3 concentration in ambient air, ppb
2599 : ! NOTE: ExtState%O3 units are now kg/kg dry air (ewl, 9/11/15)
2600 0 : VARS(3) = ExtState%O3%Arr%Val(I,J,1) &
2601 : * HcoState%Phys%AIRMW / Inst%MW_O3 &
2602 0 : * 1.e9_sp
2603 : ! end new (ewl)
2604 :
2605 : ! Solar elevation angle, degree
2606 : ! SEA0 = SEA when emitted from ship, 5-h before current model time
2607 : ! SEA5 = SEA at current model time, 5-h after emission from ship
2608 : ! Note: Since SEA = 90 - SZA, then cos(SZA) = sin(SEA) and
2609 : ! thus SEA = arcsin( cos( SZA ) )
2610 : !VARS(4) = ASIND( SC5(I,J) )
2611 : !VARS(5) = ASIND( ExtState%SUNCOS%Arr%Val(I,J) )
2612 0 : VARS(4) = ASIN( Inst%SC5(I,J) ) / HcoState%Phys%PI_180
2613 0 : VARS(5) = ASIN( ExtState%SUNCOS%Arr%Val(I,J) ) / HcoState%Phys%PI_180
2614 :
2615 : ! J(OH)/J(NO2), unitless
2616 : ! Note J(OH) is the loss rate (1/s) of O3 to OH, which accounts for
2617 : ! the temperature, pressure and water vapor dependence of these reactions
2618 0 : VARS(6) = 0.0_sp
2619 0 : IF ( JNO2 /= 0.0_sp ) THEN
2620 0 : IF ( (EXPONENT(JOH)-EXPONENT(JNO2)) < MAXEXPONENT(JOH) ) THEN
2621 0 : VARS(6) = JOH / JNO2
2622 : ENDIF
2623 : ENDIF
2624 :
2625 : ! old
2626 : ! ! NOx concetration in ambient air, ppt
2627 : ! VARS(7) = ( ( ExtState%NO%Arr%Val(I,J,1) / AIR &
2628 : ! * HcoState%Phys%AIRMW / MW_NO ) &
2629 : ! + ( ExtState%NO2%Arr%Val(I,J,1) / AIR &
2630 : ! * HcoState%Phys%AIRMW / MW_NO2 ) ) &
2631 : ! * 1.e12_sp
2632 : ! new
2633 : ! NOx concetration in ambient air, ppt
2634 : ! NOTE: ExtState vars NO and NO2 units are now kg/kg dry air (ewl, 9/11/15)
2635 0 : VARS(7) = ( ( ExtState%NO%Arr%Val(I,J,1) &
2636 : * HcoState%Phys%AIRMW / Inst%MW_NO ) &
2637 0 : + ( ExtState%NO2%Arr%Val(I,J,1) &
2638 : * HcoState%Phys%AIRMW / Inst%MW_NO2 ) ) &
2639 0 : * 1.e12_sp
2640 : ! end new (ewl)
2641 :
2642 : ! Wind speed, m/s
2643 0 : VARS(8) = SQRT( ExtState%U10M%Arr%Val(I,J)**2 &
2644 0 : + ExtState%V10M%Arr%Val(I,J)**2 )
2645 :
2646 : ! ! for debugging only
2647 : ! if(I==1.and.J==35)then
2648 : ! write(*,*) 'VARS(1) : ', VARS(1)
2649 : ! write(*,*) 'VARS(2) : ', VARS(2)
2650 : ! write(*,*) 'VARS(3) : ', VARS(3)
2651 : ! write(*,*) 'VARS(4) : ', VARS(4)
2652 : ! write(*,*) 'VARS(5) : ', VARS(5)
2653 : ! write(*,*) 'VARS(6) : ', VARS(6)
2654 : ! write(*,*) 'VARS(7) : ', VARS(7)
2655 : ! write(*,*) 'VARS(8) : ', VARS(8)
2656 : ! write(*,*) 'AIR : ', ExtState%AIR%Arr%Val(I,J,1)
2657 : ! write(*,*) 'AIRMW : ', HcoState%Phys%AIRMW
2658 : ! write(*,*) 'MWNO, NO2, O3: ', Inst%MW_NO, Inst%MW_NO2, Inst%MW_O3
2659 : ! write(*,*) 'O3conc : ', ExtState%O3%Arr%Val(I,J,1)
2660 : ! write(*,*) 'NO,NO2 conc : ', ExtState%NO%Arr%Val(I,J,1), &
2661 : ! ExtState%NO2%Arr%Val(I,J,1)
2662 : ! write(*,*) 'U, V : ', ExtState%U10M%Arr%Val(I,J), &
2663 : ! ExtState%V10M%Arr%Val(I,J)
2664 : ! endif
2665 :
2666 : !*****************************************************
2667 : ! Restoring the following lines reproduces the behavior of
2668 : ! GEOS-Chem v9-01-03 through v9-02
2669 : ! the LUT was indexed with the ratio J(O1D)/J(NO2) (cdh, 3/27/2014)
2670 : ! VARS(6) = SAFE_DIV( JO1D, JNO2, 0D0 )
2671 : ! JRATIOlev = (/ 5.e-4, 0.0015, 0.0025, 0.0055 /)
2672 : !*****************************************************
2673 :
2674 : !========================================================================
2675 : ! Find the indices of nodes and their corresponding weights for the
2676 : ! interpolation
2677 : !========================================================================
2678 :
2679 : ! Temperature:
2680 0 : CALL INTERPOL_LINWEIGHTS( Inst%Tlev, VARS(1), INDX(1,:), WTS(1,:) )
2681 :
2682 : ! J(NO2):
2683 0 : CALL INTERPOL_LINWEIGHTS( Inst%JNO2lev, VARS(2), INDX(2,:), WTS(2,:) )
2684 :
2685 : ! [O3]:
2686 0 : CALL INTERPOL_LINWEIGHTS( Inst%O3lev, VARS(3), INDX(3,:), WTS(3,:) )
2687 :
2688 : ! SEA0:
2689 0 : CALL INTERPOL_LINWEIGHTS( Inst%SEA0lev, VARS(4), INDX(4,:), WTS(4,:) )
2690 :
2691 : ! SEA5:
2692 0 : CALL INTERPOL_LINWEIGHTS( Inst%SEA5lev, VARS(5), INDX(5,:), WTS(5,:) )
2693 :
2694 : ! JRATIO:
2695 0 : CALL INTERPOL_LINWEIGHTS( Inst%JRATIOlev, VARS(6), INDX(6,:), WTS(6,:))
2696 :
2697 : ! [NOx]:
2698 0 : CALL INTERPOL_LINWEIGHTS( Inst%NOXlev, VARS(7), INDX(7,:), WTS(7,:) )
2699 :
2700 : ! Wind speed:
2701 0 : CALL INTERPOL_LINWEIGHTS( Inst%WSlev, VARS(8), INDX(8,:), WTS(8,:) )
2702 :
2703 : !========================================================================
2704 : ! Piecewise linear interpolation
2705 : !========================================================================
2706 :
2707 : ! Initialize
2708 0 : FNOX = 0.0d0
2709 0 : DNOx = 0.0d0
2710 0 : OPE = 0.0d0
2711 0 : MOE = 0.0d0
2712 :
2713 : ! Loop over wind speed
2714 0 : DO I8=1,2
2715 :
2716 : ! Point at the LUT for this wind speed
2717 : ! Last two digits in fortran variable names indicate wind speed in m/s
2718 0 : SELECT CASE ( NINT( Inst%WSlev(INDX(8,I8)) ) )
2719 : CASE ( 2 )
2720 0 : FRACNOX_LUT => Inst%FRACNOX_LUT02
2721 0 : DNOx_LUT => Inst%DNOx_LUT02
2722 0 : OPE_LUT => Inst%OPE_LUT02
2723 0 : MOE_LUT => Inst%MOE_LUT02
2724 : CASE ( 6 )
2725 0 : FRACNOX_LUT => Inst%FRACNOX_LUT06
2726 0 : DNOx_LUT => Inst%DNOx_LUT06
2727 0 : OPE_LUT => Inst%OPE_LUT06
2728 0 : MOE_LUT => Inst%MOE_LUT06
2729 : CASE ( 10 )
2730 0 : FRACNOX_LUT => Inst%FRACNOX_LUT10
2731 0 : DNOx_LUT => Inst%DNOx_LUT10
2732 0 : OPE_LUT => Inst%OPE_LUT10
2733 0 : MOE_LUT => Inst%MOE_LUT10
2734 : CASE ( 14 )
2735 0 : FRACNOX_LUT => Inst%FRACNOX_LUT14
2736 0 : DNOx_LUT => Inst%DNOx_LUT14
2737 0 : OPE_LUT => Inst%OPE_LUT14
2738 0 : MOE_LUT => Inst%MOE_LUT14
2739 : CASE ( 18 )
2740 0 : FRACNOX_LUT => Inst%FRACNOX_LUT18
2741 0 : DNOx_LUT => Inst%DNOx_LUT18
2742 0 : OPE_LUT => Inst%OPE_LUT18
2743 0 : MOE_LUT => Inst%MOE_LUT18
2744 : CASE DEFAULT
2745 0 : MSG = 'LUT error: Wind speed interpolation error!'
2746 0 : CALL HCO_ERROR(MSG, RC, THISLOC=LOC )
2747 0 : RETURN
2748 : END SELECT
2749 :
2750 : !*****************************************************
2751 : ! Restoring the following lines reproduces the behavior of
2752 : ! GEOS-Chem v9-01-03 through v9-02 in which wind speed
2753 : ! effects on FNOx and OPE are neglected (cdh, 3/25/2014)
2754 : ! FRACNOX_LUT => FRACNOX_LUT_v913
2755 : ! OPE_LUT => OPE_LUT_v913
2756 : !*****************************************************
2757 :
2758 : ! loop over all other variables
2759 0 : DO I7=1,2
2760 0 : DO I6=1,2
2761 0 : DO I5=1,2
2762 0 : DO I4=1,2
2763 0 : DO I3=1,2
2764 0 : DO I2=1,2
2765 0 : DO I1=1,2
2766 :
2767 : !------------------------------------------
2768 : ! Nodes and weights used in the interpolation
2769 : !------------------------------------------
2770 :
2771 : ! Fraction NOx from the LUT
2772 0 : FNOX_TMP = FRACNOX_LUT( INDX(1,I1), INDX(2,I2), INDX(3,I3), &
2773 0 : INDX(4,I4), INDX(5,I5), &
2774 0 : INDX(6,I6), INDX(7,I7) )
2775 :
2776 0 : DNOX_TMP = DNOx_LUT( INDX(1,I1), INDX(2,I2), INDX(3,I3), &
2777 : INDX(4,I4), INDX(5,I5), &
2778 0 : INDX(6,I6), INDX(7,I7) )
2779 :
2780 : ! OPE from the LUT
2781 0 : OPE_TMP = OPE_LUT( INDX(1,I1), INDX(2,I2), INDX(3,I3), &
2782 : INDX(4,I4), INDX(5,I5), &
2783 0 : INDX(6,I6), INDX(7,I7) )
2784 :
2785 : ! MOE from the LUT
2786 0 : MOE_TMP = MOE_LUT( INDX(1,I1), INDX(2,I2), INDX(3,I3), &
2787 : INDX(4,I4), INDX(5,I5), &
2788 0 : INDX(6,I6), INDX(7,I7) )
2789 :
2790 : ! Interpolation weight for this element
2791 : WEIGHT = WTS(1,I1) * WTS(2,I2) * WTS(3,I3) * WTS(4,I4) * &
2792 0 : WTS(5,I5) * WTS(6,I6) * WTS(7,I7) * WTS(8,I8)
2793 :
2794 : !-----------------------------------
2795 : ! Error Check
2796 : !-----------------------------------
2797 :
2798 : !IF ENCOUNTER -999 IN THE LUT PRINT ERROR!!
2799 0 : IF ( ( FNOX_TMP < 0. ) .or. ( FNOX_TMP > 1. ) ) THEN
2800 :
2801 0 : PRINT*, 'PARANOX_LUT: fracnox = ', FNOX_TMP
2802 0 : PRINT*, 'This occured at grid box ', I, J
2803 0 : PRINT*, 'Lon/Lat: ', HcoState%Grid%XMID%Val(I,J), HcoState%Grid%YMID%Val(I,J)
2804 0 : PRINT*, 'SZA 5 hours ago : ', VARS(4)
2805 0 : PRINT*, 'SZA at this time: ', VARS(5)
2806 0 : PRINT*, 'The two SZAs should not be more than 75 deg apart!'
2807 0 : PRINT*, 'If they are, your restart SZA might be wrong.'
2808 0 : PRINT*, 'You can try to coldstart PARANOx by commenting'
2809 0 : PRINT*, 'all PARANOX_SUNCOS entries in your HEMCO'
2810 0 : PRINT*, 'configuration file.'
2811 :
2812 : !print*, I1, I2, I3, I4, I5, I6, I7, I8
2813 : !print*, INDX(1,I1), INDX(2,I2), INDX(3,I3), INDX(4,I4), &
2814 : ! INDX(5,I5), INDX(6,I6), INDX(7,I7), INDX(8,I8)
2815 : !print*, VARS
2816 :
2817 0 : MSG = 'LUT error: Fracnox should be between 0 and 1!'
2818 0 : CALL HCO_ERROR(MSG, RC, THISLOC=LOC )
2819 0 : RETURN
2820 : ENDIF
2821 :
2822 : !-----------------------------------
2823 : ! Final interpolated values
2824 : !-----------------------------------
2825 :
2826 : ! Weighted sum of FNOx from the LUT
2827 0 : FNOx = FNOx + FNOX_TMP * WEIGHT
2828 :
2829 : ! Weighted sum of DNOx from the LUT
2830 0 : DNOx = DNOx + DNOX_TMP * WEIGHT
2831 :
2832 : ! Weighted sum of OPE from the LUT
2833 0 : OPE = OPE + OPE_TMP * WEIGHT
2834 :
2835 : ! Weighted sum of MOE from the LUT
2836 0 : MOE = MOE + MOE_TMP * WEIGHT
2837 :
2838 : END DO
2839 : END DO
2840 : END DO
2841 : END DO
2842 : END DO
2843 : END DO
2844 : END DO
2845 :
2846 : ! Free pointers
2847 0 : FRACNOX_LUT => NULL()
2848 0 : DNOx_LUT => NULL()
2849 0 : OPE_LUT => NULL()
2850 0 : MOE_LUT => NULL()
2851 : END DO
2852 :
2853 : ! Transfer MOE if optional output parameter is present
2854 0 : IF ( PRESENT( MOE_OUT ) ) MOE_OUT = MOE
2855 :
2856 : ! Nullify pointers
2857 0 : NULLIFY( FRACNOX_LUT )
2858 0 : NULLIFY( DNOx_LUT )
2859 0 : NULLIFY( OPE_LUT )
2860 0 : NULLIFY( MOE_LUT )
2861 :
2862 : ! Return w/ success
2863 0 : RC = HCO_SUCCESS
2864 :
2865 0 : END SUBROUTINE PARANOX_LUT
2866 : !EOC
2867 : !------------------------------------------------------------------------------
2868 : ! Harmonized Emissions Component (HEMCO) !
2869 : !------------------------------------------------------------------------------
2870 : !BOP
2871 : !
2872 : ! !IROUTINE: InstGet
2873 : !
2874 : ! !DESCRIPTION: Subroutine InstGet returns a poiner to the desired instance.
2875 : !\\
2876 : !\\
2877 : ! !INTERFACE:
2878 : !
2879 0 : SUBROUTINE InstGet ( Instance, Inst, RC, PrevInst )
2880 : !
2881 : ! !INPUT PARAMETERS:
2882 : !
2883 : INTEGER :: Instance
2884 : TYPE(MyInst), POINTER :: Inst
2885 : INTEGER :: RC
2886 : TYPE(MyInst), POINTER, OPTIONAL :: PrevInst
2887 : !
2888 : ! !REVISION HISTORY:
2889 : ! 18 Feb 2016 - C. Keller - Initial version
2890 : ! See https://github.com/geoschem/hemco for complete history
2891 : !EOP
2892 : !------------------------------------------------------------------------------
2893 : !BOC
2894 : TYPE(MyInst), POINTER :: PrvInst
2895 :
2896 : !=================================================================
2897 : ! InstGet begins here!
2898 : !=================================================================
2899 :
2900 : ! Get instance. Also archive previous instance.
2901 0 : PrvInst => NULL()
2902 0 : Inst => AllInst
2903 0 : DO WHILE ( ASSOCIATED(Inst) )
2904 0 : IF ( Inst%Instance == Instance ) EXIT
2905 0 : PrvInst => Inst
2906 0 : Inst => Inst%NextInst
2907 : END DO
2908 0 : IF ( .NOT. ASSOCIATED( Inst ) ) THEN
2909 0 : RC = HCO_FAIL
2910 0 : RETURN
2911 : ENDIF
2912 :
2913 : ! Pass output arguments
2914 0 : IF ( PRESENT(PrevInst) ) PrevInst => PrvInst
2915 :
2916 : ! Cleanup & Return
2917 0 : PrvInst => NULL()
2918 0 : RC = HCO_SUCCESS
2919 :
2920 : END SUBROUTINE InstGet
2921 : !EOC
2922 : !------------------------------------------------------------------------------
2923 : ! Harmonized Emissions Component (HEMCO) !
2924 : !------------------------------------------------------------------------------
2925 : !BOP
2926 : !
2927 : ! !IROUTINE: InstCreate
2928 : !
2929 : ! !DESCRIPTION: Subroutine InstCreate creates a new instance.
2930 : !\\
2931 : !\\
2932 : ! !INTERFACE:
2933 : !
2934 0 : SUBROUTINE InstCreate ( ExtNr, Instance, Inst, RC )
2935 : !
2936 : ! !INPUT PARAMETERS:
2937 : !
2938 : INTEGER, INTENT(IN) :: ExtNr
2939 : !
2940 : ! !OUTPUT PARAMETERS:
2941 : !
2942 : INTEGER, INTENT( OUT) :: Instance
2943 : TYPE(MyInst), POINTER :: Inst
2944 : !
2945 : ! !INPUT/OUTPUT PARAMETERS:
2946 : !
2947 : INTEGER, INTENT(INOUT) :: RC
2948 : !
2949 : ! !REVISION HISTORY:
2950 : ! 18 Feb 2016 - C. Keller - Initial version
2951 : ! See https://github.com/geoschem/hemco for complete history
2952 : !EOP
2953 : !------------------------------------------------------------------------------
2954 : !BOC
2955 : TYPE(MyInst), POINTER :: TmpInst
2956 : INTEGER :: nnInst
2957 :
2958 : !=================================================================
2959 : ! InstCreate begins here!
2960 : !=================================================================
2961 :
2962 : ! ----------------------------------------------------------------
2963 : ! Generic instance initialization
2964 : ! ----------------------------------------------------------------
2965 :
2966 : ! Initialize
2967 0 : Inst => NULL()
2968 :
2969 : ! Get number of already existing instances
2970 0 : TmpInst => AllInst
2971 0 : nnInst = 0
2972 0 : DO WHILE ( ASSOCIATED(TmpInst) )
2973 0 : nnInst = nnInst + 1
2974 0 : TmpInst => TmpInst%NextInst
2975 : END DO
2976 :
2977 : ! Create new instance
2978 0 : ALLOCATE(Inst)
2979 0 : Inst%Instance = nnInst + 1
2980 0 : Inst%ExtNr = ExtNr
2981 :
2982 : ! Attach to instance list
2983 0 : Inst%NextInst => AllInst
2984 0 : AllInst => Inst
2985 :
2986 : ! Update output instance
2987 0 : Instance = Inst%Instance
2988 :
2989 : ! ----------------------------------------------------------------
2990 : ! Type specific initialization statements follow below
2991 : ! ----------------------------------------------------------------
2992 :
2993 : ! Return w/ success
2994 0 : RC = HCO_SUCCESS
2995 :
2996 0 : END SUBROUTINE InstCreate
2997 : !EOC
2998 : !------------------------------------------------------------------------------
2999 : ! Harmonized Emissions Component (HEMCO) !
3000 : !------------------------------------------------------------------------------
3001 : !BOP
3002 : !
3003 : ! !IROUTINE: InstRemove
3004 : !
3005 : ! !DESCRIPTION: Subroutine InstRemove creates a new instance.
3006 : !\\
3007 : !\\
3008 : ! !INTERFACE:
3009 : !
3010 0 : SUBROUTINE InstRemove ( Instance )
3011 : !
3012 : ! !INPUT PARAMETERS:
3013 : !
3014 : INTEGER :: Instance
3015 : !
3016 : ! !REVISION HISTORY:
3017 : ! 18 Feb 2016 - C. Keller - Initial version
3018 : ! See https://github.com/geoschem/hemco for complete history
3019 : !EOP
3020 : !------------------------------------------------------------------------------
3021 : !BOC
3022 : INTEGER :: RC
3023 : TYPE(MyInst), POINTER :: PrevInst
3024 : TYPE(MyInst), POINTER :: Inst
3025 :
3026 : !=================================================================
3027 : ! InstRemove begins here!
3028 : !=================================================================
3029 :
3030 : ! Init
3031 0 : PrevInst => NULL()
3032 0 : Inst => NULL()
3033 :
3034 : ! Get instance. Also archive previous instance.
3035 0 : CALL InstGet ( Instance, Inst, RC, PrevInst=PrevInst )
3036 :
3037 : ! Instance-specific deallocation
3038 0 : IF ( ASSOCIATED(Inst) ) THEN
3039 :
3040 : !---------------------------------------------------------------------
3041 : ! Deallocate fields of Inst before popping off from the list
3042 : ! in order to avoid memory leaks (Bob Yantosca (17 Aug 2022)
3043 : !---------------------------------------------------------------------
3044 0 : IF ( ASSOCIATED( Inst%ShipNO ) ) THEN
3045 0 : DEALLOCATE ( Inst%ShipNO )
3046 : ENDIF
3047 0 : Inst%ShipNO => NULL()
3048 :
3049 0 : IF ( ASSOCIATED( Inst%SC5 ) ) THEN
3050 0 : DEALLOCATE ( Inst%SC5 )
3051 : ENDIF
3052 0 : Inst%SC5 => NULL()
3053 :
3054 0 : IF ( ASSOCIATED( Inst%FRACNOX_LUT02 ) ) THEN
3055 0 : DEALLOCATE( Inst%FRACNOX_LUT02 )
3056 : ENDIF
3057 0 : Inst%FRACNOX_LUT02 => NULL()
3058 :
3059 0 : IF ( ASSOCIATED( Inst%FRACNOX_LUT06 ) ) THEN
3060 0 : DEALLOCATE( Inst%FRACNOX_LUT06 )
3061 : ENDIF
3062 0 : Inst%FRACNOX_LUT06 => NULL()
3063 :
3064 0 : IF ( ASSOCIATED( Inst%FRACNOX_LUT10 ) ) THEN
3065 0 : DEALLOCATE( Inst%FRACNOX_LUT10 )
3066 : ENDIF
3067 0 : Inst%FRACNOX_LUT10 => NULL()
3068 :
3069 0 : IF ( ASSOCIATED( Inst%FRACNOX_LUT14 ) ) THEN
3070 0 : DEALLOCATE( Inst%FRACNOX_LUT14 )
3071 : ENDIF
3072 0 : Inst%FRACNOX_LUT14 => NULL()
3073 :
3074 0 : IF ( ASSOCIATED( Inst%FRACNOX_LUT18 ) ) THEN
3075 0 : DEALLOCATE( Inst%FRACNOX_LUT18 )
3076 : ENDIF
3077 0 : Inst%FRACNOX_LUT18 => NULL()
3078 :
3079 0 : IF ( ASSOCIATED( Inst%OPE_LUT02 ) ) THEN
3080 0 : DEALLOCATE( Inst%OPE_LUT02 )
3081 : ENDIF
3082 0 : Inst%OPE_LUT02 => NULL()
3083 :
3084 0 : IF ( ASSOCIATED( Inst%OPE_LUT06 ) ) THEN
3085 0 : DEALLOCATE( Inst%OPE_LUT06 )
3086 : ENDIF
3087 0 : Inst%OPE_LUT06 => NULL()
3088 :
3089 0 : IF ( ASSOCIATED( Inst%OPE_LUT10 ) ) THEN
3090 0 : DEALLOCATE( Inst%OPE_LUT10 )
3091 : ENDIF
3092 0 : Inst%OPE_LUT10 => NULL()
3093 :
3094 0 : IF ( ASSOCIATED( Inst%OPE_LUT14 ) ) THEN
3095 0 : DEALLOCATE( Inst%OPE_LUT14 )
3096 : ENDIF
3097 0 : Inst%OPE_LUT14 => NULL()
3098 :
3099 0 : IF ( ASSOCIATED( Inst%OPE_LUT18 ) ) THEN
3100 0 : DEALLOCATE( Inst%OPE_LUT18 )
3101 : ENDIF
3102 0 : Inst%OPE_LUT18 => NULL()
3103 :
3104 0 : IF ( ASSOCIATED( Inst%MOE_LUT02 ) ) THEN
3105 0 : DEALLOCATE( Inst%MOE_LUT02 )
3106 : ENDIF
3107 0 : Inst%MOE_LUT02 => NULL()
3108 :
3109 0 : IF ( ASSOCIATED( Inst%MOE_LUT06 ) ) THEN
3110 0 : DEALLOCATE( Inst%MOE_LUT06 )
3111 : ENDIF
3112 0 : Inst%MOE_LUT06 => NULL()
3113 :
3114 0 : IF ( ASSOCIATED( Inst%MOE_LUT10 ) ) THEN
3115 0 : DEALLOCATE( Inst%MOE_LUT10 )
3116 : ENDIF
3117 0 : Inst%MOE_LUT10 => NULL()
3118 :
3119 0 : IF ( ASSOCIATED( Inst%MOE_LUT14 ) ) THEN
3120 0 : DEALLOCATE( Inst%MOE_LUT14 )
3121 : ENDIF
3122 0 : Inst%MOE_LUT14 => NULL()
3123 :
3124 0 : IF ( ASSOCIATED( Inst%MOE_LUT18 ) ) THEN
3125 0 : DEALLOCATE( Inst%MOE_LUT18 )
3126 : ENDIF
3127 0 : Inst%MOE_LUT18 => NULL()
3128 :
3129 0 : IF ( ASSOCIATED( Inst%DNOx_LUT02 ) ) THEN
3130 0 : DEALLOCATE( Inst%DNOx_LUT02 )
3131 : ENDIF
3132 0 : Inst%DNOx_LUT02 => NULL()
3133 :
3134 0 : IF ( ASSOCIATED( Inst%DNOx_LUT06 ) ) THEN
3135 0 : DEALLOCATE( Inst%DNOx_LUT06 )
3136 : ENDIF
3137 0 : Inst%DNOx_LUT06 => NULL()
3138 :
3139 0 : IF ( ASSOCIATED( Inst%DNOx_LUT10 ) ) THEN
3140 0 : DEALLOCATE( Inst%DNOx_LUT10 )
3141 : ENDIF
3142 0 : Inst%DNOx_LUT10 => NULL()
3143 :
3144 0 : IF ( ASSOCIATED( Inst%DNOx_LUT14 ) ) THEN
3145 0 : DEALLOCATE( Inst%DNOx_LUT14 )
3146 : ENDIF
3147 0 : Inst%DNOx_LUT14 => NULL()
3148 :
3149 0 : IF ( ASSOCIATED( Inst%DNOx_LUT18 ) ) THEN
3150 0 : DEALLOCATE( Inst%DNOx_LUT18 )
3151 : ENDIF
3152 0 : Inst%DNOx_LUT18 => NULL()
3153 :
3154 0 : IF ( ASSOCIATED( Inst%DEPO3 ) ) THEN
3155 0 : DEALLOCATE( Inst%DEPO3 )
3156 : ENDIF
3157 0 : Inst%DEPO3 => NULL()
3158 :
3159 0 : IF ( ASSOCIATED( Inst%DEPHNO3 ) ) THEN
3160 0 : DEALLOCATE( Inst%DEPHNO3 )
3161 : ENDIF
3162 0 : Inst%DEPHNO3 => NULL()
3163 :
3164 : !---------------------------------------------------------------------
3165 : ! Pop off instance from list
3166 : !---------------------------------------------------------------------
3167 0 : IF ( ASSOCIATED(PrevInst) ) THEN
3168 0 : PrevInst%NextInst => Inst%NextInst
3169 : ELSE
3170 0 : AllInst => Inst%NextInst
3171 : ENDIF
3172 0 : DEALLOCATE(Inst)
3173 : ENDIF
3174 :
3175 : ! Free pointers before exiting
3176 0 : PrevInst => NULL()
3177 0 : Inst => NULL()
3178 :
3179 0 : END SUBROUTINE InstRemove
3180 : !EOC
3181 0 : END MODULE HCOX_ParaNOx_mod
|
