#define TRACEBACK write (gol,'("in ",a," (",a,i6,")")') rname, __FILE__, __LINE__ ; call goErr #define IF_NOTOK_RETURN(action) if (status/=0) then; TRACEBACK; action; return; end if #define IF_ERROR_RETURN(action) if (status> 0) then; TRACEBACK; action; return; end if ! #include "tm5.inc" ! !------------------------------------------------------------------------------ ! TM5 ! !------------------------------------------------------------------------------ !BOP ! ! !MODULE: SOURCES_SINKS ! ! !DESCRIPTION: Perform all calculations needed for CBM4 chemistry simulation ! in TM5: this is mainly emissions, process updates after ! changes in meteo, boundary, sedimentation, photolysis,... ! !FD: all emission are converted to kg X (fmw) /month exceptions are mentioned in the code ! !\\ !\\ ! !INTERFACE: ! MODULE SOURCES_SINKS ! ! !USES: ! use GO, only : gol, goErr, goPr, goBug, goLabel implicit none private ! ! !PUBLIC MEMBER FUNCTIONS: ! PUBLIC :: SOURCES_SINKS_INIT, SOURCES_SINKS_DONE ! Init and Done methods PUBLIC :: SS_MONTHLY_UPDATE ! monthly initialization (photolysis,..) PUBLIC :: SS_AFTER_READ_METEO_UPDATE ! Update SS after met fields are updated. Called from modelIntegration/Proces_update PUBLIC :: SOURCES_SINKS_APPLY ! apply SS ! ! !PRIVATE DATA MEMBERS: ! character(len=*), parameter :: mname = 'sources_sinks' ! ! !REVISION HISTORY: ! 19 Mar 2012 - P. Le Sager - adapted for lon-lat MPI domain decomposition ! ! !REMARKS: ! !EOP !------------------------------------------------------------------------------ contains !------------------------------------------------------------------------------ ! TM5 ! !------------------------------------------------------------------------------ !BOP ! ! !IROUTINE: SOURCES_SINKS_INIT ! ! !DESCRIPTION: switch ON required meteo; init emissions !\\ !\\ ! !INTERFACE: ! SUBROUTINE SOURCES_SINKS_INIT( status ) ! ! !USES: ! use meteo, only : Set use meteodata, only : temper_dat, humid_dat, oro_dat, gph_dat use Meteodata, only : cp_dat, lsp_dat use Meteodata, only : cvl_dat, cvh_dat, tv_dat use Meteodata, only : ci_dat, sd_dat, swvl1_dat use Meteodata, only : t2m_dat, d2m_dat use Meteodata, only : lsmask_dat, ci_dat, sst_dat #if defined(with_online_bvoc) || defined(with_online_nox) use Meteodata, only : skt_dat #endif #ifdef with_online_bvoc use emission_bvoc_data, only : online_bvoc_skt #endif #ifdef with_online_nox use online_nox_data, only : online_nox_skt use Meteodata, only : src_dat, lsmask_dat #endif use Meteodata, only : ssr_dat, sshf_dat, slhf_dat, ewss_dat, nsss_dat use Meteodata, only : src_dat, albedo_dat, nveg use chem_rates, only : rates #ifndef without_emission use emission, only : Emission_Init #endif #ifndef without_sedimentation use sedimentation , only : Sedimentation_Init #endif #ifndef without_photolysis use photolysis , only : Photolysis_Init #endif #ifndef without_boundary use boundary , only : Boundary_Init, o3du, use_o3du #endif #ifdef with_m7 use emission_dust, only : emission_dust_init use emission_ss, only : emission_ss_init #endif use GO, only : TrcFile, Init, Done, ReadRc use global_data, only : rcfile use dims, only : iglbsfc, nregions ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !REVISION HISTORY: ! 3 Oct 2012 - P. Le Sager - get Henry coeff (call rates) ! ! !REMARKS: ! !EOP !------------------------------------------------------------------------------ !BOC integer, parameter :: kr=31 ! standard unit to read auxiliary files character(len=*), parameter :: rname = mname//'/Sources_Sinks_Init' type(TrcFile) :: rcF integer :: region, iveg ! --- begin --------------------------------- !-------------------------------------------------- ! ** select meteo (cases not accounted for in the **_init procedures) !-------------------------------------------------- do region = 1, nregions #ifndef without_emission call Set( temper_dat(region), status, used=.true. ) call Set( humid_dat(region), status, used=.true. ) call Set( oro_dat(region), status, used=.true. ) call Set( gph_dat(region), status, used=.true. ) #endif ! other call Set( cvl_dat(region), status, used=.true. ) call Set( cvh_dat(region), status, used=.true. ) do iveg=1,nveg call Set( tv_dat(region,nveg), status, used=.true. ) enddo call Set( ci_dat(region), status, used=.true. ) call Set( sd_dat(region), status, used=.true. ) call Set( swvl1_dat(region), status, used=.true. ) call Set( t2m_dat(region), status, used=.true. ) call Set( d2m_dat(region), status, used=.true. ) call Set( ssr_dat(region), status, used=.true. ) call Set( sshf_dat(region), status, used=.true. ) call Set( slhf_dat(region), status, used=.true. ) call Set( ewss_dat(region), status, used=.true. ) call Set( nsss_dat(region), status, used=.true. ) call Set( src_dat(region), status, used=.true. ) call Set( albedo_dat(region), status, used=.true. ) enddo !-------------------------------------------------- ! ** Henry coefficients (must be before sedimentation_init) !-------------------------------------------------- call rates(status) IF_NOTOK_RETURN(status=1) !-------------------------------------------------- ! ** Sedimentation !-------------------------------------------------- #ifndef without_sedimentation call Sedimentation_Init( status ) IF_NOTOK_RETURN(status=1) #endif !-------------------------------------------------- ! ** Stratospheric boundary (must be before photolysis) !-------------------------------------------------- #ifndef without_boundary call Boundary_Init( .true., status ) IF_NOTOK_RETURN(status=1) #endif !-------------------------------------------------- ! ** Photolysis !-------------------------------------------------- #ifndef without_photolysis #ifdef without_boundary call photolysis_init(.true., kr ) #else if (use_o3du) then call Photolysis_Init(.true., kr, o3du ) else call Photolysis_Init(.true., kr ) end if #endif #endif !-------------------------------------------------- ! ** Emissions !-------------------------------------------------- #ifndef without_emission call Emission_Init( status ) IF_NOTOK_RETURN(status=1) #ifdef with_online_bvoc call Init( rcF, rcfile, status ) IF_NOTOK_RETURN(status=1) call ReadRc( rcF, 'online.bvoc.skt', online_bvoc_skt, status ) IF_NOTOK_RETURN(status=1) call Done( rcF, status ) IF_NOTOK_RETURN(status=1) if (online_bvoc_skt) then call Set( skt_dat(iglbsfc), status, used=.true. ) else call Set( t2m_dat(iglbsfc), status, used=.true. ) endif call Set( ssr_dat(iglbsfc), status, used=.true. ) #endif #ifdef with_online_nox call Init( rcF, rcfile, status ) IF_NOTOK_RETURN(status=1) call ReadRc( rcF, 'online.nox.skt', online_nox_skt, status ) IF_NOTOK_RETURN(status=1) call Done( rcF, status ) IF_NOTOK_RETURN(status=1) if (online_nox_skt) then call Set( skt_dat(iglbsfc), status, used=.true. ) else call Set( t2m_dat(iglbsfc), status, used=.true. ) endif call Set( lsp_dat(iglbsfc), status, used=.true. ) call Set( cp_dat(iglbsfc), status, used=.true. ) call Set( src_dat(iglbsfc), status, used=.true. ) call Set( lsmask_dat(iglbsfc), status, used=.true. ) #endif #endif /* EMISSIONS */ #ifdef with_m7 call emission_dust_init( status ) call emission_ss_init( status ) #endif !-------------------------------------------------- ! ** Done !-------------------------------------------------- status = 0 END SUBROUTINE SOURCES_SINKS_INIT !EOC !------------------------------------------------------------------------------ ! TM5 ! !------------------------------------------------------------------------------ !BOP ! ! !IROUTINE: SOURCES_SINKS_DONE ! ! !DESCRIPTION: !\\ !\\ ! !INTERFACE: ! SUBROUTINE SOURCES_SINKS_DONE( status ) ! ! !USES: ! #ifndef without_photolysis use photolysis, only: photolysis_done #endif #ifndef without_sedimentation use sedimentation, only: Sedimentation_Done #endif #ifndef without_boundary use Boundary, only: Boundary_Done #endif #ifndef without_emission use emission, only: Emission_Done #endif ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !REVISION HISTORY: ! !EOP !------------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/Sources_Sinks_Done' ! --- begin -------------------------------- #ifndef without_photolysis call photolysis_done ( status ) IF_NOTOK_RETURN(status=1) #endif #ifndef without_boundary call Boundary_Done( status ) IF_NOTOK_RETURN(status=1) #endif #ifndef without_sedimentation call Sedimentation_Done( status ) IF_NOTOK_RETURN(status=1) #endif #ifndef without_emission call Emission_Done( status ) IF_NOTOK_RETURN(status=1) #endif status = 0 END SUBROUTINE SOURCES_SINKS_DONE !EOC !------------------------------------------------------------------------------ ! TM5 ! !------------------------------------------------------------------------------ !BOP ! ! !IROUTINE: SS_MONTHLY_UPDATE ! ! !DESCRIPTION: monthly (re)initialisation of sources/sinks !\\ !\\ ! !INTERFACE: ! SUBROUTINE SS_MONTHLY_UPDATE( status ) ! ! !USES: ! use dims, only : mlen, sec_day, sec_month, sec_year use datetime, only : calc_sm #ifndef without_photolysis use photolysis , only : photolysis_init #endif #ifndef without_emission use emission, only : declare_emission #endif #ifndef without_boundary use boundary, only : Boundary_Init, o3du, use_o3du #endif ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !REVISION HISTORY: ! 19 Mar 2012 - P. Le Sager - adapted for lon-lat MPI domain decomposition ! ! !REMARKS: ! (1) routine is called at start and at beginning of each month ! !EOP !------------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/ss_monthly_update' integer, parameter :: kr=31 ! standard unit to read auxiliary files ! --- begin ------------------------------------ ! calculate some conversion factors related to time... call calc_sm( mlen, sec_day, sec_month, sec_year ) ! Read monthly emissions #ifndef without_emission call declare_emission( status ) IF_NOTOK_RETURN(status=1) #endif ! Monthly update for stratospheric boundary #ifndef without_boundary call Boundary_Init( .false., status ) IF_NOTOK_RETURN(status=1) #endif ! Monthly update for photolysis #ifndef without_photolysis #ifndef without_boundary if (use_o3du) then call Photolysis_Init( .false., kr, o3du ) end if #endif #endif status = 0 END SUBROUTINE SS_MONTHLY_UPDATE !EOC !------------------------------------------------------------------------------ ! TM5 ! !------------------------------------------------------------------------------ !BOP ! ! !IROUTINE: SS_AFTER_READ_METEO_UPDATE ! ! !DESCRIPTION: subroutine that is called after reading new met fields (clouds, ! surface winds, etc.). ! In this routine, 'chemistry' fields that depend on these ! data are calculated. Called from modelIntegration/Proces_update. !\\ !\\ ! !INTERFACE: ! SUBROUTINE SS_AFTER_READ_METEO_UPDATE( status ) ! ! !USES: ! use dims, only : nregions, sec_month use tm5_distgrid, only : dgrid, Get_DistGrid #ifndef without_photolysis use photolysis, only : ozone_info_online, slingo, aerosol_info, update_csqy #endif #ifndef without_emission use emission_nox, only : eminox_lightning #ifndef without_convection use emission_nox, only : lightningNOX #endif use emission_dms, only : getDMS #if defined(with_online_bvoc) || defined(with_online_nox) use dims, only : itau, ndyn_max #endif #ifdef with_online_bvoc use emission_bvoc, only : getBVOC #endif #ifdef with_online_nox use online_nox, only : getNOx use emission_nox, only : nat_nox #endif #ifdef with_m7 use emission_dust, only : calc_emission_dust, read_emission_dust use emission_ss, only : calc_emission_ss #endif #endif /* EMISSIONS */ #ifndef without_sedimentation use sedimentation, only : calculate_rh #endif ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !REVISION HISTORY: ! 19 Mar 2012 - P. Le Sager - adapted for lon-lat MPI domain decomposition ! ! !REMARKS: ! !EOP !------------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = 'ss_after_read_meteo_update' integer :: region, i1, j1 ! --- begin ------------------------------------ call goLabel() #ifndef without_emission #ifdef with_online_nox if (mod(itau, ndyn_max) == 0) then call getNOx(status) ! returns nat_nox in kg(N)/month IF_NOTOK_RETURN(status=1) endif #endif !get dms_emissions dms_sea, flagged by new_surface meteo fields get in kgS/month call getDMS( status ) IF_NOTOK_RETURN(status=1) #ifdef with_online_bvoc if (mod(itau, ndyn_max) == 0) then ! TvN: added fixed time step ndyn_max (1 hour) ! because the averages over the past 24 hours are calculated from hourly averages call getBVOC (status) IF_NOTOK_RETURN(status=1) endif #endif ! Lightning NOx (defined only if convection is turned on) #ifdef without_convection do region = 1, nregions eminox_lightning(region)%d3=0. end do #else do region = 1, nregions call Get_DistGrid( dgrid(region), I_STRT=i1, J_STRT=j1 ) call lightningNOX(region, i1, j1, eminox_lightning(region)%d3, status) IF_NOTOK_RETURN(status=1) eminox_lightning(region)%d3(:,:,:) = eminox_lightning(region)%d3(:,:,:)*sec_month !from kg N/s ----> kg N/month end do #endif #ifdef with_m7 call calc_emission_ss( status ) IF_NOTOK_RETURN(status=1) call read_emission_dust( status ) ! this is active if (input.emis.dust /= ONLINE) IF_NOTOK_RETURN(status=1) call calc_emission_dust( status ) ! this is active if (input.emis.dust == ONLINE) IF_NOTOK_RETURN(status=1) #endif #endif /* EMISSIONS */ #ifndef without_photolysis do region = 1, nregions ! cloud optical depth call slingo(region) ! calculate optical depth ozone from current ozone field ! note: this routine does not depend on clouds but on rm ... call ozone_info_online(region) call update_csqy( region ) ! t/p-dependent cross-sections and quantum yields end do #endif #ifndef without_sedimentation call calculate_rh #endif ! ok call goLabel() status = 0 END SUBROUTINE SS_AFTER_READ_METEO_UPDATE !EOC !------------------------------------------------------------------------------ ! TM5 ! !------------------------------------------------------------------------------ !BOP ! ! !IROUTINE: SOURCES_SINKS_APPLY ! ! !DESCRIPTION: this subroutine changes the tracer mass and its ! slopes by chemical sources. !\\ !\\ ! !INTERFACE: ! SUBROUTINE SOURCES_SINKS_APPLY( region, tr, status ) ! ! !USES: ! use GO, only : TDate #ifndef without_emission use emission, only: emission_apply #endif #ifndef without_sedimentation use sedimentation, only: Sedimentation_Apply #endif #ifndef without_boundary use boundary, only : Boundary_Apply #endif ! ! !INPUT PARAMETERS: ! integer, intent(in) :: region type(TDate) :: tr(2) ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !REVISION HISTORY: ! ! !REMARKS: ! - called each time step, during "source" step, by modelIntegration/do_steps ! !EOP !------------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/Sources_sinks_apply' ! --- begin ---------------------------------- #ifndef without_sedimentation call Sedimentation_Apply ( region, status ) IF_NOTOK_RETURN(status=1) #endif #ifndef without_emission ! note dust/ss emissions are ported to sedimentation routine call emission_apply( region, status ) IF_NOTOK_RETURN(status=1) #endif ! Apply boundary conditions to selected tracers #ifndef without_chemistry #ifndef without_boundary call Boundary_Apply( region, status ) IF_NOTOK_RETURN(status=1) #endif #endif status = 0 END SUBROUTINE SOURCES_SINKS_APPLY !EOC END MODULE SOURCES_SINKS