!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!! !!
!! 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_co2_mod
!
! Richard D. Slater
!
!
! none yet
!
!
!
!
!
!
use fms_mod, only : file_exist, write_version_number, &
mpp_pe, mpp_root_pe, &
close_file, stdlog, stdout, &
check_nml_error, error_mesg, &
open_namelist_file, FATAL, NOTE, WARNING
use tracer_manager_mod, only : get_tracer_index, tracer_manager_init
use field_manager_mod, only : MODEL_ATMOS
use diag_manager_mod, only : register_diag_field, send_data
use time_manager_mod, only : time_type
use constants_mod, only : WTMCO2, WTMAIR
use data_override_mod, only : data_override
use mpp_mod, only : mpp_pe, mpp_root_pe
implicit none
private
!-----------------------------------------------------------------------
!----- interfaces -------
public atmos_co2_sourcesink
public atmos_co2_rad
public atmos_co2_gather_data
public atmos_co2_flux_init
public atmos_co2_init
public atmos_co2_end
!-----------------------------------------------------------------------
!----------- namelist -------------------
!-----------------------------------------------------------------------
character(len=48), parameter :: mod_name = 'atmos_co2_mod'
integer, save :: ind_co2_flux = 0
integer, save :: ind_co2 = 0
integer, save :: ind_sphum = 0
integer :: id_co2restore, id_pwt
real :: restore_co2_dvmr = -1
real :: radiation_co2_dvmr = -1
!---------------------------------------------------------------------
!-------- namelist ---------
real :: restore_tscale = -1
integer :: restore_klimit = -1
logical :: do_co2_restore = .false.
logical :: co2_radiation_override = .false.
namelist /atmos_co2_nml/ &
do_co2_restore, restore_tscale, restore_klimit, &
co2_radiation_override
!-----------------------------------------------------------------------
!
! When initializing additional tracers, the user needs to make the
! following changes.
!
! Add an integer variable below for each additional tracer. This should
! be initialized to zero.
!
! Add id_tracername for each additional tracer. These are used in
! initializing and outputting the tracer fields.
!
!-----------------------------------------------------------------------
!PUBLIC VARIABLES
public :: co2_radiation_override
logical :: module_is_initialized = .FALSE.
logical :: used
!---- version number -----
character(len=128) :: version = '$$'
character(len=128) :: tagname = '$$'
!-----------------------------------------------------------------------
contains
!#######################################################################
!
!
! A subroutine to calculate the internal sources and sinks of carbon dioxide.
!
! do_co2_restore = logical to turn co2_restore on/off: default = .false.
! restore_co2_dvmr = partial pressure of co2 to which to restore (mol/mol)
! restore_klimit = atmospheric level to which to restore starting from top
! restore_tscale = timescale in seconds with which to restore
!
!
!
! A routine to calculate the sources and sinks of carbon dixoide.
!
!
!call atmos_co2_sourcesink (Time, dt, pwt, co2, sphum, co2_restore)
!
!
!
! Model time.
!
!
! Model timestep.
!
!
! The pressure weighting array. = dP/grav (kg/m2)
!
!
! The array of the carbon dioxide mixing ratio (kg co2/kg moist air)
!
!
! The array of the specific humidity mixing ratio (kg/kg)
!
!
! The array of the restoring tendency of the carbon dioxide mixing ratio.
!
!
subroutine atmos_co2_sourcesink(Time, dt, pwt, co2, sphum, co2_restore)
type (time_type), intent(in) :: Time
real, intent(in) :: dt
real, intent(in), dimension(:,:,:) :: pwt ! kg/m2
real, intent(in), dimension(:,:,:) :: co2 ! moist mmr
real, intent(in), dimension(:,:,:) :: sphum
real, intent(out), dimension(:,:,:) :: co2_restore
!
!-----------------------------------------------------------------------
! local parameters
!-----------------------------------------------------------------------
!
character(len=64), parameter :: sub_name = 'atmos_co2_sourcesink'
character(len=256), parameter :: error_header = &
'==>Error from ' // trim(mod_name) // '(' // trim(sub_name) // '):'
character(len=256), parameter :: warn_header = &
'==>Warning from ' // trim(mod_name) // '(' // trim(sub_name) // '):'
character(len=256), parameter :: note_header = &
'==>Note from ' // trim(mod_name) // '(' // trim(sub_name) // '):'
integer :: i,j,k,id,jd,kd, logunit
logical :: sent
!-----------------------------------------------------------------------
id=size(co2,1); jd=size(co2,2); kd=min(size(co2,3),restore_klimit)
co2_restore(:,:,:)=0.0
logunit=stdlog()
if (ind_co2 > 0 .and. do_co2_restore) then
! input is in dvmr (mol/mol)
call data_override('ATM', 'co2_dvmr_restore', restore_co2_dvmr, Time, override=used)
if (.not. used) then
call error_mesg (trim(error_header), ' data override needed for co2_dvmr_restore ', FATAL)
endif
! if (mpp_pe() == mpp_root_pe() ) &
! write (logunit,*)' atmos_co2_sourcesink: mean restore co2_dvmr = ', restore_co2_dvmr
if (restore_tscale .gt. 0 .and. restore_co2_dvmr .ge. 0.0) then
! co2mmr = (wco2/wair) * co2vmr; wet_mmr is approximated as dry_mmr * (1-Q)
do k=1,kd
do j=1,jd
do i=1,id
! convert restore_co2_dvmr to wet mmr and get tendency
co2_restore(i,j,k) = (restore_co2_dvmr * (WTMCO2/WTMAIR) * (1.0 - &
sphum(i,j,k)) - co2(i,j,k))/restore_tscale
enddo
enddo
enddo
! restoring diagnostic in moles co2/m2/sec
! pwt is moist air, so no need to divide by 1-sphum here
if (id_co2restore > 0) sent = send_data (id_co2restore, co2_restore * &
pwt / (WTMCO2*1.e-3), Time)
endif
!else
! if (mpp_pe() == mpp_root_pe() ) &
! write (logunit,*)' atmos_co2_sourcesink: CO2 restoring not active: ',do_co2_restore
endif
!! add pwt as a diagnostic
if (id_pwt > 0) sent = send_data (id_pwt, pwt, Time)
end subroutine atmos_co2_sourcesink
!
!#######################################################################
!
!
! Subroutine to get global avg co2 to be used in radiation.
! input co2 field is from data override
!
subroutine atmos_co2_rad(Time, radiation_co2_dvmr)
!
!-----------------------------------------------------------------------
! arguments
!-----------------------------------------------------------------------
!
type (time_type), intent(in) :: Time
real, intent(inout) :: radiation_co2_dvmr
!
!-----------------------------------------------------------------------
! local variables
!-----------------------------------------------------------------------
!
!
!-----------------------------------------------------------------------
! local parameters
!-----------------------------------------------------------------------
!
character(len=64), parameter :: sub_name = 'atmos_co2_rad'
character(len=256), parameter :: error_header = &
'==>Error from ' // trim(mod_name) // '(' // trim(sub_name) // '):'
character(len=256), parameter :: warn_header = &
'==>Warning from ' // trim(mod_name) // '(' // trim(sub_name) // '):'
character(len=256), parameter :: note_header = &
'==>Note from ' // trim(mod_name) // '(' // trim(sub_name) // '):'
integer :: logunit
!-----------------------------------------------------------------------
logunit=stdlog()
if (ind_co2 > 0 .and. co2_radiation_override) then
! input is in dvmr (mol/mol)
call data_override('ATM', 'co2_dvmr_rad', radiation_co2_dvmr, Time, override=used)
if (.not. used) then
call error_mesg (trim(error_header), ' data override needed for co2_dvmr_rad ', FATAL)
endif
! if (mpp_pe() == mpp_root_pe() ) &
! write (logunit,*)' atmos_co2_rad : mean radiation co2_dvmr = ', radiation_co2_dvmr
!else
! if (mpp_pe() == mpp_root_pe() ) &
! write (logunit,*)' atmos_co2_rad: CO2 radiation override not active: ',co2_radiation_override
endif
!-----------------------------------------------------------------------
end subroutine atmos_co2_rad
!
!#######################################################################
!
!
! A subroutine to gather fields needed for calculating the CO2 gas flux
!
!
subroutine atmos_co2_gather_data (gas_fields, tr_bot)
use coupler_types_mod, only: coupler_2d_bc_type, ind_pcair
implicit none
!-----------------------------------------------------------------------
type(coupler_2d_bc_type), intent(inout) :: gas_fields
real, dimension(:,:,:), intent(in) :: tr_bot
!
!-----------------------------------------------------------------------
! local parameters
!-----------------------------------------------------------------------
!
character(len=64), parameter :: sub_name = 'atmos_co2_gather_data'
character(len=256), parameter :: error_header = &
'==>Error from ' // trim(mod_name) // '(' // trim(sub_name) // '):'
character(len=256), parameter :: warn_header = &
'==>Warning from ' // trim(mod_name) // '(' // trim(sub_name) // '):'
character(len=256), parameter :: note_header = &
'==>Note from ' // trim(mod_name) // '(' // trim(sub_name) // '):'
!-----------------------------------------------------------------------
!-----------------------------------------------------------------------
! 2008/06/17 JPD/jgj: OCMIP calculation expects pco2 in dry vmr (mol/mol) units
! atm co2 is in moist mass mixing ratio (kg co2/kg moist air)
! tr_bot: co2 bottom layer moist mass mixing ratio
! convert to dry_mmr and then to dry_vmr for ocean model.
! dry_mmr = wet_mmr / (1-Q); co2vmr = (wair/wco2) * co2mmr
if (ind_co2_flux .gt. 0) then
gas_fields%bc(ind_co2_flux)%field(ind_pcair)%values(:,:) = (tr_bot(:,:,ind_co2) / &
(1.0 - tr_bot(:,:,ind_sphum))) * (WTMAIR/gas_fields%bc(ind_co2_flux)%mol_wt)
endif
end subroutine atmos_co2_gather_data
!
!#######################################################################
!#######################################################################
!
!
! Subroutine to initialize the carbon dioxide flux
!
subroutine atmos_co2_flux_init
use atmos_ocean_fluxes_mod, only: aof_set_coupler_flux
!
!-----------------------------------------------------------------------
! arguments
!-----------------------------------------------------------------------
!
!
!-----------------------------------------------------------------------
! local variables
!-----------------------------------------------------------------------
!
integer :: n
!
!-----------------------------------------------------------------------
! local parameters
!-----------------------------------------------------------------------
!
character(len=64), parameter :: sub_name = 'atmos_co2_flux_init'
character(len=256), parameter :: error_header = &
'==>Error from ' // trim(mod_name) // '(' // trim(sub_name) // '):'
character(len=256), parameter :: warn_header = &
'==>Warning from ' // trim(mod_name) // '(' // trim(sub_name) // '):'
character(len=256), parameter :: note_header = &
'==>Note from ' // trim(mod_name) // '(' // trim(sub_name) // '):'
integer :: logunit
if ( .not. module_is_initialized) then
!----- set initial value of carbon ------------
call tracer_manager_init ! need to call here since the ocean pes never call it
n = get_tracer_index(MODEL_ATMOS,'co2')
if (n > 0) then
ind_co2 = n
if (ind_co2 > 0) then
logunit=stdout()
write (logunit,*) trim(note_header), ' CO2 was initialized as tracer number ', ind_co2
logunit=stdlog()
write (logunit,*) trim(note_header), ' CO2 was initialized as tracer number ', ind_co2
endif
endif
module_is_initialized = .TRUE.
endif
!
! initialize coupler flux
!
if (ind_co2 > 0) then
ind_co2_flux = aof_set_coupler_flux('co2_flux', &
flux_type = 'air_sea_gas_flux', implementation = 'ocmip2', &
atm_tr_index = ind_co2, &
mol_wt = WTMCO2, param = (/ 9.36e-07, 9.7561e-06 /), &
caller = trim(mod_name) // '(' // trim(sub_name) // ')')
endif
!-----------------------------------------------------------------------
end subroutine atmos_co2_flux_init
!
!#######################################################################
!
!
! Subroutine to initialize the carbon dioxide module.
!
subroutine atmos_co2_init (Time, axes)
!
!-----------------------------------------------------------------------
! arguments
!-----------------------------------------------------------------------
!
type(time_type), intent(in) :: Time
integer, dimension(3), intent(in) :: axes
!
!-----------------------------------------------------------------------
! local variables
! unit io unit number used to read namelist file
! ierr error code
! io error status returned from io operation
!-----------------------------------------------------------------------
!
integer :: ierr, unit, io, logunit
integer :: n
real :: missing_value = -1.e10
character(len=64) :: desc
!
!-----------------------------------------------------------------------
! local parameters
!-----------------------------------------------------------------------
!
character(len=64), parameter :: sub_name = 'atmos_co2_init'
character(len=256), parameter :: error_header = &
'==>Error from ' // trim(mod_name) // '(' // trim(sub_name) // '):'
character(len=256), parameter :: warn_header = &
'==>Warning from ' // trim(mod_name) // '(' // trim(sub_name) // '):'
character(len=256), parameter :: note_header = &
'==>Note from ' // trim(mod_name) // '(' // trim(sub_name) // '):'
if (module_is_initialized) return
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_co2_nml, iostat=io, end=10)
ierr = check_nml_error(io,'atmos_co2_nml')
end do
10 call close_file (unit)
endif
!---------------------------------------------------------------------
! write namelist to logfile.
!---------------------------------------------------------------------
logunit=stdlog()
if (mpp_pe() == mpp_root_pe() ) &
write (logunit, nml=atmos_co2_nml)
!----- set initial value of carbon ------------
n = get_tracer_index(MODEL_ATMOS,'co2')
if (n > 0) then
ind_co2 = n
logunit=stdout()
write (logunit,*) trim(note_header), ' CO2 was initialized as tracer number ', ind_co2
logunit=stdlog()
write (logunit,*) trim(note_header), ' CO2 was initialized as tracer number ', ind_co2
! initialize diagnostics
desc = ' restoring tendency'
id_co2restore = register_diag_field ('atmos_co2_restoring', 'co2_restore', axes, Time, &
'CO2'//trim(desc), 'moles co2/m2/s',missing_value=missing_value)
desc = ' pressure weighting array = dP/grav'
id_pwt = register_diag_field ('atmos_co2', 'pwt', axes, Time, &
trim(desc), 'kg/m2',missing_value=missing_value)
!
! get the index for sphum
!
ind_sphum = get_tracer_index(MODEL_ATMOS,'sphum')
if (ind_sphum .le. 0) then
call error_mesg (trim(error_header), ' Could not find index for sphum', FATAL)
endif
endif
logunit=stdlog()
if (.not.(ind_co2 > 0 .and. do_co2_restore)) then
if (mpp_pe() == mpp_root_pe() ) &
write (logunit,*)' CO2 restoring not active: ',do_co2_restore
endif
if (.not.(ind_co2 > 0 .and. co2_radiation_override)) then
if (mpp_pe() == mpp_root_pe() ) &
write (logunit,*)' CO2 radiation override not active: ',co2_radiation_override
endif
call write_version_number (version, tagname)
module_is_initialized = .TRUE.
!-----------------------------------------------------------------------
end subroutine atmos_co2_init
!
!
subroutine atmos_co2_end
!
!-----------------------------------------------------------------------
! local parameters
!-----------------------------------------------------------------------
!
character(len=64), parameter :: sub_name = 'atmos_co2_end'
character(len=256), parameter :: error_header = &
'==>Error from ' // trim(mod_name) // '(' // trim(sub_name) // '):'
character(len=256), parameter :: warn_header = &
'==>Warning from ' // trim(mod_name) // '(' // trim(sub_name) // '):'
character(len=256), parameter :: note_header = &
'==>Note from ' // trim(mod_name) // '(' // trim(sub_name) // '):'
module_is_initialized = .FALSE.
end subroutine atmos_co2_end
!
end module atmos_co2_mod