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