!### 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 emission !BOP !----------------------------------------------------------------------- ! purpose ! ------- ! perform emissions needed for TM5 CBM4 version ! ! interface ! --------- ! call declare_emission ! call emission_apply ! call free_emission ! ! method ! ------ ! subroutine declare_emission is called from trace0 ! subroutine emission_apply is called from source1 ! subroutine free_emission is called from trace_end !----------------------------------------------------------------------- ! Table of Contents: ! ! module emission ! subroutine Emis_Init(status) ! end subroutine Emis_Init ! subroutine emission_apply(region) ! end subroutine emission_apply ! subroutine Emis_Done(budget_file,status) ! end subroutine Emis_Done ! end module emission !EOP use GO, only: goLoCase,gol, goErr, goPr, goLabel #ifdef with_budgets use emission_data, only: budemi_dat #endif use emission_data, only: plandr, emis2D use emission_data, only: sum_emission use emission_common, only: declare_emission_averagers, free_emission_averagers, calc_emission_averagers, save_emission_averagers use emission_common, only: EmisInputDir, ParamInputFile, C13InputDir,GetParameterMaps use emission_co2_ff, only: declare_emission_co2_ff, free_emission_co2_ff, calc_emission_co2_ff, emission_apply_co2_ff use emission_co2_ocean, only: declare_emission_co2_ocean, free_emission_co2_ocean, calc_emission_co2_ocean, emission_apply_co2_ocean use emission_co2_bio, only: declare_emission_co2_bio, free_emission_co2_bio, calc_emission_co2_bio, emission_apply_co2_bio use emission_co2_fires, only: declare_emission_co2_fires, free_emission_co2_fires, calc_emission_co2_fires, emission_apply_co2_fires use emission_co2_all, only: declare_emission_co2_all, free_emission_co2_all, emission_apply_co2_all, calc_emission_co2_all use emission_c13_bio, only: declare_emission_c13_bio, free_emission_c13_bio, calc_emission_c13_bio, emission_apply_c13_bio use emission_c13_ocean, only: declare_emission_c13_ocean, free_emission_c13_ocean, calc_emission_c13_ocean, emission_apply_c13_ocean use emission_c13_ff, only: declare_emission_c13_ff, free_emission_c13_ff, calc_emission_c13_ff, emission_apply_c13_ff use emission_c13_fires, only: declare_emission_c13_fires, free_emission_c13_fires, calc_emission_c13_fires, emission_apply_c13_fires use emission_c13_biodis, only: declare_emission_c13_biodis, free_emission_c13_biodis, calc_emission_c13_biodis, emission_apply_c13_biodis use emission_c13_oceandis, only: declare_emission_c13_oceandis, free_emission_c13_oceandis, calc_emission_c13_oceandis, emission_apply_c13_oceandis use emission_c13_all, only: declare_emission_c13_all, free_emission_c13_all, emission_apply_c13_all, calc_emission_c13_all private ! The default inverse behavior is to use emission_averagers to ! accumulate weekly average fluxes, to be written to the ! savestate.hdf file. If the rc file specifies flux1x1 output, ! then the averaging window may be different, and the output ! is to flux1x1*hdf files. However, both facilities use ! the emission_averagers and they cannot both be turned on. If ! flux1x1 output is requested, this logical will be .true. and ! the weekly average output to savestate will be suppressed. ! Note that we cannot use value of flux1x1_data from user_output ! to determine this override, because it is manipulated in ! main__inv.F90. logical :: flux1x1_override = .False. ! public routines: public :: Emis_Init, emission_apply, Emis_Done ! --- const -------------------------------------- character(len=*), parameter :: mname = 'emission' contains subroutine Emis_Init(status) !BOP ! !ROUTINE: Emis_Init ! !DESCRIPTION: ! ! purpose ! ------- ! allocate space needed to handle the emissions ! ! interface ! --------- ! call declare_emission is called from trace0 (start and newmonth) ! ! method ! ------ ! dynamic allocation, depending on dimensions problem ! opening of daily input files ! reading monthly/yearly emissions ! !--------------------------------------------------------------------- !EOP use dims, only: im,jm,lm,newsrun,nregions, nlon360, nlat180 use chem_param, only: names,ntracet, ico2_ff, ico2_fires, ico2_oce, ico2_bio, nmembersloc #ifdef with_budgets use budget_global, only: nbud_vg #endif use emission_common, only: flux_means_are_priors, enkf_set_inverse, fluxmultiplier,fluxmultiplier_biodis,fluxmultiplier_epsbio use GO, only: goLoCase,TrcFile, Init, Done, ReadRc use global_data, only: rcfile use ParTools, only: ntracetloc,ntracet_ar,myid implicit none ! --- in/out -------------------------------- integer, intent(out) :: status ! --- const ------------------------------ character(len=*), parameter :: rname = mname//'/Emis_Init' ! local integer :: region, imr, jmr, lmr,offsetn type(TrcFile) :: rcF !WP! Set the run mode to inverse, use local tracer numbers as local ensemble !WP! size call enkf_set_inverse !print*,'enkf_set_inverse',enkf_set_inverse(),'emission_inv' nmembersloc = ntracetloc !/2 ! nmembersloc = ntracetloc/2 write(gol,*) 'CarbonTracker run mode set to INVERSE' ; call goPr allocate(fluxmultiplier(nlon360,nlat180,nmembersloc)) allocate(fluxmultiplier_biodis(nlon360,nlat180,nmembersloc)) allocate(fluxmultiplier_epsbio(nlon360,nlat180,nmembersloc)) fluxmultiplier = 1.0 fluxmultiplier_biodis = 1.0 fluxmultiplier_epsbio = 1.0 call Init( rcF, rcfile, status ) IF_NOTOK_RETURN(status=1) ! This ReadRC call commented out but not removed in case someone cares to reinstate it later. ! Today, changed enkf_set_inverse routine (in emission_common.F90) to set the default ! value of flux_means_are_priors. For the inverse code, the default value in .TRUE.; ! posterior fluxes are available from a follow-on forward run or by scaling the prior ! fluxes by the parameter values (presuming that post = prior * scalefac). ! -Andy Jacobson, 29 Sep 2007. ! ! Reinstated! 15 Nov 2008, ARJ. Defaults work just fine, no need to change the rc file. ! If you want the inverse run to have posterior fluxes in its savestate output, then ! set flux_means_are_priors : F in the rc file. Otherwise, if that key is missing in the ! rc file, or if it is set to T, then the "standard" scheme will prevail. call ReadRc( rcF, 'flux_means_are_priors', flux_means_are_priors ,status ,default=.True.) IF_NOTOK_RETURN(status=1) call ReadRc( rcF, 'emis.input.dir', EmisInputDir ,status ,default='dir.not.specified') IF_NOTOK_RETURN(status=1) call ReadRc( rcF, 'c13.input.dir', C13InputDir ,status ,default='dir.not.specified') IF_NOTOK_RETURN(status=1) call ReadRc( rcF, 'ct.params.input.file', ParamInputFile, status ,default='./') IF_NOTOK_RETURN(status=1) call Done( rcF, status ) IF_NOTOK_RETURN(status=1) !WP! Get the values for the fluxmultipliers from file call GetParameterMaps(status) IF_NOTOK_RETURN(status=1) ! (maybe) initialize the average flux accumulators for weekly output if(.not. flux1x1_override) then call declare_emission_averagers endif #ifdef with_budgets ! on all PEs to be made available: do region=1,nregions imr = im(region) ;jmr = jm(region) ; lmr = lm(region) allocate(budemi_dat(region)%emi(im(region),jm(region),nbud_vg,ntracet)) budemi_dat(region)%emi = 0.0 sum_emission(region) = 0.0 enddo #endif !offsetn = sum(ntracet_ar(0:myid-1)) !if (goLoCase(trim(names(offsetn+1)))=='co2'.or.goLoCase(trim(names(offsetn+1)))=='co2_bg') then call declare_emission_co2_ff call declare_emission_co2_fires call declare_emission_co2_ocean call declare_emission_co2_bio call declare_emission_co2_all !endif !if (goLoCase(trim(names(offsetn+1)))=='co2c13'.or.goLoCase(trim(names(offsetn+1)))=='co2c13_bg') then call declare_emission_c13_bio call declare_emission_c13_ocean call declare_emission_c13_ff call declare_emission_c13_fires call declare_emission_c13_biodis call declare_emission_c13_oceandis call declare_emission_c13_all !endif status=0 end subroutine Emis_Init ! --------------------------------- subroutine emission_apply(region, status) use Dims, only: okdebug use chem_param, only: ntracet, ico2_ff, ico2_fires, ico2_oce, ico2_bio,names use ParTools, only: ntracetloc,ntracet_ar, tracer_active, myid, root_t implicit none ! --- const --------------------------------------- character(len=*), parameter :: rname = mname//'/emission_apply' integer :: offsetn integer, intent(in) :: region integer, intent(out) :: status if(okdebug) write(*,*) 'start of emission_apply' if(ntracetloc == 0) return !if no tracer on this processor do nothing if(region == 1) then !WP! only needed once, on 1x1 grid, coarsened array used after that call calc_emission_co2_ff call calc_emission_co2_bio call calc_emission_co2_fires call calc_emission_co2_ocean call calc_emission_c13_bio call calc_emission_c13_ocean call calc_emission_c13_ff call calc_emission_c13_fires call calc_emission_c13_biodis call calc_emission_c13_oceandis endif call calc_emission_co2_all call emission_apply_co2_all(region) call calc_emission_c13_all call emission_apply_c13_all(region) if(okdebug) write(*,*) 'End of adding emission ' !ok status = 0 end subroutine emission_apply subroutine Emis_Done(budget_file,status) ! ! purpose ! ------- ! deallocate space needed to handle the emissions ! ! interface ! --------- ! call free_emission is called from trace_end ! ! method ! ------ ! collecting budgets on coarser budget regions ! writing budget file ! deallocation of memory ! closing of input files !--------------------------------------------------------------------- use dims, only : im,jm, nregions use chem_param, only: names,ntracet, ico2_ff, ico2_fires, ico2_oce, ico2_bio use emission_common, only: fluxmultiplier,fluxmultiplier_biodis, fluxmultiplier_epsbio #ifdef with_budgets use budget_global, only : nbudg, nbud_vg, budg_dat, apply_budget_global, budget_hdf_names_unique #endif use io_hdf, only : DFACC_WRITE, io_write #ifdef MPI use mpi_const, only: ntracet_ar, lmloc, tracer_active, ntracetloc, myid, root, my_real, mpi_sum, com_trac, ierr, root_t #else use ParTools, only: lmloc, ntracet_ar,tracer_active, myid, root, ntracetloc, root_t #endif implicit none ! --- in/out -------------------------------- integer, intent(out) :: status character(len=*) :: budget_file ! --- const ------------------------------ character(len=*), parameter :: rname = mname//'/Emis_Done' ! local integer :: region, imr, jmr integer :: i,j,n, nzone,offsetn #ifdef with_budgets real, dimension(nbudg, nbud_vg, ntracet) :: budemig real, dimension(:,:,:,:),allocatable :: budget4d #endif integer :: ifail,sfend, nsend, io, sfstart #ifdef with_budgets if(apply_budget_global) then !handle budget emission: if(myid == root) then io = sfstart(budget_file, DFACC_WRITE) endif budemig(:,:,:) = 0.0 do region = 1, nregions allocate(budget4d(im(region),jm(region),nbud_vg,ntracet)) #ifdef MPI nsend = im(region)*jm(region)*nbud_vg*ntracet call mpi_allreduce(budemi_dat(region)%emi,budget4d, & nsend,my_real,mpi_sum,com_trac,ierr) #else budget4d = budemi_dat(region)%emi #endif call io_write(io,im(region),'im',jm(region),'jm', nbud_vg, 'nbud_vg', ntracet,'ntracet', & budget4d ,'budemi_dat') do n=1,ntracet do j=1,jm(region) do i=1,im(region) nzone = budg_dat(region)%nzong(i,j) budemig(nzone,:,n) = budemig(nzone,:,n) + & budget4d(i,j,:,n) end do end do!j end do ! deallocate(budget4d) enddo ! regions if(myid == root) then call io_write(io,nbudg,'nbud',nbud_vg,'nbud_v', & ntracet,'ntracet',budemig,'budemi') print *, 'Close budget file:', sfend(io) endif endif ! apply_budget_global do region = 1, nregions deallocate(budemi_dat(region)%emi) enddo #endif if(.not. flux1x1_override) then call calc_emission_averagers call save_emission_averagers(status) IF_NOTOK_RETURN(status=1) endif call free_emission_co2_ff call free_emission_co2_fires call free_emission_co2_ocean call free_emission_co2_bio call free_emission_co2_all !endif call free_emission_c13_bio call free_emission_c13_ocean call free_emission_c13_ff call free_emission_c13_fires call free_emission_c13_biodis call free_emission_c13_oceandis call free_emission_c13_all ! endif if(.not. flux1x1_override) then call free_emission_averagers endif deallocate(fluxmultiplier) deallocate(fluxmultiplier_biodis) deallocate(fluxmultiplier_epsbio) status=0 end subroutine Emis_Done end module emission