Line data Source code
1 : !-------------------------------------------------------------------
2 : ! manages reading and interpolation of prescribed ghg tracers
3 : ! Created by: Francis Vitt
4 : !-------------------------------------------------------------------
5 : module prescribed_ghg
6 :
7 : use shr_kind_mod, only : r8 => shr_kind_r8
8 : use cam_abortutils, only : endrun
9 : use spmd_utils, only : masterproc
10 : use tracer_data, only : trfld, trfile
11 : use cam_logfile, only : iulog
12 :
13 : implicit none
14 : private
15 : save
16 :
17 : type(trfld), pointer :: fields(:)
18 : type(trfile) :: file
19 :
20 : public :: prescribed_ghg_init
21 : public :: prescribed_ghg_adv
22 : public :: write_prescribed_ghg_restart
23 : public :: read_prescribed_ghg_restart
24 : public :: has_prescribed_ghg
25 : public :: prescribed_ghg_register
26 : public :: init_prescribed_ghg_restart
27 : public :: prescribed_ghg_readnl
28 :
29 : logical :: has_prescribed_ghg = .false.
30 : integer, parameter, public :: N_GHG = 5
31 : integer :: number_flds
32 :
33 : character(len=256) :: filename = 'NONE'
34 : character(len=256) :: filelist = ''
35 : character(len=256) :: datapath = ''
36 : character(len=32) :: datatype = 'SERIAL'
37 : logical :: rmv_file = .false.
38 : integer :: cycle_yr = 0
39 : integer :: fixed_ymd = 0
40 : integer :: fixed_tod = 0
41 : character(len=16) :: specifier(N_GHG) = ''
42 :
43 : character(len=8) :: ghg_names(N_GHG) = (/ 'prsd_co2', 'prsd_ch4', 'prsd_n2o', 'prsd_f11', 'prsd_f12' /)
44 : real(r8), parameter :: molmass(N_GHG) = (/ 44.00980_r8, 16.04060_r8, 44.01288_r8, 137.3675_r8, 120.9132_r8 /)
45 :
46 : integer :: index_map(N_GHG)
47 :
48 : contains
49 :
50 :
51 : !-------------------------------------------------------------------
52 : !-------------------------------------------------------------------
53 1536 : subroutine prescribed_ghg_readnl(nlfile)
54 :
55 : use namelist_utils, only: find_group_name
56 : use units, only: getunit, freeunit
57 : use mpishorthand
58 :
59 : character(len=*), intent(in) :: nlfile ! filepath for file containing namelist input
60 :
61 : ! Local variables
62 : integer :: unitn, ierr
63 : character(len=*), parameter :: subname = 'prescribed_ghg_readnl'
64 :
65 : character(len=16) :: prescribed_ghg_specifier(N_GHG)
66 : character(len=256) :: prescribed_ghg_file
67 : character(len=256) :: prescribed_ghg_filelist
68 : character(len=256) :: prescribed_ghg_datapath
69 : character(len=32) :: prescribed_ghg_type
70 : logical :: prescribed_ghg_rmfile
71 : integer :: prescribed_ghg_cycle_yr
72 : integer :: prescribed_ghg_fixed_ymd
73 : integer :: prescribed_ghg_fixed_tod
74 :
75 : namelist /prescribed_ghg_nl/ &
76 : prescribed_ghg_specifier, &
77 : prescribed_ghg_file, &
78 : prescribed_ghg_filelist, &
79 : prescribed_ghg_datapath, &
80 : prescribed_ghg_type, &
81 : prescribed_ghg_rmfile, &
82 : prescribed_ghg_cycle_yr, &
83 : prescribed_ghg_fixed_ymd, &
84 : prescribed_ghg_fixed_tod
85 : !-----------------------------------------------------------------------------
86 :
87 : ! Initialize namelist variables from local module variables.
88 9216 : prescribed_ghg_specifier= specifier
89 1536 : prescribed_ghg_file = filename
90 1536 : prescribed_ghg_filelist = filelist
91 1536 : prescribed_ghg_datapath = datapath
92 1536 : prescribed_ghg_type = datatype
93 1536 : prescribed_ghg_rmfile = rmv_file
94 1536 : prescribed_ghg_cycle_yr = cycle_yr
95 1536 : prescribed_ghg_fixed_ymd= fixed_ymd
96 1536 : prescribed_ghg_fixed_tod= fixed_tod
97 :
98 : ! Read namelist
99 1536 : if (masterproc) then
100 2 : unitn = getunit()
101 2 : open( unitn, file=trim(nlfile), status='old' )
102 2 : call find_group_name(unitn, 'prescribed_ghg_nl', status=ierr)
103 2 : if (ierr == 0) then
104 0 : read(unitn, prescribed_ghg_nl, iostat=ierr)
105 0 : if (ierr /= 0) then
106 0 : call endrun(subname // ':: ERROR reading namelist')
107 : end if
108 : end if
109 2 : close(unitn)
110 2 : call freeunit(unitn)
111 : end if
112 :
113 : #ifdef SPMD
114 : ! Broadcast namelist variables
115 1536 : call mpibcast(prescribed_ghg_specifier,len(prescribed_ghg_specifier(1))*N_GHG, mpichar, 0, mpicom)
116 1536 : call mpibcast(prescribed_ghg_file, len(prescribed_ghg_file), mpichar, 0, mpicom)
117 1536 : call mpibcast(prescribed_ghg_filelist, len(prescribed_ghg_filelist), mpichar, 0, mpicom)
118 1536 : call mpibcast(prescribed_ghg_datapath, len(prescribed_ghg_datapath), mpichar, 0, mpicom)
119 1536 : call mpibcast(prescribed_ghg_type, len(prescribed_ghg_type), mpichar, 0, mpicom)
120 1536 : call mpibcast(prescribed_ghg_rmfile, 1, mpilog, 0, mpicom)
121 1536 : call mpibcast(prescribed_ghg_cycle_yr, 1, mpiint, 0, mpicom)
122 1536 : call mpibcast(prescribed_ghg_fixed_ymd,1, mpiint, 0, mpicom)
123 1536 : call mpibcast(prescribed_ghg_fixed_tod,1, mpiint, 0, mpicom)
124 : #endif
125 :
126 : ! Update module variables with user settings.
127 9216 : specifier = prescribed_ghg_specifier
128 1536 : filename = prescribed_ghg_file
129 1536 : filelist = prescribed_ghg_filelist
130 1536 : datapath = prescribed_ghg_datapath
131 1536 : datatype = prescribed_ghg_type
132 1536 : rmv_file = prescribed_ghg_rmfile
133 1536 : cycle_yr = prescribed_ghg_cycle_yr
134 1536 : fixed_ymd = prescribed_ghg_fixed_ymd
135 1536 : fixed_tod = prescribed_ghg_fixed_tod
136 :
137 : ! Turn on prescribed volcanics if user has specified an input dataset.
138 1536 : if (len_trim(filename) > 0 .and. filename.ne.'NONE') has_prescribed_ghg = .true.
139 :
140 1536 : end subroutine prescribed_ghg_readnl
141 :
142 : !-------------------------------------------------------------------
143 : !-------------------------------------------------------------------
144 1536 : subroutine prescribed_ghg_register()
145 : use ppgrid, only: pver, pcols
146 : use physics_buffer, only : pbuf_add_field, dtype_r8
147 :
148 : integer :: i,idx
149 :
150 1536 : if (has_prescribed_ghg) then
151 0 : do i = 1,N_GHG
152 0 : call pbuf_add_field(ghg_names(i),'physpkg',dtype_r8,(/pcols,pver/),idx)
153 : enddo
154 : endif
155 :
156 1536 : endsubroutine prescribed_ghg_register
157 : !-------------------------------------------------------------------
158 : !-------------------------------------------------------------------
159 1536 : subroutine prescribed_ghg_init()
160 :
161 1536 : use tracer_data, only : trcdata_init
162 : use cam_history, only : addfld
163 :
164 : implicit none
165 :
166 : integer :: ndx, istat, i
167 :
168 1536 : if ( has_prescribed_ghg ) then
169 0 : if ( masterproc ) then
170 0 : write(iulog,*) 'ghg is prescribed in :'//trim(filename)
171 : endif
172 : else
173 : return
174 : endif
175 :
176 0 : allocate(file%in_pbuf(size(specifier)))
177 0 : file%in_pbuf(:) = .true.
178 : call trcdata_init( specifier, filename, filelist, datapath, fields, file, &
179 0 : rmv_file, cycle_yr, fixed_ymd, fixed_tod, datatype)
180 :
181 0 : number_flds = 0
182 0 : if (associated(fields)) number_flds = size( fields )
183 :
184 0 : if( number_flds < 1 ) then
185 0 : if ( masterproc ) then
186 0 : write(iulog,*) 'There are no prescribed ghg tracers'
187 0 : write(iulog,*) ' '
188 : endif
189 0 : return
190 : end if
191 :
192 0 : do i = 1,number_flds
193 0 : ndx = get_ndx( fields(i)%fldnam )
194 0 : index_map(i) = ndx
195 :
196 0 : if (ndx < 1) then
197 0 : call endrun('prescribed_ghg_init: '//trim(fields(i)%fldnam)//' is not one of the named ghg fields in pbuf2d')
198 : endif
199 0 : call addfld( fields(i)%fldnam, (/ 'lev' /), 'I','kg/kg', 'prescribed ghg' )
200 : enddo
201 :
202 1536 : end subroutine prescribed_ghg_init
203 :
204 : !-------------------------------------------------------------------
205 : !-------------------------------------------------------------------
206 741888 : subroutine prescribed_ghg_adv( state, pbuf2d)
207 :
208 1536 : use tracer_data, only : advance_trcdata
209 : use physics_types,only : physics_state
210 : use ppgrid, only : begchunk, endchunk
211 : use ppgrid, only : pcols, pver
212 : use string_utils, only : to_lower, GLC
213 : use cam_history, only : outfld
214 : use physconst, only : mwdry ! molecular weight dry air ~ kg/kmole
215 : use physconst, only : boltz ! J/K/molecule
216 :
217 : use physics_buffer, only : physics_buffer_desc, pbuf_get_field, pbuf_set_field, pbuf_get_chunk
218 :
219 : implicit none
220 :
221 : type(physics_state), intent(in) :: state(begchunk:endchunk)
222 :
223 : type(physics_buffer_desc), pointer :: pbuf2d(:,:)
224 :
225 370944 : type(physics_buffer_desc), pointer :: pbuf_chnk(:)
226 : integer :: ind,c,ncol,i
227 : real(r8) :: to_mmr(pcols,pver)
228 : real(r8) :: outdata(pcols,pver)
229 370944 : real(r8),pointer :: tmpptr(:,:)
230 :
231 : character(len=32) :: units_str
232 :
233 370944 : if( .not. has_prescribed_ghg ) return
234 :
235 0 : call advance_trcdata( fields, file, state, pbuf2d )
236 :
237 : ! set the correct units and invoke history outfld
238 0 : do i = 1,number_flds
239 0 : ind = index_map(i)
240 :
241 0 : units_str = trim(to_lower(trim(fields(i)%units(:GLC(fields(i)%units)))))
242 :
243 : !$OMP PARALLEL DO PRIVATE (C, NCOL, OUTDATA, TO_MMR, tmpptr, pbuf_chnk)
244 0 : do c = begchunk,endchunk
245 0 : ncol = state(c)%ncol
246 :
247 : select case ( units_str )
248 : case ("molec/cm3","/cm3","molecules/cm3","cm^-3","cm**-3")
249 0 : to_mmr(:ncol,:) = (molmass(ind)*1.e6_r8*boltz*state(c)%t(:ncol,:))/(mwdry*state(c)%pmiddry(:ncol,:))
250 : case ('kg/kg','mmr')
251 0 : to_mmr(:ncol,:) = 1._r8
252 : case ('mol/mol','mole/mole','vmr','fraction')
253 0 : to_mmr(:ncol,:) = molmass(ind)/mwdry
254 : case default
255 0 : print*, 'prescribed_ghg_adv: units = ',trim(fields(i)%units) ,' are not recognized'
256 0 : call endrun('prescribed_ghg_adv: units are not recognized')
257 : end select
258 :
259 0 : pbuf_chnk => pbuf_get_chunk(pbuf2d, c)
260 0 : call pbuf_get_field(pbuf_chnk, fields(i)%pbuf_ndx, tmpptr )
261 :
262 0 : tmpptr(:ncol,:) = tmpptr(:ncol,:)*to_mmr(:ncol,:)
263 :
264 0 : outdata(:ncol,:) = tmpptr(:ncol,:)
265 0 : call outfld( fields(1)%fldnam, outdata(:ncol,:), ncol, state(c)%lchnk )
266 :
267 : enddo
268 : enddo
269 :
270 370944 : end subroutine prescribed_ghg_adv
271 :
272 : !-------------------------------------------------------------------
273 :
274 : !-------------------------------------------------------------------
275 1536 : subroutine init_prescribed_ghg_restart( piofile )
276 370944 : use pio, only : file_desc_t
277 : use tracer_data, only : init_trc_restart
278 : implicit none
279 : type(file_desc_t),intent(inout) :: pioFile ! pio File pointer
280 :
281 1536 : call init_trc_restart( 'prescribed_ghg', piofile, file )
282 :
283 1536 : end subroutine init_prescribed_ghg_restart
284 : !-------------------------------------------------------------------
285 1536 : subroutine write_prescribed_ghg_restart( piofile )
286 1536 : use tracer_data, only : write_trc_restart
287 : use pio, only : file_desc_t
288 : implicit none
289 :
290 : type(file_desc_t) :: piofile
291 :
292 1536 : call write_trc_restart( piofile, file )
293 :
294 1536 : end subroutine write_prescribed_ghg_restart
295 :
296 : !-------------------------------------------------------------------
297 : !-------------------------------------------------------------------
298 768 : subroutine read_prescribed_ghg_restart( pioFile )
299 1536 : use tracer_data, only : read_trc_restart
300 : use pio, only : file_desc_t
301 : implicit none
302 :
303 : type(file_desc_t) :: piofile
304 :
305 768 : call read_trc_restart( 'prescribed_ghg', piofile, file )
306 :
307 768 : end subroutine read_prescribed_ghg_restart
308 : !-------------------------------------------------------------------
309 0 : integer function get_ndx( name )
310 :
311 : implicit none
312 : character(len=*), intent(in) :: name
313 :
314 : integer :: i
315 :
316 0 : get_ndx = 0
317 0 : do i = 1,N_GHG
318 0 : if ( trim(name) == trim(ghg_names(i)) ) then
319 0 : get_ndx = i
320 0 : return
321 : endif
322 : enddo
323 :
324 768 : end function get_ndx
325 :
326 : end module prescribed_ghg
|
