Line data Source code
1 : module elevated_emissions_mod
2 : !---------------------------------------------------------------
3 : ! ... elevalted emissions module
4 : !---------------------------------------------------------------
5 :
6 : use shr_kind_mod, only : r8 => shr_kind_r8, shr_kind_cl
7 : use spmd_utils, only : masterproc
8 : use cam_abortutils,only : endrun
9 : use ioFileMod, only : getfil
10 : use cam_logfile, only : iulog
11 : use tracer_data, only : trfld,trfile
12 : use infnan, only : nan, assignment(=)
13 : use cam_history, only : addfld, outfld, add_default, fieldname_len
14 :
15 : implicit none
16 :
17 : type :: emission
18 : integer :: bufndx
19 : real(r8) :: scalefactor
20 : character(len=256):: filename
21 : character(len=16) :: species
22 : character(len=32) :: units
23 : integer :: nsectors
24 : character(len=32),pointer :: sectors(:)
25 : type(trfld), pointer :: fields(:)
26 : type(trfile) :: file
27 : end type emission
28 :
29 : private
30 :
31 : public :: elevated_emissions_readnl
32 : public :: elevated_emissions_reg
33 : public :: elevated_emissions_init
34 : public :: elevated_emissions_adv
35 : public :: elevated_emissions_set
36 :
37 : integer, parameter :: NMAX=50
38 :
39 : type(emission), allocatable :: elev_emis(:)
40 : integer :: n_emis_files = 0
41 : integer :: n_pbuf_flds = 0
42 :
43 : character(len=shr_kind_cl) :: elev_emis_specifier(NMAX) = ' '
44 : character(len=24) :: elev_emis_type
45 : integer :: elev_emis_cycle_yr
46 : integer :: elev_emis_fixed_ymd
47 : integer :: elev_emis_fixed_tod
48 :
49 : character(len=fieldname_len) :: names(NMAX) = ' '
50 : character(len=32) :: units(NMAX) = ' '
51 : integer :: indexes(NMAX) = -1
52 : integer :: n_diags = 0
53 :
54 : contains
55 :
56 : !-----------------------------------------------------------------------
57 : !-----------------------------------------------------------------------
58 1536 : subroutine elevated_emissions_readnl(nlfile)
59 :
60 : use namelist_utils, only : find_group_name
61 : use spmd_utils, only : mpicom, masterprocid, mpi_integer, mpi_character
62 :
63 : character(len=*), intent(in) :: nlfile ! filepath for file containing namelist input
64 :
65 : ! Local variables
66 : integer :: unitn, ierr, i
67 : logical :: logmsg
68 : character(len=*), parameter :: subname = 'elevated_emissions_readnl'
69 :
70 : namelist /elevated_emissions_opts/ elev_emis_specifier, elev_emis_type, elev_emis_cycle_yr, &
71 : elev_emis_fixed_ymd, elev_emis_fixed_tod
72 :
73 : ! Read namelist
74 1536 : if (masterproc) then
75 2 : open( newunit=unitn, file=trim(nlfile), status='old' )
76 2 : call find_group_name(unitn, 'elevated_emissions_opts', status=ierr)
77 2 : if (ierr == 0) then
78 0 : read(unitn, elevated_emissions_opts, iostat=ierr)
79 0 : if (ierr /= 0) then
80 0 : call endrun(subname // ':: ERROR reading namelist')
81 : end if
82 : end if
83 2 : close(unitn)
84 :
85 : end if
86 :
87 : ! Broadcast namelist variables
88 1536 : call mpi_bcast(elev_emis_specifier,len(elev_emis_specifier(1))*NMAX, mpi_character, masterprocid, mpicom, ierr)
89 1536 : call mpi_bcast(elev_emis_type, len(elev_emis_type), mpi_character, masterprocid, mpicom, ierr)
90 1536 : call mpi_bcast(elev_emis_cycle_yr, 1, mpi_integer, masterprocid, mpicom, ierr)
91 1536 : call mpi_bcast(elev_emis_fixed_ymd, 1, mpi_integer, masterprocid, mpicom, ierr)
92 1536 : call mpi_bcast(elev_emis_fixed_tod, 1, mpi_integer, masterprocid, mpicom, ierr)
93 :
94 1536 : logmsg = .false.
95 1536 : if (masterproc) then
96 102 : do i = 1,NMAX
97 102 : if (len_trim(elev_emis_specifier(i))>0) then
98 0 : logmsg = .true.
99 0 : write(iulog,'(2a)') subname,': elev_emis_specifier: ',trim(elev_emis_specifier(i))
100 : endif
101 : enddo
102 2 : if (logmsg) then
103 0 : write(iulog,*) subname,': elev_emis_type: ',elev_emis_type
104 0 : write(iulog,*) subname,': elev_emis_cycle_yr: ',elev_emis_cycle_yr
105 0 : write(iulog,*) subname,': elev_emis_fixed_ymd: ',elev_emis_fixed_ymd
106 0 : write(iulog,*) subname,': elev_emis_fixed_tod: ',elev_emis_fixed_tod
107 : endif
108 : endif
109 :
110 1536 : end subroutine elevated_emissions_readnl
111 :
112 : !-----------------------------------------------------------------------
113 : !-----------------------------------------------------------------------
114 1536 : subroutine elevated_emissions_reg( )
115 : use m_MergeSorts, only : IndexSort
116 : use physics_buffer, only : pbuf_add_field, dtype_r8, pbuf_get_index
117 : use ppgrid, only : pcols, pver
118 :
119 : !-----------------------------------------------------------------------
120 : ! ... local variables
121 : !-----------------------------------------------------------------------
122 : integer :: astat
123 : integer :: j, l, m, n, i, nn, kk ! Indices
124 : character(len=16) :: spc_name
125 : character(len=256) :: filename
126 :
127 : character(len=16) :: emis_species(NMAX)
128 : character(len=256) :: emis_filenam(NMAX)
129 : integer :: emis_indexes(NMAX)
130 : integer :: indx(NMAX)
131 : real(r8) :: emis_scalefactor(NMAX)
132 :
133 : character(len=256) :: tmp_string = ' '
134 : character(len=32) :: xchr = ' '
135 : real(r8) :: xdbl
136 :
137 :
138 : integer :: err
139 : character(len=32) :: bname
140 :
141 : character(len=*), parameter :: prefix = 'elevated_emissions_reg: '
142 1536 : kk = 0
143 1536 : nn = 0
144 1536 : indx(:) = 0
145 78336 : emis_species = ' '
146 78336 : emis_indexes = -1
147 78336 : emis_filenam = 'NONE'
148 :
149 1536 : count_emis: do n=1,size(elev_emis_specifier)
150 1536 : if ( len_trim(elev_emis_specifier(n) ) == 0 ) then
151 : exit count_emis
152 : endif
153 :
154 0 : i = scan(elev_emis_specifier(n),'->')
155 0 : spc_name = trim(adjustl(elev_emis_specifier(n)(:i-1)))
156 :
157 : ! need to parse out scalefactor ...
158 0 : tmp_string = adjustl(elev_emis_specifier(n)(i+2:))
159 0 : j = scan( tmp_string, '*' )
160 0 : if (j>0) then
161 0 : xchr = tmp_string(1:j-1) ! get the multipler (left of the '*')
162 0 : read( xchr, * ) xdbl ! convert the string to a real
163 0 : tmp_string = adjustl(tmp_string(j+1:)) ! get the filepath name (right of the '*')
164 : else
165 0 : xdbl = 1._r8
166 : endif
167 0 : filename = trim(tmp_string)
168 :
169 0 : bname = trim(spc_name)//'_elevemis'
170 :
171 0 : m = pbuf_get_index(bname,errcode=err)
172 0 : if (m<1) then
173 0 : call pbuf_add_field(bname, 'physpkg', dtype_r8, (/pcols,pver/), m)
174 0 : kk = kk+1
175 0 : names(kk) = bname
176 0 : indexes(kk) = m
177 : endif
178 :
179 0 : nn = nn+1
180 0 : emis_species(nn) = spc_name
181 0 : emis_filenam(nn) = filename
182 0 : emis_indexes(nn) = m
183 0 : emis_scalefactor(nn) = xdbl
184 :
185 1536 : indx(n)=n
186 : enddo count_emis
187 :
188 1536 : n_diags = kk
189 1536 : n_emis_files = nn
190 :
191 1536 : if (masterproc) write(iulog,*) prefix,' n_emis_files = ',n_emis_files
192 :
193 9216 : allocate( elev_emis(n_emis_files), stat=astat )
194 1536 : if( astat/= 0 ) then
195 0 : write(iulog,*) 'elev_emis_inti: failed to allocate emissions array; error = ',astat
196 0 : call endrun('elev_emis_inti: failed to allocate emissions array')
197 : end if
198 :
199 : !-----------------------------------------------------------------------
200 : ! Sort the input files so that the emissions sources are summed in the
201 : ! same order regardless of the order of the input files in the namelist
202 : !-----------------------------------------------------------------------
203 1536 : if (n_emis_files > 0) then
204 0 : call IndexSort(n_emis_files, indx, emis_filenam)
205 : end if
206 :
207 : !-----------------------------------------------------------------------
208 : ! ... setup the emission type array
209 : !-----------------------------------------------------------------------
210 1536 : do m=1,n_emis_files
211 0 : elev_emis(m)%bufndx = emis_indexes(indx(m))
212 0 : elev_emis(m)%species = emis_species(indx(m))
213 0 : elev_emis(m)%filename = emis_filenam(indx(m))
214 1536 : elev_emis(m)%scalefactor = emis_scalefactor(indx(m))
215 : enddo
216 1536 : end subroutine elevated_emissions_reg
217 :
218 : !-----------------------------------------------------------------------
219 : !-----------------------------------------------------------------------
220 1536 : subroutine elevated_emissions_init(pbuf2d)
221 1536 : use tracer_data, only : trcdata_init
222 : use cam_pio_utils, only : cam_pio_openfile
223 : use pio, only : pio_inquire, pio_nowrite, pio_closefile, pio_inq_varndims
224 : use pio, only : pio_inq_varname, pio_inq_vardimid, pio_inq_dimid
225 : use pio, only : file_desc_t, pio_get_att, PIO_NOERR, PIO_GLOBAL
226 : use pio, only : pio_seterrorhandling, PIO_BCAST_ERROR,PIO_INTERNAL_ERROR
227 : use string_utils, only : GLC
228 : use physics_buffer,only : physics_buffer_desc, pbuf_set_field
229 :
230 : !--------------------------------------------------------
231 : ! ... Dummy arguments
232 : !--------------------------------------------------------
233 : type(physics_buffer_desc), pointer :: pbuf2d(:,:)
234 :
235 : !-----------------------------------------------------------------------
236 : ! ... local variables
237 : !-----------------------------------------------------------------------
238 : integer :: ierr, astat, l, m, n
239 : logical :: unstructured
240 : integer :: vid, nvars, isec, num_dims_emis
241 : integer :: vndims
242 1536 : logical, allocatable :: is_sector(:)
243 : type(file_desc_t) :: ncid
244 : character(len=32) :: varname
245 : character(len=256) :: locfn
246 : character(len=80) :: file_interp_type = ' '
247 1536 : integer, allocatable :: dimids(:)
248 : integer :: time_dimid, ncol_dimid
249 :
250 : character(len=32) :: emis_type = ' '
251 : character(len=1), parameter :: filelist = ''
252 : character(len=1), parameter :: datapath = ''
253 : logical , parameter :: rmv_file = .false.
254 : real(r8) :: xnan
255 :
256 1536 : xnan = nan
257 : !-----------------------------------------------------------------------
258 : ! read emis files to determine number of sectors
259 : !-----------------------------------------------------------------------
260 1536 : files_loop: do m = 1, n_emis_files
261 :
262 0 : elev_emis(m)%nsectors = 0
263 0 : call getfil (elev_emis(m)%filename, locfn, 0)
264 0 : call cam_pio_openfile ( ncid, trim(locfn), PIO_NOWRITE)
265 0 : ierr = pio_inquire (ncid, nVariables=nvars)
266 :
267 0 : call pio_seterrorhandling(ncid, PIO_BCAST_ERROR)
268 0 : ierr = pio_inq_dimid( ncid, 'ncol', ncol_dimid )
269 0 : unstructured = ierr==PIO_NOERR
270 0 : call pio_seterrorhandling(ncid, PIO_INTERNAL_ERROR)
271 :
272 0 : allocate(is_sector(nvars))
273 0 : is_sector(:) = .false.
274 :
275 0 : if (unstructured) then
276 0 : ierr = pio_inq_dimid( ncid, 'time', time_dimid )
277 : end if
278 :
279 0 : do vid = 1,nvars
280 :
281 0 : ierr = pio_inq_varndims (ncid, vid, vndims)
282 :
283 0 : if (unstructured) then
284 : num_dims_emis = 3
285 : else
286 0 : num_dims_emis = 4
287 : endif
288 :
289 0 : if( vndims < num_dims_emis ) then
290 : cycle
291 0 : elseif( vndims > num_dims_emis ) then
292 0 : ierr = pio_inq_varname (ncid, vid, varname)
293 0 : write(iulog,*) 'elev_emis_inti: Skipping variable ', trim(varname),', ndims = ',vndims, &
294 0 : ' , species=',trim(elev_emis(m)%species)
295 0 : cycle
296 : end if
297 :
298 0 : if (unstructured) then
299 0 : allocate( dimids(vndims) )
300 0 : ierr = pio_inq_vardimid( ncid, vid, dimids )
301 0 : if ( any(dimids(:)==ncol_dimid) .and. any(dimids(:)==time_dimid) ) then
302 0 : elev_emis(m)%nsectors = elev_emis(m)%nsectors+1
303 0 : is_sector(vid)=.true.
304 : endif
305 0 : deallocate(dimids)
306 : else
307 0 : elev_emis(m)%nsectors = elev_emis(m)%nsectors+1
308 0 : is_sector(vid)=.true.
309 : end if
310 :
311 : enddo
312 :
313 0 : allocate( elev_emis(m)%sectors(elev_emis(m)%nsectors), stat=astat )
314 0 : if( astat/= 0 ) then
315 0 : write(iulog,*) 'elev_emis_inti: failed to allocate elev_emis(m)%sectors array; error = ',astat
316 0 : call endrun
317 : end if
318 :
319 0 : isec = 1
320 :
321 0 : do vid = 1,nvars
322 0 : if( is_sector(vid) ) then
323 0 : ierr = pio_inq_varname(ncid, vid, elev_emis(m)%sectors(isec))
324 0 : isec = isec+1
325 : endif
326 : enddo
327 0 : deallocate(is_sector)
328 :
329 : ! Global attribute 'input_method' overrides the srf_emis_type namelist setting on
330 : ! a file-by-file basis. If the emis file does not contain the 'input_method'
331 : ! attribute then the srf_emis_type namelist setting is used.
332 0 : call pio_seterrorhandling(ncid, PIO_BCAST_ERROR)
333 0 : ierr = pio_get_att(ncid, PIO_GLOBAL, 'input_method', file_interp_type)
334 0 : call pio_seterrorhandling(ncid, PIO_INTERNAL_ERROR)
335 0 : if ( ierr == PIO_NOERR) then
336 0 : l = GLC(file_interp_type)
337 0 : emis_type(1:l) = file_interp_type(1:l)
338 0 : emis_type(l+1:) = ' '
339 : else
340 0 : emis_type = trim(elev_emis_type)
341 : endif
342 :
343 0 : call pio_closefile (ncid)
344 :
345 0 : allocate(elev_emis(m)%file%in_pbuf(size(elev_emis(m)%sectors)))
346 0 : elev_emis(m)%file%in_pbuf(:) = .false.
347 :
348 0 : call trcdata_init( elev_emis(m)%sectors, &
349 : elev_emis(m)%filename, filelist, datapath, &
350 : elev_emis(m)%fields, &
351 : elev_emis(m)%file, &
352 : rmv_file, elev_emis_cycle_yr, &
353 0 : elev_emis_fixed_ymd, elev_emis_fixed_tod, trim(emis_type) )
354 :
355 0 : elev_emis(m)%units = elev_emis(m)%fields(1)%units
356 :
357 0 : call pbuf_set_field(pbuf2d, elev_emis(m)%bufndx, xnan)
358 :
359 1536 : set_units: do n = 1,n_diags
360 0 : if (trim(elev_emis(m)%species)//'_elevemis'==names(n)) then
361 0 : units(n) = elev_emis(m)%fields(1)%units
362 0 : exit set_units
363 : end if
364 : end do set_units
365 :
366 : enddo files_loop
367 :
368 1536 : do n = 1, n_diags
369 1536 : call addfld(names(n), (/ 'lev' /), 'A', units(n), 'pbuf elev emis '//trim(names(n)))
370 : end do
371 :
372 3072 : end subroutine elevated_emissions_init
373 :
374 : !-----------------------------------------------------------------------
375 : !-----------------------------------------------------------------------
376 32256 : subroutine elevated_emissions_adv( pbuf2d, state )
377 : !-----------------------------------------------------------------------
378 : ! ... check serial case for time span
379 : !-----------------------------------------------------------------------
380 :
381 1536 : use physics_types,only : physics_state
382 : use ppgrid, only : begchunk, endchunk
383 : use tracer_data, only : advance_trcdata
384 : use physics_buffer, only : physics_buffer_desc, pbuf_set_field
385 :
386 : !--------------------------------------------------------
387 : ! ... Dummy arguments
388 : !--------------------------------------------------------
389 : type(physics_state), intent(in):: state(begchunk:endchunk)
390 : type(physics_buffer_desc), pointer :: pbuf2d(:,:)
391 :
392 : !-----------------------------------------------------------------------
393 : ! ... local variables
394 : !-----------------------------------------------------------------------
395 : integer :: m
396 :
397 16128 : do m = 1,n_emis_files
398 0 : call advance_trcdata( elev_emis(m)%fields, elev_emis(m)%file, state, pbuf2d )
399 16128 : call pbuf_set_field(pbuf2d, elev_emis(m)%bufndx, 0._r8)
400 : end do
401 :
402 16128 : end subroutine elevated_emissions_adv
403 :
404 : !-----------------------------------------------------------------------
405 : !-----------------------------------------------------------------------
406 80640 : subroutine elevated_emissions_set( lchnk, ncol, pbuf )
407 16128 : use physics_buffer, only : physics_buffer_desc, pbuf_get_field
408 :
409 : !--------------------------------------------------------
410 : ! ... Dummy arguments
411 : !--------------------------------------------------------
412 : integer, intent(in) :: ncol
413 : integer, intent(in) :: lchnk
414 : type(physics_buffer_desc), pointer :: pbuf(:)
415 :
416 : !--------------------------------------------------------
417 : ! ... local variables
418 : !--------------------------------------------------------
419 : integer :: isec, m, n
420 80640 : real(r8), pointer :: flux(:,:)
421 :
422 : !--------------------------------------------------------
423 : ! ... set non-zero emissions
424 : !--------------------------------------------------------
425 80640 : do m = 1,n_emis_files
426 0 : call pbuf_get_field(pbuf, elev_emis(m)%bufndx, flux)
427 80640 : do isec = 1,elev_emis(m)%nsectors
428 0 : flux(:ncol,:) = flux(:ncol,:) + elev_emis(m)%scalefactor*elev_emis(m)%fields(isec)%data(:ncol,:,lchnk)
429 : enddo
430 : end do
431 :
432 80640 : do n = 1, n_diags
433 0 : call pbuf_get_field(pbuf, indexes(n), flux)
434 80640 : call outfld(names(n), flux(:ncol,:), ncol, lchnk)
435 : end do
436 :
437 161280 : end subroutine elevated_emissions_set
438 :
439 80640 : end module elevated_emissions_mod
|
