PartMC  2.2.1
gas_data.F90
Go to the documentation of this file.
00001 ! Copyright (C) 2005-2011 Nicole Riemer and Matthew West
00002 ! Licensed under the GNU General Public License version 2 or (at your
00003 ! option) any later version. See the file COPYING for details.
00004 
00005 !> \file
00006 !> The pmc_gas_data module.
00007 
00008 !> The gas_data_t structure and associated subroutines.
00009 module pmc_gas_data
00010 
00011   use pmc_spec_file
00012   use pmc_mpi
00013   use pmc_util
00014   use pmc_netcdf
00015 #ifdef PMC_USE_MPI
00016   use mpi
00017 #endif
00018 
00019   !> Maximum length of the name of a gas.
00020   integer, parameter :: GAS_NAME_LEN = 100
00021 
00022   !> Constant gas data.
00023   !!
00024   !! Each gas species is identified by an integer \c i between 1 and
00025   !! \c n_spec. Species \c i has name \c gas_data%%name(i). The
00026   !! variable gas data describing the current mixing ratios is stored
00027   !! in the gas_state_t structure, so the mixing ratio of species \c i
00028   !! is gas_state%%mix_rat(i).
00029   type gas_data_t
00030      !> Number of species.
00031      integer :: n_spec
00032      !> Species name [length \c n_spec].
00033      character(len=GAS_NAME_LEN), pointer :: name(:)
00034      !> Index of the corresponding MOSAIC species [length \c
00035      !> n_spec]. \c to_mosaic(i) is the mosaic index of species \c i,
00036      !> or 0 if there is no match.
00037      integer, pointer :: mosaic_index(:)
00038   end type gas_data_t
00039 
00040 contains
00041 
00042 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00043 
00044   !> Allocate storage for gas species.
00045   subroutine gas_data_allocate(gas_data)
00046 
00047     !> Gas data.
00048     type(gas_data_t), intent(out) :: gas_data
00049 
00050     gas_data%n_spec = 0
00051     allocate(gas_data%name(0))
00052     allocate(gas_data%mosaic_index(0))
00053 
00054   end subroutine gas_data_allocate
00055 
00056 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00057 
00058   !> Allocate storage for gas species with the given size.
00059   subroutine gas_data_allocate_size(gas_data, n_spec)
00060 
00061     !> Gas data.
00062     type(gas_data_t), intent(out) :: gas_data
00063     !> Number of species.
00064     integer, intent(in) :: n_spec
00065 
00066     gas_data%n_spec = n_spec
00067     allocate(gas_data%name(n_spec))
00068     allocate(gas_data%mosaic_index(n_spec))
00069 
00070   end subroutine gas_data_allocate_size
00071 
00072 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00073 
00074   !> Free all storage.
00075   subroutine gas_data_deallocate(gas_data)
00076 
00077     !> Gas data.
00078     type(gas_data_t), intent(inout) :: gas_data
00079 
00080     deallocate(gas_data%name)
00081     deallocate(gas_data%mosaic_index)
00082 
00083   end subroutine gas_data_deallocate
00084 
00085 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00086 
00087   !> Returns the number of the species in gas with the given name, or
00088   !> returns 0 if there is no such species.
00089   integer function gas_data_spec_by_name(gas_data, name)
00090 
00091     !> Gas data.
00092     type(gas_data_t), intent(in) :: gas_data
00093     !> Name of species to find.
00094     character(len=*), intent(in) :: name
00095 
00096     integer i
00097     logical found
00098 
00099     found = .false.
00100     do i = 1,gas_data%n_spec
00101        if (name == gas_data%name(i)) then
00102           found = .true.
00103           exit
00104        end if
00105     end do
00106     if (found) then
00107        gas_data_spec_by_name = i
00108     else
00109        gas_data_spec_by_name = 0
00110     end if
00111 
00112   end function gas_data_spec_by_name
00113 
00114 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00115 
00116   !> Fills in gas_data%mosaic_index.
00117   subroutine gas_data_set_mosaic_map(gas_data)
00118 
00119     !> Gas data.
00120     type(gas_data_t), intent(inout) :: gas_data
00121 
00122     integer, parameter :: n_mosaic_species = 77
00123     character(len=GAS_NAME_LEN), parameter, dimension(n_mosaic_species) 
00124          :: mosaic_species = [ 
00125          "H2SO4        ", "HNO3         ", "HCl          ", 
00126          "NH3          ", "NO           ", "NO2          ", 
00127          "NO3          ", "N2O5         ", "HONO         ", 
00128          "HNO4         ", "O3           ", "O1D          ", 
00129          "O3P          ", "OH           ", "HO2          ", 
00130          "H2O2         ", "CO           ", "SO2          ", 
00131          "CH4          ", "C2H6         ", "CH3O2        ", 
00132          "ETHP         ", "HCHO         ", "CH3OH        ", 
00133          "ANOL         ", "CH3OOH       ", "ETHOOH       ", 
00134          "ALD2         ", "HCOOH        ", "RCOOH        ", 
00135          "C2O3         ", "PAN          ", "ARO1         ", 
00136          "ARO2         ", "ALK1         ", "OLE1         ", 
00137          "API1         ", "API2         ", "LIM1         ", 
00138          "LIM2         ", "PAR          ", "AONE         ", 
00139          "MGLY         ", "ETH          ", "OLET         ", 
00140          "OLEI         ", "TOL          ", "XYL          ", 
00141          "CRES         ", "TO2          ", "CRO          ", 
00142          "OPEN         ", "ONIT         ", "ROOH         ", 
00143          "RO2          ", "ANO2         ", "NAP          ", 
00144          "XO2          ", "XPAR         ", "ISOP         ", 
00145          "ISOPRD       ", "ISOPP        ", "ISOPN        ", 
00146          "ISOPO2       ", "API          ", "LIM          ", 
00147          "DMS          ", "MSA          ", "DMSO         ", 
00148          "DMSO2        ", "CH3SO2H      ", "CH3SCH2OO    ", 
00149          "CH3SO2       ", "CH3SO3       ", "CH3SO2OO     ", 
00150          "CH3SO2CH2OO  ", "SULFHOX      "]
00151 
00152     integer spec, mosaic_spec, i
00153 
00154     gas_data%mosaic_index = 0
00155     do spec = 1,gas_data%n_spec
00156        mosaic_spec = 0
00157        do i = 1,n_mosaic_species
00158           if (gas_data%name(spec) == mosaic_species(i)) then
00159              mosaic_spec = i
00160           end if
00161        end do
00162        if (mosaic_spec > 0) then
00163           gas_data%mosaic_index(spec) = mosaic_spec
00164        end if
00165     end do
00166 
00167   end subroutine gas_data_set_mosaic_map
00168 
00169 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00170 
00171   !> Read gas data from a .spec file.
00172   subroutine spec_file_read_gas_data(file, gas_data)
00173 
00174     !> Spec file to read data from.
00175     type(spec_file_t), intent(inout) :: file
00176     !> Gas data.
00177     type(gas_data_t), intent(inout) :: gas_data
00178 
00179     integer :: n_species, species, i
00180     character(len=SPEC_LINE_MAX_VAR_LEN), pointer :: species_name(:)
00181     real(kind=dp), pointer :: species_data(:,:)
00182 
00183     !> \page input_format_gas_data Input File Format: Gas Material Data
00184     !!
00185     !! A gas material data file must consist of one line per gas
00186     !! species, with each line having the species name. This specifies
00187     !! which species are to be recognized as gases. For example, a \c
00188     !! gas_data file could contain:
00189     !! <pre>
00190     !! H2SO4
00191     !! HNO3
00192     !! HCl
00193     !! NH3
00194     !! </pre>
00195     !!
00196     !! See also:
00197     !!   - \ref spec_file_format --- the input file text format
00198     !!   - \ref output_format_gas_data --- the corresponding output format
00199 
00200     ! read the gas data from the specified file
00201     allocate(species_name(0))
00202     allocate(species_data(0,0))
00203     call spec_file_read_real_named_array(file, 0, species_name, &
00204          species_data)
00205 
00206     ! check the data size
00207     if (size(species_data, 2) /= 0) then
00208        call die_msg(614290516, 'each line in ' // trim(file%name) &
00209             // ' must only contain the species name')
00210     end if
00211 
00212     ! allocate and copy over the data
00213     n_species = size(species_data, 1)
00214     call gas_data_deallocate(gas_data)
00215     call gas_data_allocate_size(gas_data, n_species)
00216     do i = 1,n_species
00217        gas_data%name(i) = species_name(i)(1:GAS_NAME_LEN)
00218     end do
00219     deallocate(species_name)
00220     deallocate(species_data)
00221     
00222     call gas_data_set_mosaic_map(gas_data)
00223 
00224   end subroutine spec_file_read_gas_data
00225 
00226 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00227 
00228   !> Determines the number of bytes required to pack the given value.
00229   integer function pmc_mpi_pack_size_gas_data(val)
00230 
00231     !> Value to pack.
00232     type(gas_data_t), intent(in) :: val
00233 
00234     pmc_mpi_pack_size_gas_data = &
00235          pmc_mpi_pack_size_integer(val%n_spec) &
00236          + pmc_mpi_pack_size_string_array(val%name) &
00237          + pmc_mpi_pack_size_integer_array(val%mosaic_index)
00238 
00239   end function pmc_mpi_pack_size_gas_data
00240 
00241 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00242 
00243   !> Packs the given value into the buffer, advancing position.
00244   subroutine pmc_mpi_pack_gas_data(buffer, position, val)
00245 
00246     !> Memory buffer.
00247     character, intent(inout) :: buffer(:)
00248     !> Current buffer position.
00249     integer, intent(inout) :: position
00250     !> Value to pack.
00251     type(gas_data_t), intent(in) :: val
00252 
00253 #ifdef PMC_USE_MPI
00254     integer :: prev_position
00255 
00256     prev_position = position
00257     call pmc_mpi_pack_integer(buffer, position, val%n_spec)
00258     call pmc_mpi_pack_string_array(buffer, position, val%name)
00259     call pmc_mpi_pack_integer_array(buffer, position, val%mosaic_index)
00260     call assert(449872094, &
00261          position - prev_position <= pmc_mpi_pack_size_gas_data(val))
00262 #endif
00263 
00264   end subroutine pmc_mpi_pack_gas_data
00265 
00266 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00267 
00268   !> Unpacks the given value from the buffer, advancing position.
00269   subroutine pmc_mpi_unpack_gas_data(buffer, position, val)
00270 
00271     !> Memory buffer.
00272     character, intent(inout) :: buffer(:)
00273     !> Current buffer position.
00274     integer, intent(inout) :: position
00275     !> Value to pack.
00276     type(gas_data_t), intent(inout) :: val
00277 
00278 #ifdef PMC_USE_MPI
00279     integer :: prev_position
00280 
00281     prev_position = position
00282     call pmc_mpi_unpack_integer(buffer, position, val%n_spec)
00283     call pmc_mpi_unpack_string_array(buffer, position, val%name)
00284     call pmc_mpi_unpack_integer_array(buffer, position, val%mosaic_index)
00285     call assert(137879163, &
00286          position - prev_position <= pmc_mpi_pack_size_gas_data(val))
00287 #endif
00288 
00289   end subroutine pmc_mpi_unpack_gas_data
00290 
00291 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00292 
00293   !> Write the gas species dimension to the given NetCDF file if it
00294   !> is not already present and in any case return the associated
00295   !> dimid.
00296   subroutine gas_data_netcdf_dim_gas_species(gas_data, ncid, &
00297        dimid_gas_species)
00298 
00299     !> Gas_data structure.
00300     type(gas_data_t), intent(in) :: gas_data
00301     !> NetCDF file ID, in data mode.
00302     integer, intent(in) :: ncid
00303     !> Dimid of the species dimension.
00304     integer, intent(out) :: dimid_gas_species
00305 
00306     integer :: status, i_spec
00307     integer :: varid_gas_species
00308     integer :: gas_species_centers(gas_data%n_spec)
00309     character(len=((GAS_NAME_LEN + 2) * gas_data%n_spec)) :: gas_species_names
00310 
00311     ! try to get the dimension ID
00312     status = nf90_inq_dimid(ncid, "gas_species", dimid_gas_species)
00313     if (status == NF90_NOERR) return
00314     if (status /= NF90_EBADDIM) call pmc_nc_check(status)
00315 
00316     ! dimension not defined, so define now define it
00317     call pmc_nc_check(nf90_redef(ncid))
00318 
00319     call pmc_nc_check(nf90_def_dim(ncid, "gas_species", &
00320          gas_data%n_spec, dimid_gas_species))
00321     gas_species_names = ""
00322     do i_spec = 1,gas_data%n_spec
00323        gas_species_names((len_trim(gas_species_names) + 1):) &
00324             = trim(gas_data%name(i_spec))
00325        if (i_spec < gas_data%n_spec) then
00326           gas_species_names((len_trim(gas_species_names) + 1):) = ","
00327        end if
00328     end do
00329     call pmc_nc_check(nf90_def_var(ncid, "gas_species", NF90_INT, &
00330          dimid_gas_species, varid_gas_species))
00331     call pmc_nc_check(nf90_put_att(ncid, varid_gas_species, "names", &
00332          gas_species_names))
00333     call pmc_nc_check(nf90_put_att(ncid, varid_gas_species, "description", &
00334          "dummy dimension variable (no useful value) - read species names " &
00335          // "as comma-separated values from the 'names' attribute"))
00336 
00337     call pmc_nc_check(nf90_enddef(ncid))
00338 
00339     do i_spec = 1,gas_data%n_spec
00340        gas_species_centers(i_spec) = i_spec
00341     end do
00342     call pmc_nc_check(nf90_put_var(ncid, varid_gas_species, &
00343          gas_species_centers))
00344 
00345   end subroutine gas_data_netcdf_dim_gas_species
00346 
00347 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00348 
00349   !> Write full state.
00350   subroutine gas_data_output_netcdf(gas_data, ncid)
00351     
00352     !> Gas_data to write.
00353     type(gas_data_t), intent(in) :: gas_data
00354     !> NetCDF file ID, in data mode.
00355     integer, intent(in) :: ncid
00356 
00357     integer :: dimid_gas_species
00358 
00359     !> \page output_format_gas_data Output File Format: Gas Material Data
00360     !!
00361     !! The gas material data NetCDF dimensions are:
00362     !!   - \b gas_species: number of gas species
00363     !!
00364     !! The gas material data NetCDF variables are:
00365     !!   - \b gas_species (dim \c gas_species): dummy dimension variable
00366     !!     (no useful value) - read species names as comma-separated values
00367     !!     from the 'names' attribute
00368     !!   - \b gas_mosaic_index (dim \c gas_species): MOSAIC indices of
00369     !!     gas species
00370     !!
00371     !! See also:
00372     !!   - \ref input_format_gas_data --- the corresponding input format
00373 
00374     call gas_data_netcdf_dim_gas_species(gas_data, ncid, &
00375          dimid_gas_species)
00376 
00377     call pmc_nc_write_integer_1d(ncid, gas_data%mosaic_index, &
00378          "gas_mosaic_index", (/ dimid_gas_species /), &
00379          long_name="MOSAIC indices of gas species")
00380 
00381   end subroutine gas_data_output_netcdf
00382 
00383 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00384 
00385   !> Read full state.
00386   subroutine gas_data_input_netcdf(gas_data, ncid)
00387     
00388     !> Gas_data to read.
00389     type(gas_data_t), intent(inout) :: gas_data
00390     !> NetCDF file ID, in data mode.
00391     integer, intent(in) :: ncid
00392 
00393     character(len=1000) :: name
00394     integer :: dimid_gas_species, n_spec, varid_gas_species, i_spec, i
00395     character(len=((GAS_NAME_LEN + 2) * 1000)) :: gas_species_names
00396 
00397     call pmc_nc_check(nf90_inq_dimid(ncid, "gas_species", dimid_gas_species))
00398     call pmc_nc_check(nf90_Inquire_Dimension(ncid, dimid_gas_species, name, &
00399          n_spec))
00400     call gas_data_deallocate(gas_data)
00401     call gas_data_allocate_size(gas_data, n_spec)
00402     call assert(719237193, n_spec < 1000)
00403 
00404     call pmc_nc_read_integer_1d(ncid, gas_data%mosaic_index, &
00405          "gas_mosaic_index")
00406 
00407     call pmc_nc_check(nf90_inq_varid(ncid, "gas_species", varid_gas_species))
00408     call pmc_nc_check(nf90_get_att(ncid, varid_gas_species, "names", &
00409          gas_species_names))
00410     ! gas_species_names are comma-separated, so unpack them
00411     do i_spec = 1,gas_data%n_spec
00412        i = 1
00413        do while ((gas_species_names(i:i) /= " ") &
00414             .and. (gas_species_names(i:i) /= ","))
00415           i = i + 1
00416        end do
00417        call assert(173021381, i > 1)
00418        gas_data%name(i_spec) = gas_species_names(1:(i-1))
00419        gas_species_names = gas_species_names((i+1):)
00420     end do
00421     call assert(469721220, gas_species_names == "")
00422 
00423   end subroutine gas_data_input_netcdf
00424 
00425 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00426 
00427 end module pmc_gas_data