!----------------------------------------------------------------------- ! purpose ! ------- ! perform all calculations needed for TM5 test run ! with 1 ppb concentration ! ! interface ! --------- ! call trace0 ! call trace1 ! call source ! call chemie ! ! method ! ------ ! subroutine trace0 is called from start, and is used ! to input tracer specific values. ! subroutine trace1 is called from start when istart = 2, and provides ! a way to do some initialization, that has to be done only ! at the beginning of a run. ! subroutine source is called every nsrce seconds and provides ! the source/sink calculation. ! subroutine chem is called every nchem seconds and provides the ! chemistry calculation on the tracer fields. ! in zoom mode, recursive calling of the children is performed ! !### macro's ##################################################### ! #define TRACEBACK write (gol,'("in ",a," (",a,", line",i5,")")') 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" ! !################################################################# module sources_sinks use GO , only : gol, goErr, goPr use emission , only : Emis_Init, Emis_Done use global_data , only : outdir implicit none ! --- in/out -------------------------- private public :: Sources_sinks_Init, Sources_sinks_Done public :: Trace0, Trace1, Trace_End, trace_after_read public :: Source1 public :: Source2 ! --- const -------------------------------------- character(len=*), parameter :: mname = 'sources_sinks' contains ! ================================================================ subroutine Sources_Sinks_Init( status ) ! --- in/out -------------------------------- integer, intent(out) :: status ! --- const ------------------------------ character(len=*), parameter :: rname = mname//'/Sources_Sinks_Init' ! --- begin -------------------------------- ! nothing to be init ! ok call Emis_Init(status) IF_NOTOK_RETURN(status=1) status = 0 end subroutine Sources_Sinks_Init ! *** subroutine Sources_Sinks_Done( status ) ! --- in/out -------------------------------- integer, intent(out) :: status ! --- const ------------------------------ character(len=*), parameter :: rname = mname//'/Sources_Sinks_Done' ! --- begin -------------------------------- ! nothing to be done ! ok call Emis_Done(trim(outdir)//'/budget_global.hdf',status) status = 0 end subroutine Sources_Sinks_Done ! ================================================================ !----------------------------------------------------------------------- ! ! trace0 ! ! ! purpose ! ------- ! initialise data needed for chemistry, ! e.g. read emissions and boundary conditions ! routine is calculated at start and at beginning each month... ! ! interface ! --------- ! call trace0 ! ! method ! ------ ! subroutine trace0 is called from module input, and provides a way ! to input tracer specific values. ! ! external ! -------- ! calc_sm ! ini_zone ! ini_zoneg ! call_diagbud ! call_diagbudg ! ! reference ! --------- ! see above !--------------------------------------------------------------------- subroutine trace0( status ) #ifndef without_photolysis use photolysis , only : init_photolysis #endif use datetime, only : calc_sm use dims, only : mlen, sec_month, sec_year, sec_day use dims, only : okdebug ! --- in/out -------------------------- integer, intent(out) :: status ! --- local -------------------------- ! --- const ------------------------------ character(len=*), parameter :: rname = mname//'/Trace0' ! --- begin --------------------------- #ifndef without_photolysis call init_photolysis !( kr, o3du ) !if ( okdebug ) then ! write (gol,*) 'On-line photolysis rate calculation initialized', myid; call goPr !end if #endif write (gol,'(a,": time parameters ...")') rname; call goPr call calc_sm( mlen, sec_day,sec_month,sec_year ) if ( okdebug ) then write (gol,*) 'seconds in a day, month and year',sec_day,sec_month,sec_year; call goPr end if ! ok status = 0 end subroutine trace0 subroutine trace1 !----------------------------------------------------------------------- ! ! purpose ! ------- ! this subroutine determines the initial values for tracer mass and ! and its slopes. it is called from input when istart = 2. ! ! interface ! --------- ! call trace0 ! ! method ! ------ ! ! reference ! --------- ! see above ! !----------------------------------------------------------------------- use dims, only : im, jm, lm use dims, only : nregions use toolbox, only : print_totalmass, escape_tm use global_data, only : mass_dat use meteodata , only : m_dat use tracer_data, only : chem_dat use zoom_tools, only : update_parent use chem_param, only : fscale, ntracet, ntrace, ntrace_chem use ParTools , only : myid use partools , only : lmar, lmloc, offsetl use partools , only : ntracetloc, ntracet_ar !_______________________start_____________________- character(len=*), parameter :: rname = mname//'/trace1' ! local real,dimension(:,:,:,:),pointer :: rm #ifdef slopes real,dimension(:,:,:,:),pointer :: rxm, rym, rzm #endif real,dimension(:,:,:), pointer :: m real,dimension(:,:,:,:),pointer :: rmc integer :: region, n, l, offsetn ! start if ( ntracetloc == 0 ) return !WP! only PE's with nonzero ntracetloc do region=1,nregions m => m_dat(region)%data rm => mass_dat(region)%rm_t #ifdef slopes rxm => mass_dat(region)%rxm_t rym => mass_dat(region)%rym_t rzm => mass_dat(region)%rzm_t #endif offsetn = 0 if ( myid > 0 ) offsetn = sum(ntracet_ar(0:myid-1)) do n=1,ntracetloc rm(1:im(region),1:jm(region),1:lm(region),n) = & (n+offsetn)*1e-9*m(1:im(region),1:jm(region),1:lm(region))/ & fscale(n+offsetn) if((n+offsetn) == 5) rm(1:im(region),1:jm(region),1:lm(region),n) = 0.0 ! rm(1:10,1:10,1:10,n) = 40e-9*m(1:10,1:10,1:10)/fscale(n) #ifdef slopes rxm(:,:,:,n) = 0.0 rym(:,:,:,n) = 0.0 rzm(:,:,:,n) = 0.0 #endif end do ! non-transported tracers ? if ( ntrace_chem > 0 ) then ! for safety ... if ( .not. associated(chem_dat(region)%rm_k) ) then write (*,'("ERROR - array chem_dat(region)%rm_k not allocated ...")') write (*,'("ERROR in ",a)') 'sourcec_sinks/trace1' call escape_tm('Error in trace1') !allocate( chem_dat(region)%rm_k(im(region),jm(region),lmloc, & ! ntracet+1:ntracet+ntrace_chem) ) end if rmc => chem_dat(region)%rm_k do n=ntracet+1,ntrace do l=1,lmloc rmc(1:im(region),1:jm(region),l,n) = & n*1e-9*m(1:im(region),1:jm(region),l+offsetl)/fscale(n) end do end do end if ! non transported tracers nullify(m) nullify(rm) #ifdef slopes nullify(rxm) nullify(rym) nullify(rzm) #endif nullify(rmc) end do write (gol,*) 'trace1: rm initialized at mixing ratio :', (myid+1)*1e-9; call goPr do region = nregions,2,-1 call update_parent(region) end do call print_totalmass(1) end subroutine trace1 !--------------------------------------------------------------------- ! ! subroutine that is called after reading new fields ! (clouds, surface winds, etc.) ! In this routine, 'chemistry' fields that depend on these data are calculated ! Routine is called from tracer.f90 (after reading data) ! !--------------------------------------------------------------------- subroutine trace_after_read( status ) ! --- in/out --------------------------------- integer, intent(out) :: status ! --- const --------------------------------- character(len=*), parameter :: rname = 'trace_after_read' ! --- begin --------------------------------- ! nothing to be done ! ok status = 0 end subroutine trace_after_read !--------------------------------------------------------------------- ! ! source1 ! ! ! purpose ! ------- ! this subroutine changes the tracer mass and its slopes ! by chemical sources ! ! interface ! --------- ! call source1 ! ! method ! ------ ! ! ! ! reference ! -------- ! see above ! !--------------------------------------------------------------------- subroutine Source1( region, tr, status ) use GO , only : TDate use emission, only: emission_apply ! use GO, only : gol, goPr, goErr, goBug, goLabel ! --- in/out --------------------------------- integer, intent(in) :: region type(TDate), intent(in) :: tr(2) ! time range integer, intent(out) :: status ! --- const --------------------------------------- character(len=*), parameter :: rname = mname//'/Source1' ! --- begin ---------------------------------- call emission_apply(region, status) IF_NOTOK_RETURN(status=1) ! ok ! status = 0 end subroutine source1 ! --------------------------------------------------------------------- ! ! source2 ! ! ! purpose ! ------- ! apply boundary conditions to selected tracers ! ! interface ! --------- ! call source2 ! ! method ! ------ ! ! ! external ! -------- ! none ! ! reference ! -------- ! see above ! ! !--------------------------------------------------------------------- subroutine source2( region, tr, status ) use GO , only : TDate use dims, only : dx, dy, xref, yref, tref, nsrce use toolbox, only : dumpfieldi,dumpfield use global_data, only : mass_dat,region_dat ! --- in/out -------------------------------- integer, intent(in) :: region type(TDate), intent(in) :: tr(2) ! time range integer, intent(out) :: status ! --- const -------------------------------- character(len=*), parameter :: rname = mname//'/source2' ! --- local --------------------------------- ! --- begin ---------------------------------- ! ok status = 0 end subroutine source2 ! ! deallocate the memory and do everything needed to finalise chemistry ! subroutine trace_end( status ) #ifndef without_photolysis use photolysis, only : end_photolysis #endif ! --- in/out ------------------------------------ integer, intent(out) :: status ! --- local ------------------------------------- ! --- begin -------------------------------------- #ifndef without_photolysis call end_photolysis #endif ! ok status = 0 end subroutine trace_end end module sources_sinks