!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !! !! !! GNU General Public License !! !! !! !! This file is part of the Flexible Modeling System (FMS). !! !! !! !! FMS is free software; you can redistribute it and/or modify !! !! it and are expected to follow the terms of the GNU General Public !! !! License as published by the Free Software Foundation. !! !! !! !! FMS is distributed in the hope that it will be useful, !! !! but WITHOUT ANY WARRANTY; without even the implied warranty of !! !! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the !! !! GNU General Public License for more details. !! !! !! !! You should have received a copy of the GNU General Public License !! !! along with FMS; if not, write to: !! !! Free Software Foundation, Inc. !! !! 59 Temple Place, Suite 330 !! !! Boston, MA 02111-1307 USA !! !! or see: !! !! http://www.gnu.org/licenses/gpl.txt !! !! !! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! module atmos_ch3i_mod use fms_mod, only : file_exist, & write_version_number, & error_mesg, & FATAL, & NOTE, & mpp_pe, & mpp_root_pe, & lowercase, & open_namelist_file, & check_nml_error, & close_file, & stdlog use field_manager_mod, only : MODEL_ATMOS, & parse use tracer_manager_mod, only : get_tracer_index, & get_tracer_names, & query_method use time_manager_mod, only : time_type use diag_manager_mod, only : send_data, & register_diag_field use interpolator_mod, only : interpolate_type, & interpolator_init, & obtain_interpolator_time_slices, & unset_interpolator_time_flag, & interpolator_end, & interpolator, & query_interpolator, & CONSTANT, & INTERP_WEIGHTED_P use constants_mod, only : WTMAIR, & AVOGNO, & SECONDS_PER_DAY implicit none private public :: atmos_ch3i_init, atmos_ch3i_time_vary, atmos_ch3i, & atmos_ch3i_endts, atmos_ch3i_end !----------------------------------------------------------------------- ! ... namelist !----------------------------------------------------------------------- character(len=128) :: conc_filename = '' namelist /atmos_ch3i_nml/ & conc_filename logical :: has_emissions character(len=128) :: emis_filename type(interpolate_type), save :: ch3i_emissions, input_conc character(len=128), allocatable :: field_names(:) integer :: id_emissions, id_loss, id_j_ch3i, ind_ch3i real, parameter :: g_to_kg = 1.e-3, & !conversion factor (kg/g) m2_to_cm2 = 1.e4 !conversion factor (cm2/m2) real, parameter :: emis_cons = WTMAIR * g_to_kg * m2_to_cm2 / AVOGNO real, parameter :: boltz = 1.38044e-16 ! Boltzmann's Constant (erg/K) character(len=7), parameter :: module_name = 'tracers' !---- version number ----- character(len=128) :: version = '$Id: atmos_ch3i.F90,v 15.0.4.1.2.1.2.1 2009/07/30 10:52:46 rsh Exp $' character(len=128) :: tagname = '$Name: mom4p1_pubrel_dec2009_nnz $' logical :: module_is_initialized = .FALSE. contains subroutine atmos_ch3i_init( lonb_mod, latb_mod, axes, Time, mask ) !----------------------------------------------------------------------- ! ! mask = optional mask (0. or 1.) that designates which grid points ! are above (=1.) or below (=0.) the ground dimensioned as ! (nlon,nlat,nlev). ! !----------------------------------------------------------------------- real, intent(in), dimension(:,:) :: lonb_mod real, intent(in), dimension(:,:) :: latb_mod type(time_type), intent(in) :: Time integer , intent(in) :: axes(4) real, intent(in), dimension(:,:,:), optional :: mask !----------------------------------------------------------------------- integer :: unit, nfields, flag_file integer :: ierr, io, logunit character(len=128) :: tracer_name, tracer_units, name, control if (module_is_initialized) then return end if ind_ch3i = get_tracer_index(MODEL_ATMOS, 'ch3i') !----------------------------------------------------------------------- ! ... write version number !----------------------------------------------------------------------- call write_version_number(version, tagname) !----------------------------------------------------------------------- ! ... read namelist !----------------------------------------------------------------------- if ( file_exist('input.nml')) then unit = open_namelist_file () ierr=1; do while (ierr /= 0) read (unit, nml=atmos_ch3i_nml, iostat=io, end=10) ierr = check_nml_error(io, 'atmos_ch3i_nml') enddo 10 call close_file (unit) endif logunit=stdlog() if(mpp_pe() == mpp_root_pe()) then write(logunit, nml=atmos_ch3i_nml) endif !----------------------------------------------------------------------- ! ... set up emissions !----------------------------------------------------------------------- has_emissions = .false. if( query_method('emissions',MODEL_ATMOS,ind_ch3i,name,control) ) then if( trim(name) == 'file' ) then flag_file = parse(control, 'file',emis_filename) if (flag_file > 0) then has_emissions = .true. call interpolator_init( ch3i_emissions, trim(emis_filename), & lonb_mod, latb_mod, & data_out_of_bounds=(/CONSTANT/) ) call query_interpolator(ch3i_emissions,nfields=nfields) allocate( field_names(nfields) ) call query_interpolator(ch3i_emissions,field_names=field_names) id_emissions = register_diag_field( module_name, 'ch3i_emis', axes(1:2), & Time, 'ch3i_emis', 'molec/cm2/s' ) else call error_mesg('atmos_ch3i_init','CH3I emission file not specified in field table',FATAL) end if else call error_mesg('atmos_ch3i_init','CH3I emission file not specified in field table',FATAL) end if else call error_mesg('atmos_ch3i_init','No emissions specified for CH3I in field table',NOTE) end if !----------------------------------------------------------------------- ! ... set up tracer concentrations !----------------------------------------------------------------------- if (conc_filename /= '') then call interpolator_init( input_conc,trim(conc_filename), lonb_mod, latb_mod,& data_out_of_bounds=(/CONSTANT/), & vert_interp=(/INTERP_WEIGHTED_P/) ) end if call get_tracer_names( MODEL_ATMOS,ind_ch3i,tracer_name,units=tracer_units) id_loss =register_diag_field( module_name, 'ch3i_loss', axes(1:3), Time, 'ch3i_loss', TRIM(tracer_units)//'/s' ) id_j_ch3i =register_diag_field( module_name, 'j_ch3i', axes(1:3), Time, 'j_ch3i', '/s' ) module_is_initialized = .true. end subroutine atmos_ch3i_init !##################################################################### subroutine atmos_ch3i_time_vary (Time) type(time_type), intent(in) :: Time if (has_emissions) then call obtain_interpolator_time_slices (ch3i_emissions, Time) endif if (conc_filename /= '') then call obtain_interpolator_time_slices (input_conc, Time) endif end subroutine atmos_ch3i_time_vary !###################################################################### subroutine atmos_ch3i_endts if (has_emissions) then call unset_interpolator_time_flag (ch3i_emissions) endif if (conc_filename /= '') then call unset_interpolator_time_flag (input_conc) endif end subroutine atmos_ch3i_endts !----------------------------------------------------------------------- subroutine atmos_ch3i( lon, lat, land, pwt, ch3i, ch3i_dt, & Time, phalf, pfull, t, is, js, dt, & z_half, z_full, q, tsurf, albedo, coszen, & Time_next, kbot) real, intent(in), dimension(:,:) :: lon, lat real, intent(in), dimension(:,:) :: land real, intent(in), dimension(:,:,:) :: pwt real, intent(in), dimension(:,:,:) :: ch3i real, intent(out), dimension(:,:,:) :: ch3i_dt type(time_type), intent(in) :: Time, Time_next integer, intent(in) :: is,js real, intent(in), dimension(:,:,:) :: phalf,pfull,t real, intent(in) :: dt !to be passed into chemdr real, intent(in), dimension(:,:,:) :: z_half !height in meters at half levels real, intent(in), dimension(:,:,:) :: z_full !height in meters at full levels real, intent(in), dimension(:,:,:) :: q !specific humidity at current time step in kg/kg real, intent(in), dimension(:,:) :: tsurf !surface temperature real, intent(in), dimension(:,:) :: albedo real, intent(in), dimension(:,:) :: coszen integer, intent(in), dimension(:,:), optional :: kbot !----------------------------------------------------------------------- logical :: used integer :: i, j, k, id, jd, kd, kb real, dimension(size(ch3i,1),size(ch3i,2),size(ch3i,3)) :: & conc_o3, conc_oh, j_ch3i, k_ch3i_oh, & emis_source, ch3i_loss, & air_dens real, dimension(size(ch3i,1),size(ch3i,2)) :: emis,temp_data ! ! Tracer_driver_init needs to be called before tracer_driver. ! if (.not. module_is_initialized) & call error_mesg ('atmos_ch3i','atmos_ch3i_init must be called first1', FATAL) id=size(ch3i,1); jd=size(ch3i,2); kd=size(ch3i,3) !----------------------------------------------------------------------- ! ... Get emissions !----------------------------------------------------------------------- emis_source(:,:,:) = 0. emis(:,:) = 0. if (has_emissions) then do k = 1,size(field_names) call interpolator( ch3i_emissions, Time, temp_data, field_names(k), is, js ) emis(:,:) = emis(:,:) + temp_data(:,:) end do if (present(kbot)) then do j=1,jd do i=1,id kb=kbot(i,j) emis_source(i,j,kb) = emis(i+is,j+js)/pwt(i,j,kb) * emis_cons enddo enddo else emis_source(:,:,kd) = emis(:,:)/pwt(:,:,kd) * emis_cons end if used = send_data(id_emissions,emis,Time,is_in=is,js_in=js) end if !----------------------------------------------------------------------- ! ... Get tracer concentrations !----------------------------------------------------------------------- call interpolator(input_conc, Time, phalf, conc_o3, 'ox') call interpolator(input_conc, Time, phalf, conc_oh, 'oh') air_dens(:,:,:) = 10. * pfull(:,:,:) / (boltz*t(:,:,:)) ! molec/cm3 conc_oh(:,:,:) = conc_oh(:,:,:) * air_dens(:,:,:) ! convert VMR to molec/cm3 !----------------------------------------------------------------------- ! ... reaction rates !----------------------------------------------------------------------- k_ch3i_oh(:,:,:) = 2.9e-12 * exp( -1100./t(:,:,:) ) ! s^-1 do k = 1,kd where( coszen(:,:) > 0. ) j_ch3i(:,:,k) = 2./(4.*SECONDS_PER_DAY) ! s^-1 elsewhere j_ch3i(:,:,k) = 0. endwhere end do ch3i_loss(:,:,:) = ( j_ch3i(:,:,:) + k_ch3i_oh(:,:,:)*conc_oh(:,:,:) ) & * ch3i(:,:,:) ! VMR/s used = send_data( id_loss, ch3i_loss, Time, is_in=is, js_in=js ) used = send_data( id_j_ch3i, j_ch3i, Time, is_in=is, js_in=js ) ch3i_dt(:,:,:) = emis_source(:,:,:) - ch3i_loss(:,:,:) end subroutine atmos_ch3i !----------------------------------------------------------------------- subroutine atmos_ch3i_end deallocate( field_names ) module_is_initialized = .FALSE. end subroutine atmos_ch3i_end end module atmos_ch3i_mod