!### 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" ! !----------------------------------------------------------------------------- ! TM5 ! !----------------------------------------------------------------------------- !BOP ! ! !MODULE: REDGRIDZOOM ! ! !DESCRIPTION: !\\ !\\ ! !INTERFACE: ! MODULE REDGRIDZOOM ! ! !USES: ! use GO , only : gol, goErr, goPr use binas, only : pi use dims, only : nregions, im, jm, okdebug implicit none private ! ! !PUBLIC MEMBER FUNCTIONS: ! public :: uni2red_mf, redgrid_init, redgrid_done, calc_pdiff public :: uni2red, red2uni, red2uni_em #ifdef secmom public :: uni2red_2nd, red2uni_em_2nd #endif ! ! !PUBLIC DATA MEMBERS: ! public :: nred, nredmax, clustsize, grid_reduced, idle_xadv public :: jred, imredj public :: mfl_alt1, m_alt1, rm_alt1, rxm_alt1, rym_alt1, rzm_alt1 integer, parameter :: nredmax=60 logical :: grid_reduced=.true. logical :: idle_xadv=.false. ! proc is idle during x-advection (possible if reduced grid with npe_lon/=1) ! the number of reduced latitude circles per region on the current proc integer, dimension(nregions) :: nred=0 ! number of joined cells per latitude circle and region integer, allocatable :: clustsize(:,:) ! j1:j2,nregions ! reduced dimension...... integer, allocatable :: imredj(:,:) ! j1:j2,nregions ! latitaude numbers where the reduction applies integer, dimension(nredmax,nregions) :: jred ! alternate mass array (used on row_roots, to deal with reduced grid). With and w/o halo. real, dimension(:,:,:), target, allocatable :: mfl_alt1, m_alt1 real, dimension(:,:,:,:), target, allocatable :: rm_alt1, rxm_alt1, rym_alt1, rzm_alt1 real, dimension(:,:,:), allocatable :: mfl_alt2, m_alt2 real, dimension(:,:,:,:), allocatable :: rm_alt2 ! ! !PRIVATE DATA MEMBERS: ! character(len=*), parameter :: mname = 'redgridZoom' real, parameter :: dl = 2.0*pi/im(1) real, parameter :: dp = pi/jm(1) real, parameter :: epsx = epsilon(0.0) ! ! !REVISION HISTORY: ! 20 Dec 2012 - Ph Le Sager - TM5-MP version. ! ! !REMARKS: ! !EOP !------------------------------------------------------------------------ CONTAINS ! Fortran code for Split-rg (reduced grid) scheme for advection on the ! globe. ! ! VERSION: 1.1 ! DATE: November 12, 1998 ! ! Version 1.0 Dated October 19, 1998. ! Version 1.1 Corrected red2uni routine and added new uni2red_m ! routine. ! (Orography effects were erroneously not taken into ! account in the reduced to uniform grid conversion) ! ! This code solves the discrete advection equation for the tracer ! mass, with the wind field specified in air mass fluxes: ! ! dtracm/dt = ( ! mu chi (i-half) - mu chi (i+half) ! mv chi (j-half) - mv chi (j+half) ! mw chi (k-half) - mw chi (k+half) ) ! ! where ! ! tracm = tracer mass ! mu, mv, mw = air mass fluxes ! chi = tracm / mass = tracer mixing ratio ! mass = air mass ! ! The integer values of indices i, j, and k correspond to cell centers. ! ! See: ! ! Petersen, A.C., E.J. Spee, H. van Dop, and W. Hundsdorfer, ! An evaluation and intercomparison of four new advection schemes ! for use in global chemistry models, ! Journal of Geophysical Research, 103, 19253-19269, 1998. ! ! Written by Edwin Spee, CWI, Amsterdam, The Netherlands ! ! Parts of the code are written by ! Joke Blom and Willem Hundsdorfer (CWI) and ! ! The advection routine is called in the following way: ! ! call advect_split(tracm, mu, mv, mw, mass, dt) ! ! tracm = tracer mass real(im,jm,lm,ntracet) [kg] ! ! mu = longitudinal mass flux real(im,jm,lm) [kg/s] ! (defined on eastward side of grid ! cells; for all grid cells) ! ! mv = latitudinal mass flux real(im,jm,lm) [kg/s] ! (defined on southward side of grid ! cells; only j=2,..,jm are used) ! ! mw = vertical massflux real(im,jm,lm) [kg/s] ! (defined on upper side of grid cell; ! only l=1,..,lm-1 are used since it is ! assumed that there is no flux across ! the top of the model) ! ! The number of latitudinal elements actually used for mu, mv, and mw ! depends on the advection grid; for mv this number is the determined ! by the latitude row with the largest number of cells bordering a ! certain cell-bounding latitude circle. ! ! mass = air mass real(im,jm,lm) [kg] ! ! dt = timestep real [s] ! ! The physical units were chosen for compatibility with the TM3 model ! but can be changed for other global tracer models, provided that ! the calculation of the fluxes is done with an equivalent of tracer ! mixing ratios and not of tracer concentrations (make for other tracer ! models, if necessary, the appropriate changes in the routines ! split_lab, split_phi, and split_eta, where now tracer mixing ratios ! are calculated as tracm / mass; also the grid transformations as done ! in advect_split are different for tracer mass and concentration). ! ! The longitude is counted from west to east, the latitude from south ! to north, and the height from surface to top. ! ! The variables tracm, mu, mv, mw and mass should be declared in ! the calling program under any name; dt is a parameter. The ! parameters im, jm, and lm denote the number of grid cells in ! longitudinal, latitudinal, and vertical direction, respectively. ! ntracet is the number of transported tracers. The letters 'u', 'v', ! and 'w' for the three dimensions are in the code often replaced by ! 'l', 'p', and 'e', denoting the coordinates lambda, phi, and eta, ! respectively (see Petersen et al. 1998). The advection scheme works ! in principle with any grid and coordinate system of the main model. ! However, in the current implementation (made for the TM3 model) the ! variables tracm and mass are assumed to be defined in the main ! program on a uniform grid with polar cap, and the variables mu, mv, ! and mw on the advection grid (i.e., reduced grid with polar cap). ! Locally in the advection routines a transformation of tracm and mass, ! if necessary, is made to the advection grid. Please contact us for ! information on the changes that have to be made to work with other ! main model grids. Depending on the interests of users, we intend to ! generalize the code with respect to the main model grid in the next ! version. ! ! The following parameters and variables, related to the specification ! of the advection grid are used by the routines of the advection ! scheme: ! ! integer, parameter :: nredmax = 10 ! real, parameter :: dl = 2.0 * pi / im ! real, parameter :: dp = pi / jm ! integer nred ! integer clustsize(0:nredmax), ! jred(-nredmax:nredmax+1),imredj(jm) ! real maxtau ! integer advsub ! ! nredmax = maximum number of grid reductions per hemisphere ! (in this implementation the minimum number of grid ! reductions per hemisphere equals one, since use is ! made of polar caps) ! dl = longitudinal uniform grid cell size in lambda ! coordinates ! dp = latitudinal grid cell size in phi coordinates ! nred = number of grid reductions ! clustsize = array with the ratio im / imred, as a function of the ! absolute value of the latitude zone index, where imred ! is the actually used number of cells in the longitudinal ! direction for the given latitude zone ! jred = latitude number where grid reduction starts (seen from ! SP), as a function of the reduction zone index ! imredj(j) = number of cells for latitude #j ! maxtau = time step such that CFL <= 1 in all directions (given ! the directional splitting, see routine advect_split) ! advsub = largest integer such that dtadv / advsub <= maxtau, ! with dtadv the advection time step in operator splitting ! ! The advection grid variables for the reduced grid are automatically ! set by the following line: ! ! call redgrid_init(clustsize, jred, imredj, nred, nredmax, im, jm, ! . pi, .false.) ! ! The variables maxtau and advsub should be set each time new mass ! fluxes (wind fields) are available in the model through ! ! call split_maxtau(mu, mv, mw, mass, dtadv) ! ! where dtadv is the advection time step in the main model. ! ! Furthermore, some routines in the advection scheme need information ! on the machine precision. The following variable is used for this ! purpose: ! ! real, parameter :: epsx = epsilon(0.0) ! ! The value of epsx is dependent on the computing platform that is used. ! epsilon is a Fortran 90 intrinsic function. For example, for our Cray ! C90 eps is about 1e-14. ! ! The above-mentioned parameters and variables are assumed to be ! declared in the module global_data (Fortran 90) or in some common ! block (Fortran 77). The current implementation (Fortran 90) uses the ! statements ! ! use global_data ! implicit none ! ! If Fortran 77 is preferred one should define a common block or ! use an existing common block to include the above-given parameters ! and variables, e.g. ! ! implicit none ! (declarations as given in the above) ! common/advgrid/nredmax,dl,dp,nred,clustsize,jred,maxtau,advsub ! ! Also the parameters im, jm, lm, ntracet, pi, twopi (= 2 * pi), and ! eps are assumed to be declared in the module global_data. ! ! -------------------------- ! Example of a reduced grid ! -------------------------- ! ! NP=North Pole ! SP=South Pole ! EQ=Equator ! | grid boundary in longitudal direction ! o cell center ! j imredj clustsize latitude zone ! NP| o | 13 1 12 3 ! | o | o | o | 12 3 4 2 ! | o | o | o | o | o | o | 11 6 2 1 ! | o | o | o | o | o | o | 10 6 2 1 ! |o|o|o|o|o|o|o|o|o|o|o|o| 9 12 1 0 ! |o|o|o|o|o|o|o|o|o|o|o|o| 8 12 1 0 ! EQ|o|o|o|o|o|o|o|o|o|o|o|o| 7 12 1 0 ! |o|o|o|o|o|o|o|o|o|o|o|o| 6 12 1 0 ! |o|o|o|o|o|o|o|o|o|o|o|o| 5 12 1 0 ! | o | o | o | o | o | o | 4 6 2 -1 ! | o | o | o | o | o | o | 3 6 2 -1 ! | o | o | o | 2 3 4 -2 ! SP| o | 1 1 12 -3 ! ! This grid would have ! im = 12, jm = 13, nred = 3 ! ! jred(-3) = 1 by definition ! jred(-2) = 2 ! jred(-1) = 3 ! jred( 0) = 5 ! jred( 1) =10 ! jred( 2) =12 ! jred( 3) =13 ! jred( 4) =14=jm+1 by definition ! ! clustsize(3) = 12 ! clustsize(2) = 4 ! clustsize(1) = 2 ! clustsize(0) = 1 !CMK total new definition to allow zoom regios.. ! nred(nregions) :: number of reduced latitudes circles < jm(region) ! nredmax :: maximum number of reduced latitude circles... ! clustsize(nredmax.nregions) :: number of joined cells for n'th circle ! jred(nredmax.nregions) :: corresponding j value in 1...jm(region) array ! imredj(nredmax.nregions) :: 'new' im value for the reduced grid... !CMK !PLS new definitions for TM5-MP ! nred(nregions) :: number of reduced latitudes circles on current tile (ie proc) and region ! nredmax :: maximum number of reduced latitude circles ! clustsize(j1:j2, nregions) :: number of joined cells for J'th circle (default=1=no_cluster=no_reduced_grid) ! jred (nredmax, nregions) :: j value for 1...nred(region) in 1..jmr ! imredj (j1:j2, nregions) :: 'new' im value for the reduced grid, default to im(region) !PLS !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: uni2red_mf ! ! !DESCRIPTION: Converts mass fluxes mfl from uniform grid to reduced grid !\\ !\\ ! !INTERFACE: ! subroutine uni2red_mf(region,i1,i2,j1,j2, status) ! ! !USES: ! use dims, only : lm, im use global_data, only : wind_dat use partools, only : isRowRoot, npe_lon use tm5_distgrid, only : dgrid, gather_j_band ! ! !INPUT PARAMETERS: ! integer,intent(in) :: region,i1,i2,j1,j2 ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !REVISION HISTORY: ! Written by Edwin Spee, CWI, Amsterdam, The Netherlands ! cmk changed..... ! 20 Dec 2012 - Ph Le Sager - TM5-MP version ! ! !REMARKS: ! - no need to check on tile bounds, zero-trip loop for jreg does the job !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/uni2red_mf' real,dimension(:,:,:),pointer :: mfl, mfl_alt integer :: clust_size, i, j, l, jreg integer :: lmr, imr, halo ! -- start mfl => wind_dat(region)%am halo=wind_dat(region)%halo imr=im(region) lmr=lm(region) ! Two scenarios if (npe_lon /= 1) then ! ! decomposition along longitudes => work on row_root ! if (nred(region)/=0) then ! some reduce grid on this tile ! array to gather (remove halo) allocate(mfl_alt( i1:i2, j1:j2, 0:lmr+1)) mfl_alt = mfl(i1:i2,j1:j2,:) CALL GATHER_J_BAND(dgrid(region), mfl_alt, mfl_alt2, status, jref=j1, rowcom=.true.) IF_NOTOK_RETURN(status=1) if (isRowRoot) then mfl_alt1(1:imr,:,:) = mfl_alt2 mfl_alt1( 0,:,:) = mfl_alt1(imr,:,:) mfl_alt1( -1,:,:) = mfl_alt1(imr-1,:,:) mfl_alt1(imr+1,:,:) = mfl_alt1(1,:,:) !PLS mfl_alt1(imr+2,:,:) = mfl_alt1(2,:,:) ! full on do l=1,lmr do jreg = 1,nred(region) ! loop over the reduced latitudes j = jred(jreg,region) ! latitude clust_size = clustsize(j,region) ! clustersize do i=1,imredj(j,region) mfl_alt1(i,j,l) = mfl_alt1(i*clust_size,j,l) ! compress mfl end do do i=imredj(j,region)+1,im(region) mfl_alt1(i,j,l) = -1.0 end do ! Update boundary condition. Only meaningful for x-cyclic regions (which is the ! only possibility in TM5-MP). mfl_alt1( 0,j,l)=mfl_alt1(imredj(j,region),j,l) mfl_alt1(-1,j,l)=mfl_alt1(imredj(j,region)-1,j,l) ! needed to limit mpi comm in advectx__slopes.F90 mfl_alt1(imredj(j,region)+1,j,l)=mfl_alt1(1,j,l) ! needed to limit mpi comm in advectx__slopes.F90 end do end do endif deallocate( mfl_alt ) endif else ! ! Reduced grid without longitudinal decomposition ! do l=1,lmr do jreg = 1,nred(region) ! loop over the reduced latitudes j = jred(jreg,region) ! latitude clust_size = clustsize(j,region) ! clustersize do i=1,imredj(j,region) mfl(i,j,l) = mfl(i*clust_size,j,l) ! compress mfl end do do i=imredj(j,region)+1,im(region) mfl(i,j,l) = -1.0 end do ! Update boundary condition. Only meaningful for x-cyclic regions (which is the ! only possibility in TM5-MP). mfl( 0,j,l)=mfl(imredj(j,region),j,l) mfl(-1,j,l)=mfl(imredj(j,region)-1,j,l) ! needed to limit mpi comm in advectx__slopes.F90 mfl(imredj(j,region)+1,j,l)=mfl(1,j,l) ! needed to limit mpi comm in advectx__slopes.F90 end do end do endif nullify(mfl) status=0 end subroutine uni2red_mf !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: REDGRID_INIT ! ! !DESCRIPTION: allocate work arrays, read settings from rc file and fill ! reduced grid parameters (clustsize, jred, and imredj) for ! current region and tile. !\\ !\\ ! !INTERFACE: ! SUBROUTINE REDGRID_INIT( region, status ) ! ! !USES: ! use GO, only : TrcFile, Init, Done, ReadRc use dims, only : lm, im use global_data, only : wind_dat use tracer_data, only : mass_dat use meteodata , only : m_dat use chem_param, only : ntracet use tm5_distgrid, only : dgrid, Get_DistGrid use partools, only : isRowRoot, npe_lon use global_data, only : rcfile ! ! !INPUT PARAMETERS: ! integer, intent(in) :: region ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !REVISION HISTORY: ! Written by Edwin Spee, CWI, Amsterdam, The Netherlands, based on work of Joke Blom. ! 20 Dec 2012 - Ph Le Sager - TM5-MP ! ! !REMARKS: ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/redgrid_init' integer :: j, i1, i2, j1, j2, imr, lmr, halo character(len=256) :: fname, line type(TrcFile) :: rcF ! --- begin ----------------------------- ! Open rcfile call Init( rcF, rcfile, status ) IF_NOTOK_RETURN(status=1) ! grid reduction ? call ReadRc( rcF, 'proces.advection.reduced', grid_reduced, status ) IF_NOTOK_RETURN(status=1) if (grid_reduced) then ! Arrays CALL GET_DISTGRID( dgrid(region), I_STRT=i1, I_STOP=i2, J_STRT=j1, J_STOP=j2 ) allocate( clustsize(j1:j2, nregions)) allocate( imredj(j1:j2, nregions)) call Read_RedGrid_rc( region, status ) IF_NOTOK_RETURN(status=1) ! Extra arrays to work on row_root (only if really needed) if (nred(region)/=0 .and. (npe_lon/=1)) then if (isRowRoot) then imr = im(region) lmr = lm(region) halo= wind_dat(1)%halo allocate(mfl_alt1( 1-halo:imr+halo, j1:j2, 0:lmr+1)) allocate(mfl_alt2( imr, j1:j2, 0:lmr+1)) halo= m_dat(1)%halo allocate( m_alt1( 1-halo:imr+halo, j1:j2, lmr)) allocate( m_alt2( imr, j1:j2, lmr)) halo= mass_dat(1)%halo allocate(rm_alt1( 1-halo:imr+halo, j1:j2, lmr, ntracet)) allocate(rm_alt2( imr, j1:j2, lmr, ntracet)) #ifdef slopes allocate(rxm_alt1( 1-halo:imr+halo, j1:j2, lmr, ntracet)) allocate(rym_alt1( 1-halo:imr+halo, j1:j2, lmr, ntracet)) allocate(rzm_alt1( 1-halo:imr+halo, j1:j2, lmr, ntracet)) #endif else idle_xadv=.true. allocate(mfl_alt2(1,1,1) ) allocate( m_alt2(1,1,1) ) allocate( rm_alt2(1,1,1,1)) endif endif endif ! close call Done( rcF, status ) IF_NOTOK_RETURN(status=1) status = 0 END SUBROUTINE REDGRID_INIT !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: REDGRID_DONE ! ! !DESCRIPTION: deallocate work arrays !\\ !\\ ! !INTERFACE: ! SUBROUTINE REDGRID_DONE( status ) ! ! !USES: ! use partools, only : isRowRoot, npe_lon ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !REMARKS: used only for region=1, in a hardcoded way ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/redgrid_done' ! --- begin ----------------------------- if (grid_reduced) then deallocate( clustsize) deallocate( imredj) if (nred(1)/=0 .and. (npe_lon/=1)) then if (isRowRoot) then deallocate(mfl_alt1) deallocate(mfl_alt2) deallocate( m_alt1) deallocate( m_alt2) deallocate( rm_alt1) deallocate( rm_alt2) #ifdef slopes deallocate(rxm_alt1) deallocate(rym_alt1) deallocate(rzm_alt1) #endif else deallocate(mfl_alt2) deallocate( m_alt2) deallocate( rm_alt2) endif endif endif status = 0 END SUBROUTINE REDGRID_DONE !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: READ_REDGRID_RC ! ! !DESCRIPTION: !\\ !\\ ! !INTERFACE: ! subroutine Read_RedGrid_rc( region, status ) ! ! !USES: ! use GO, only : TrcFile, Init, Done, ReadRc use global_data, only : rcfile use dims, only : region_name, ybeg, yend, dy, yref use tm5_distgrid, only : dgrid, Get_DistGrid use partools, only : isRoot ! ! !INPUT PARAMETERS: ! integer, intent(in) :: region ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !REVISION HISTORY: ! 20 Dec 2012 - Ph Le Sager - TM5-MP version ! 6 Jun 2013 - P. Le Sager - check for consecutive 1-cell rings ! ! !REMARKS: ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/Read_RedGrid_rc' type(TrcFile) :: rcF integer :: imr, jmr, xam, nim integer :: i,j, i1, i2, j1, j2 integer :: jband integer :: nred_nh, nred_sh integer, allocatable :: ncombs_nh(:), ncombs_sh(:) real :: y0, y1 ! --- begin ---------------------------------------- ! number of lats imr = im(region) jmr = jm(region) ! open settings: call Init( rcF, rcfile, status ) IF_NOTOK_RETURN(status=1) ! number of reduced lat bands: call ReadRc( rcF, 'region.'//trim(region_name(region))//'.redgrid.nh.n', nred_nh, status ) IF_NOTOK_RETURN(status=1) ! number of reduced lat bands: call ReadRc( rcF, 'region.'//trim(region_name(region))//'.redgrid.sh.n', nred_sh, status ) IF_NOTOK_RETURN(status=1) ! total number of reduced bands nred(region) = nred_nh + nred_sh ! local bound CALL GET_DISTGRID( dgrid(region), I_STRT=i1, I_STOP=i2, J_STRT=j1, J_STOP=j2 ) ! reduced grid ? if ( nred(region) > 0 ) then if ( nred(region) > nredmax ) then write (gol,'("problem with reduced grid. nred > nredmax")'); call goErr TRACEBACK; status=1; return end if ! test done, reset nred(region)=0 ! to read number of combined cells per lat band allocate( ncombs_sh(jmr) ) ncombs_sh=1 allocate( ncombs_nh(jmr) ) ! southern hemisphere if ( nred_sh > 0 ) then call ReadRc( rcF, 'region.'//trim(region_name(region))//'.redgrid.sh.comb', ncombs_sh(1:nred_sh), status ) IF_NOTOK_RETURN(status=1) end if ! northern hemisphere if ( nred_nh > 0 ) then call ReadRc( rcF, 'region.'//trim(region_name(region))//'.redgrid.nh.comb', ncombs_nh(1:nred_nh), status ) IF_NOTOK_RETURN(status=1) do jband= 1, nred_nh ncombs_sh(jmr-jband+1) = ncombs_nh(jband) end do end if ! fill work arrays clustsize(:,region) = ncombs_sh(j1:j2) nred(region) = count(clustsize(:,region)/=1) jred(1:nred(region), region) = pack((/(j,j=j1,j2)/), clustsize(:,region)/=1) imredj(:,region) = imr ! default do i=1,nred(region) j = jred(i, region) imredj(j,region) = imr/clustsize(j,region) end do ! Testing (global) do j = 1, jmr ! check if number of combined cells matches with grid: if ( modulo(imr,ncombs_sh(j)) /= 0 ) then write (gol,'("number of combined cells not ok:")') ; call goErr write (gol,'(" region : ",i2," ",a)') region, trim(region_name(region)) ; call goErr write (gol,'(" lat band : ",i4)') j ; call goErr write (gol,'(" combined cells : ",i4)') ncombs_sh(j) ; call goErr write (gol,'(" im : ",i4)') imr ; call goErr TRACEBACK; status=1; return end if ! check with previous ... if ( j > 1 ) then xam=max(ncombs_sh(j-1),ncombs_sh(j)) nim=min(ncombs_sh(j-1),ncombs_sh(j)) if ( modulo(xam,nim) /= 0 ) then write (gol,'("number of combined cells does match with previous:")') ; call goErr write (gol,'(" region : ",i2," ",a)') region, trim(region_name(region)) ; call goErr write (gol,'(" lat band : ",i4)') j ; call goErr write (gol,'(" combined cells : ",i4)') ncombs_sh(j) ; call goErr write (gol,'(" previous : ",i4)') ncombs_sh(j-1) ; call goErr TRACEBACK; status=1; return end if end if ! check for consecutive rings with a single reduced grid cell if ( j > 1 ) then if ( ncombs_sh(j-1)==imr .and. ncombs_sh(j)==imr ) then write (gol,'("Consecutive rings with a single reduced grid cell detected:")') ; call goErr write (gol,'(" region : ",i2," ",a)') region, trim(region_name(region)) ; call goErr write (gol,'(" lat band : ",i4)') j ; call goErr write (gol,'(" combined cells : ",i4)') ncombs_sh(j) ; call goErr write (gol,'(" previous : ",i4)') ncombs_sh(j-1) ; call goErr TRACEBACK; status=1; return end if end if end do ! Verbose if (isRoot) then write (gol,'(" [Reduced grid] Region """,a,""":")') trim(region_name(region)); call goPr write (gol,'(" [Reduced grid] - Uniform grid: no. of cells in each latitude band: ",i3,".")') imr; call goPr do j=1,jmr if ( ncombs_sh(j) /= 1 ) then y0=ybeg(region)+(dy/yref(region))*(j-1) y1=ybeg(region)+(dy/yref(region))*j write (gol,'(" [Reduced grid] - from ",f6.2," to ",f6.2," degrees latitude: ",i3," cells.")') & y0,y1,ncombs_sh(j); call goPr end if end do end if deallocate( ncombs_sh ) deallocate( ncombs_nh ) else write (gol,'("No reduced grid for region :",a)') trim(region_name(region)); call goPr nred(region) = 0 end if ! Done call Done( rcF, status ) IF_NOTOK_RETURN(status=1) status = 0 END SUBROUTINE READ_REDGRID_RC !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: CALC_PDIFF ! ! !DESCRIPTION: calculate max pressure difference, accounting for reduce ! grid if used. !\\ !\\ ! !INTERFACE: ! SUBROUTINE CALC_PDIFF( region, pdiffmax ) ! ! !USES: ! use meteodata, only : sp_dat, sp1_dat use tm5_distgrid, only : dgrid, Get_DistGrid use partools, only : isRowRoot, npe_lon ! ! !INPUT PARAMETERS: ! integer, intent(in) :: region ! ! !OUTPUT PARAMETERS: ! real, intent(out) :: pdiffmax ! max pressure difference ! ! !REVISION HISTORY: ! 20 Dec 2012 - Ph Le Sager - TM5-MP version ! ! !REMARKS: ! !EOP !------------------------------------------------------------------------ !BOC integer :: lrg, i, j, ratio, idx, iu, i1, i2, j1, j2 real :: work, work1 real, pointer :: p(:,:,:), pold(:,:,:) !---- start p => sp_dat(region)%data pold => sp1_dat(region)%data pdiffmax = 0.0 ! local bound CALL GET_DISTGRID( dgrid(region), I_STRT=i1, I_STOP=i2, J_STRT=j1, J_STOP=j2 ) ! if no reduced grid if ( .not. grid_reduced ) then pdiffmax = maxval(abs(p(i1:i2,j1:j2,1)-pold(i1:i2,j1:j2,1))) else if (npe_lon==1) then ! first maximum over the reduced grid do lrg = 1,nred(region) j = jred(lrg,region) ratio = clustsize(j,region) do i = 1,imredj(j,region) work = 0.0 work1 = 0.0 do iu = 1,ratio idx = (i-1)*ratio+iu work = work + p(idx,j,1) work1 = work1 + pold(idx,j,1) end do pdiffmax=max(pdiffmax,abs(work-work1)/ratio) end do end do ! now the pressure difference in the non-reduced part j2loop: do j=j1,j2 if(imredj(j,region).ne.im(region)) cycle j2loop !if reduced...skip i2loop: do i=i1,i2 ! should be same as 1,im(region) by design pdiffmax = max(pdiffmax,abs(p(i,j,1)-pold(i,j,1))) end do i2loop end do j2loop else if (nred(region)/=0) then jloop:do j=j1,j2 if(imredj(j,region).ne.im(region)) then ! this requires gathering of sp_dat and sp1_dat across mpi tasks - just skip it (FIXME?) cycle jloop else pdiffmax = max( pdiffmax, maxval(abs(p(i1:i2,j,1)-pold(i1:i2,j,1))) ) endif enddo jloop else pdiffmax = maxval(abs(p(i1:i2,j1:j2,1)-pold(i1:i2,j1:j2,1))) endif endif end if end subroutine calc_pdiff !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: UNI2RED ! ! !DESCRIPTION: transforms data (air mass & tracers) from uniform grid to reduced grid !\\ !\\ ! !INTERFACE: ! subroutine uni2red( region, i1, i2, j1, j2, status) ! ! !USES: ! use dims, only : im,jm,lm use meteodata , only : m_dat use tracer_data, only : mass_dat use chem_param, only : ntracet use partools, only : isRowRoot, npe_lon use tm5_distgrid, only : dgrid, gather_j_band ! ! !INPUT PARAMETERS: ! integer, intent(in) :: region,i1,i2,j1,j2 ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !REVISION HISTORY: ! written by mike botchev, march-june 1999 ! modified by Maarten Krol, dec 2002 ! 20 Dec 2012 - Ph Le Sager - TM5-MP ! ! !REMARKS: ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'_MOD_uni2red' real,dimension(:,:,:,:),pointer :: rm, rm_alt #ifdef slopes real,dimension(:,:,:,:),pointer :: rxm,rym,rzm #endif real,dimension(:,:,:), pointer :: m, m_alt integer :: i,ie,is,j,l,lrg,redfact,n, imr, lmr, halo real :: summ,sumrm ! -------------- start m => m_dat(region)%data rm => mass_dat(region)%rm #ifdef slopes rxm => mass_dat(region)%rxm rym => mass_dat(region)%rym rzm => mass_dat(region)%rzm #endif imr=im(region) lmr=lm(region) ! Two scenarios if (npe_lon /= 1) then ! ! decomposition along longitudes => work on row_root ! if (nred(region)/=0) then ! some reduce grid on this tile halo=m_dat(region)%halo ! local array to gather (remove i-halo) allocate( m_alt( i1:i2, j1:j2, lmr)) allocate(rm_alt( i1:i2, j1:j2, lmr, ntracet)) m_alt = m(i1:i2,j1:j2,:) CALL GATHER_J_BAND(dgrid(region), m_alt, m_alt2, status, jref=j1, rowcom=.true.) IF_NOTOK_RETURN(status=1) if (isRowRoot) then m_alt1( 1:imr,:,:) = m_alt2 m_alt1( 0,:,:) = m_alt1( imr,:,:) m_alt1( -1,:,:) = m_alt1(imr-1,:,:) !PLS m_alt1( imr+1,:,:) = m_alt1(1,:,:) !PLS m_alt1( imr+2,:,:) = m_alt1(2,:,:) endif rm_alt = rm(i1:i2,j1:j2,:,:) CALL GATHER_J_BAND(dgrid(region), rm_alt, rm_alt2, status, jref=j1, rowcom=.true.) IF_NOTOK_RETURN(status=1) if (isRowRoot) then rm_alt1( 1:imr,:,:,:) = rm_alt2 rm_alt1( 0,:,:,:) = rm_alt1( imr,:,:,:) rm_alt1( -1,:,:,:) = rm_alt1(imr-1,:,:,:) !PLS rm_alt1( imr+1,:,:,:) = rm_alt1( 1,:,:,:) !PLS rm_alt1( imr+2,:,:,:) = rm_alt1( 2,:,:,:) endif #ifdef slopes rm_alt = rxm(i1:i2,j1:j2,:,:) CALL GATHER_J_BAND(dgrid(region), rm_alt, rm_alt2, status, jref=j1, rowcom=.true.) IF_NOTOK_RETURN(status=1) if (isRowRoot) then rxm_alt1( 1:imr,:,:,:) = rm_alt2 rxm_alt1( 0,:,:,:) = rm_alt1( imr,:,:,:) rxm_alt1( -1,:,:,:) = rm_alt1(imr-1,:,:,:) !PLS rxm_alt1( imr+1,:,:,:) = rm_alt1( 1,:,:,:) !PLS rxm_alt1( imr+2,:,:,:) = rm_alt1( 2,:,:,:) endif rm_alt = rym(i1:i2,j1:j2,:,:) CALL GATHER_J_BAND(dgrid(region), rm_alt, rm_alt2, status, jref=j1, rowcom=.true.) IF_NOTOK_RETURN(status=1) if (isRowRoot) then rym_alt1( 1:imr,:,:,:) = rm_alt2 rym_alt1( 0,:,:,:) = rm_alt1( imr,:,:,:) rym_alt1( -1,:,:,:) = rm_alt1(imr-1,:,:,:) !PLS rym_alt1( imr+1,:,:,:) = rm_alt1( 1,:,:,:) !PLS rym_alt1( imr+2,:,:,:) = rm_alt1( 2,:,:,:) endif rm_alt = rzm(i1:i2,j1:j2,:,:) CALL GATHER_J_BAND(dgrid(region), rm_alt, rm_alt2, status, jref=j1, rowcom=.true.) IF_NOTOK_RETURN(status=1) if (isRowRoot) then rzm_alt1( 1:imr,:,:,:) = rm_alt2 rzm_alt1( 0,:,:,:) = rm_alt1( imr,:,:,:) rzm_alt1( -1,:,:,:) = rm_alt1(imr-1,:,:,:) !PLS rzm_alt1( imr+1,:,:,:) = rm_alt1( 1,:,:,:) !PLS rzm_alt1( imr+2,:,:,:) = rm_alt1( 2,:,:,:) endif #endif if (isRowRoot) then do lrg=1,nred(region) j = jred(lrg,region) redfact=clustsize(j,region) do l=1,lmr ! air mass do i = 1,imredj(j,region) ! the is:ie array section will be reduced to i is = (i-1)*redfact + 1 ie = i*redfact summ = sum(m_alt1(is:ie,j,l)) m_alt1(i,j,l) = summ end do ! tracers do n=1, ntracet do i = 1,imredj(j,region) is = (i-1)*redfact + 1 ie = i*redfact sumrm = sum(rm_alt1(is:ie,j,l,n)) rm_alt1(i,j,l,n) = sumrm #ifdef slopes sumrm = sum(rxm_alt1(is:ie,j,l,n)) rxm_alt1(i,j,l,n) = sumrm sumrm = sum(rym_alt1(is:ie,j,l,n)) rym_alt1(i,j,l,n) = sumrm sumrm = sum(rzm_alt1(is:ie,j,l,n)) rzm_alt1(i,j,l,n) = sumrm #endif end do ! JFM: set remaining masses to zero ! for consistency with adjoint do i = imredj(j,region)+1, im(region) rm_alt1(i,j,l,n) = 0. #ifdef slopes rxm_alt1(i,j,l,n) = 0. rym_alt1(i,j,l,n) = 0. rzm_alt1(i,j,l,n) = 0. #endif end do end do end do !l enddo endif deallocate(m_alt, rm_alt) endif else ! ! Reduced grid without longitudinal decomposition ! do lrg=1,nred(region) j = jred(lrg,region) redfact=clustsize(j,region) do l=1,lmr ! air mass do i = 1,imredj(j,region) ! the is:ie array section will be reduced to i is = (i-1)*redfact + 1 ie = i*redfact summ = sum(m(is:ie,j,l)) m(i,j,l) = summ end do !cmkm_uni(is:ie,j,l) = m_uni(is:ie,j,l)/summ ! use as distribution function ! when transferring back from reduced--->uniform grid... ! these summations mean that mixing ratio and the ! the slopes are averaged out within the is:ie section ! with m(is:ie,...) taken as the weights ! tracers do n=1, ntracet do i = 1,imredj(j,region) is = (i-1)*redfact + 1 ie = i*redfact sumrm = sum(rm(is:ie,j,l,n)) rm(i,j,l,n) = sumrm #ifdef slopes sumrm = sum(rxm(is:ie,j,l,n)) rxm(i,j,l,n) = sumrm sumrm = sum(rym(is:ie,j,l,n)) rym(i,j,l,n) = sumrm sumrm = sum(rzm(is:ie,j,l,n)) rzm(i,j,l,n) = sumrm #endif end do ! JFM: set remaining masses to zero ! for consistency with adjoint do i = imredj(j,region)+1, im(region) rm(i,j,l,n) = 0. #ifdef slopes rxm(i,j,l,n) = 0. rym(i,j,l,n) = 0. rzm(i,j,l,n) = 0. #endif end do end do !put periodic boundary... end do !l end do !redgrid... endif nullify(m) nullify(rm) #ifdef slopes nullify(rxm) nullify(rym) nullify(rzm) #endif status=0 end subroutine uni2red !EOC #ifdef secmom !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: uni2red_2nd ! ! !DESCRIPTION: Transforms data from uniform grid to reduced grid. @nd ! moments version. !\\ !\\ ! !INTERFACE: ! subroutine uni2red_2nd(region) ! ! !USES: ! use dims use meteodata, only : m_dat use global_data, only : mass_dat use chem_param, only : ntracet ! ! !INPUT PARAMETERS: ! integer,intent(in) :: region ! ! !REVISION HISTORY: ! written by mike botchev, march-june 1999 ! modified by Maarten Krol, dec 2002 ! 20 Dec 2012 - Ph Le Sager - ! ! !REMARKS: ! !EOP !------------------------------------------------------------------------ !BOC real, dimension(:,:,:,:), pointer :: rm, rxm,rym,rzm, rxxm, rxym, rxzm, ryym, ryzm, rzzm real, dimension(:,:,:), pointer :: m integer i,ie,is,j,l,lrg,redfact,n,lmr real summ,sumrm ! start m => m_dat(region)%data rm => mass_dat(region)%rm rxm => mass_dat(region)%rxm rym => mass_dat(region)%rym rzm => mass_dat(region)%rzm rxxm => mass_dat(region)%rxxm rxym => mass_dat(region)%rxym rxzm => mass_dat(region)%rxzm ryym => mass_dat(region)%ryym ryzm => mass_dat(region)%ryzm rzzm => mass_dat(region)%rzzm lmr=lm(region) do lrg=1,nred(region) j = jred(lrg,region) redfact=clustsize(j,region) do l=1,lmr do i = 1,imredj(j,region) ! the is:ie array section will be reduced to i is = (i-1)*redfact + 1 ie = i*redfact summ = sum(m(is:ie,j,l)) m(i,j,l) = summ !cmkm_uni(is:ie,j,l) = m_uni(is:ie,j,l)/summ ! use as distribution function ! when transferring back from reduced--->uniform grid... ! these summations mean that mixing ratio and the ! the slopes are averaged out within the is:ie section ! with m(is:ie,...) taken as the weights do n=1,ntracet !#endif sumrm = sum(rm(is:ie,j,l,n)) rm(i,j,l,n) = sumrm sumrm = sum(rxm(is:ie,j,l,n)) rxm(i,j,l,n) = sumrm sumrm = sum(rym(is:ie,j,l,n)) rym(i,j,l,n) = sumrm sumrm = sum(rzm(is:ie,j,l,n)) rzm(i,j,l,n) = sumrm sumrm = sum(rxxm(is:ie,j,l,n)) rxxm(i,j,l,n) = sumrm sumrm = sum(rxym(is:ie,j,l,n)) rxym(i,j,l,n) = sumrm sumrm = sum(rxzm(is:ie,j,l,n)) rxzm(i,j,l,n) = sumrm sumrm = sum(ryym(is:ie,j,l,n)) ryym(i,j,l,n) = sumrm sumrm = sum(ryzm(is:ie,j,l,n)) ryzm(i,j,l,n) = sumrm sumrm = sum(rzzm(is:ie,j,l,n)) rzzm(i,j,l,n) = sumrm end do !n end do !i !put periodic boundary... end do !l end do !redgrid... nullify(m) nullify(rm) nullify(rxm) nullify(rym) nullify(rzm) nullify(rxxm) nullify(rxym) nullify(rxzm) nullify(ryym) nullify(ryzm) nullify(rzzm) end subroutine uni2red_2nd !EOC #endif !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: RED2UNI ! ! !DESCRIPTION: transforms data from reduced grid back to uniform grid !\\ !\\ ! !INTERFACE: ! subroutine red2uni(region, i1, i2, j1, j2, halo, m_uni, status) ! ! !USES: ! use dims, only : im, lm use tracer_data, only : mass_dat use meteodata, only : m_dat use chem_param, only : ntracet use partools, only : isRowRoot, npe_lon use tm5_distgrid, only : dgrid, scatter_j_band, gather_j_band, update_halo_jband ! ! !INPUT PARAMETERS: ! integer, intent(in) :: region, i1, i2, j1, j2, halo real, intent(in) :: m_uni(i1-halo:,j1-halo:,1:) ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !REVISION HISTORY: ! written by mike botchev, march-june 1999 ! 20 Dec 2012 - Ph Le Sager - TM5-MP version ! ! !REMARKS: ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'_MOD_red2uni' ! local real,dimension(:,:,:,:),pointer :: rm, rm_ #ifdef slopes real,dimension(:,:,:,:),pointer :: rxm,rym,rzm #endif real,dimension(:,:,:) ,pointer :: m, m_, m_alt integer :: i, ie, is, j, l, lrg, n, redfact, imr, lmr, nt, halo_r real :: hi, mass, mass_coord, rmm, slope, m_old character(len=5) :: distr_mode ! start nt=ntracet imr=im(region) lmr=lm(region) m => m_dat(region)%data rm => mass_dat(region)%rm #ifdef slopes rxm => mass_dat(region)%rxm rym => mass_dat(region)%rym rzm => mass_dat(region)%rzm #endif ! Two scenarios if (npe_lon /= 1) then ! ! decomposition along longitudes => work on row_root ! if (nred(region)/=0) then ! gather m_uni on row_root, where we redistribute the masses allocate(m_(i1:i2,j1:j2,1:lmr)) m_ = m_uni(i1:i2,j1:j2,1:lmr) CALL GATHER_J_BAND(dgrid(region), m_, m_alt2, status, jref=j1, rowcom=.true.) IF_NOTOK_RETURN(status=1) if (isRowRoot) then do lrg=1,nred(region) j = jred(lrg,region) redfact=clustsize(j,region) !m(1:imr,j,:)= m_alt2(:,j,:) do l=1,lmr do i=imredj(j,region),1,-1 is = (i-1)*redfact + 1 ie = i*redfact mass=m_alt1(i,j,l) do n=1,nt rmm = rm_alt1(i,j,l,n) rm_alt1(is:ie,j,l,n)= m_alt2(is:ie,j,l)/mass* rmm #ifdef slopes rmm = rxm_alt1(i,j,l,n) rxm_alt1(is:ie,j,l,n)= m_alt2(is:ie,j,l)/mass * rmm ! rym and rzm are always distributed uniformly: rmm = rym_alt1(i,j,l,n) rym_alt1(is:ie,j,l,n)= m_alt2(is:ie,j,l)/mass * rmm rmm = rzm_alt1(i,j,l,n) rzm_alt1(is:ie,j,l,n)= m_alt2(is:ie,j,l)/mass * rmm #endif end do end do enddo enddo endif ! scatter results if (nred(region)< (j2-j1+1)) then if (isRowRoot) m_alt2 = m_alt1(1:imr,:,:) CALL SCATTER_J_BAND(dgrid(region), m_, m_alt2, status, jref=j1, rowcom=.true.) IF_NOTOK_RETURN(status=1) m(i1:i2,j1:j2,:) = m_ endif allocate(rm_( i1:i2, j1:j2, lmr, ntracet)) if (isRowRoot) rm_alt2 = rm_alt1(1:imr,:,:,:) CALL SCATTER_J_BAND(dgrid(region), rm_, rm_alt2, status, jref=j1, rowcom=.true.) IF_NOTOK_RETURN(status=1) rm(i1:i2,j1:j2,:,:) = rm_ #ifdef slopes if (isRowRoot) rm_alt2 = rxm_alt1(1:imr,:,:,:) CALL SCATTER_J_BAND(dgrid(region), rm_, rm_alt2, status, jref=j1, rowcom=.true.) IF_NOTOK_RETURN(status=1) rxm(i1:i2,j1:j2,:,:) = rm_ if (isRowRoot) rm_alt2 = rym_alt1(1:imr,:,:,:) CALL SCATTER_J_BAND(dgrid(region), rm_, rm_alt2, status, jref=j1, rowcom=.true.) IF_NOTOK_RETURN(status=1) rym(i1:i2,j1:j2,:,:) = rm_ if (isRowRoot) rm_alt2 = rzm_alt1(1:imr,:,:,:) CALL SCATTER_J_BAND(dgrid(region), rm_, rm_alt2, status, jref=j1, rowcom=.true.) IF_NOTOK_RETURN(status=1) rzm(i1:i2,j1:j2,:,:) = rm_ #endif ! (2) bands with reduced grid, just use m_uni do lrg=1,nred(region) j = jred(lrg,region) m(:,j,:)= m_uni(:,j,:) end do ! update halo halo_r = mass_dat(region)%halo call UPDATE_HALO_JBAND(dgrid(region), m, m_dat(region)%halo, status ) call UPDATE_HALO_JBAND(dgrid(region), rm(:,:,:,1:nt), halo_r, status ) #ifdef slopes call UPDATE_HALO_JBAND(dgrid(region), rxm, halo_r, status ) call UPDATE_HALO_JBAND(dgrid(region), rym, halo_r, status ) call UPDATE_HALO_JBAND(dgrid(region), rzm, halo_r, status ) #endif deallocate(m_,rm_) endif else ! ! Reduced grid without longitudinal decomposition ! do lrg=1,nred(region) j = jred(lrg,region) redfact=clustsize(j,region) do l=1,lmr do i = imredj(j,region),1,-1 ! the i cell will be distributed within the is:ie array section is = (i-1)*redfact + 1 ie = i*redfact !m_uni is the mass-distribution in the non-reduced grid/divided by !the reduced_grid mass. This is used as distribution function!.... mass=m(i,j,l) m(is:ie,j,l)= m_uni(is:ie,j,l) do n=1,nt rmm = rm(i,j,l,n) rm(is:ie,j,l,n)= m_uni(is:ie,j,l)/mass* rmm #ifdef slopes rmm = rxm(i,j,l,n) rxm(is:ie,j,l,n)= m_uni(is:ie,j,l)/mass * rmm ! rym and rzm are always distributed uniformly: rmm = rym(i,j,l,n) rym(is:ie,j,l,n)= m_uni(is:ie,j,l)/mass * rmm rmm = rzm(i,j,l,n) rzm(is:ie,j,l,n)= m_uni(is:ie,j,l)/mass * rmm #endif end do end do ! update cell(0,...) according to the periodic bc's: rm(0,j,l,1:nt) = rm(im(region),j,l,1:nt) #ifdef slopes rxm(0,j,l,:) = rxm(im(region),j,l,:) rym(0,j,l,:) = rym(im(region),j,l,:) rzm(0,j,l,:) = rzm(im(region),j,l,:) #endif m(0,j,l) = m(im(region),j,l) end do end do endif nullify(m) nullify(rm) #ifdef slopes nullify(rxm) nullify(rym) nullify(rzm) #endif end subroutine red2uni !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: RED2UNI_EM ! ! !DESCRIPTION: transforms data from reduced grid back to uniform grid, ! using EQUAL MASSES instead of reduced mass array (m_uni). !\\ !\\ ! !INTERFACE: ! subroutine red2uni_em(region) ! !USES: ! use dims use meteodata , only : m_dat use tracer_data, only : mass_dat use chem_param, only : ntracet ! ! !INPUT PARAMETERS: ! integer,intent(in) :: region ! ! !REVISION HISTORY: ! written by mike botchev, march-june 1999 ! 20 Dec 2012 - Ph Le Sager - TM5-MP version ! ! !REMARKS: ! !EOP !------------------------------------------------------------------------ !BOC real,dimension(:,:,:,:),pointer :: rm #ifdef slopes real,dimension(:,:,:,:),pointer :: rxm,rym,rzm #endif real,dimension(:,:,:) ,pointer :: m integer i,ie,ii,is,j,l,lrg,n,redfact,lmr real hi,mass,mass_coord,rmm,slope,m_old character(len=5) distr_mode !------------ start lmr=lm(region) m => m_dat(region)%data rm => mass_dat(region)%rm #ifdef slopes rxm => mass_dat(region)%rxm rym => mass_dat(region)%rym rzm => mass_dat(region)%rzm #endif distr_mode = 'unfrm' ! 'slope' or 'unfrm' do lrg=1,nred(region) j = jred(lrg,region) redfact=clustsize(j,region) do l=1,lmr do i = imredj(j,region),1,-1 ! the i cell will be distributed within the is:ie array section is = (i-1)*redfact + 1 ie = i*redfact !m_uni is the mass-distribution in the non-reduced grid/divided by !the reduced_grid mass. This is used as distribution function!.... mass=m(i,j,l); m(is:ie,j,l)= mass/(ie-is+1) if (distr_mode=='unfrm') then ! mixing ratio and x-slope will be UNiFoRMly distributed ! within is:ie do n=1,ntracet rmm = rm(i,j,l,n) rm(is:ie,j,l,n)= m(is:ie,j,l)/mass* rmm #ifdef slopes rmm = rxm(i,j,l,n) rxm(is:ie,j,l,n)= m(is:ie,j,l)/mass * rmm ! rym and rzm are always distributed uniformly: rmm = rym(i,j,l,n) rym(is:ie,j,l,n)= m(is:ie,j,l)/mass * rmm rmm = rzm(i,j,l,n) rzm(is:ie,j,l,n)= m(is:ie,j,l)/mass * rmm #endif end do end if !cmkelseif(distr_mode=='slope') then end do ! update cell(0,...) according to the periodic bc's: rm(0,j,l,1:ntracet) = rm(im(region),j,l,1:ntracet) #ifdef slopes rxm(0,j,l,:) = rxm(im(region),j,l,:) rym(0,j,l,:) = rym(im(region),j,l,:) rzm(0,j,l,:) = rzm(im(region),j,l,:) #endif m(0,j,l) = m(im(region),j,l) end do end do nullify(m) nullify(rm) #ifdef slopes nullify(rxm) nullify(rym) nullify(rzm) #endif end subroutine red2uni_em !EOC #ifdef secmom !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: red2uni_em_2nd ! ! !DESCRIPTION: transforms data from reduced grid back to uniform grid !\\ !\\ ! !INTERFACE: ! subroutine red2uni_em_2nd(region) ! ! !USES: ! use dims use meteodata, only : m_dat use global_data, only : mass_dat use chem_param, only : ntracet ! ! !INPUT PARAMETERS: ! integer,intent(in) :: region ! ! !REVISION HISTORY: ! written by mike botchev, march-june 1999 ! 20 Dec 2012 - Ph Le Sager - TM5-MP version ! ! !REMARKS: ! !EOP !------------------------------------------------------------------------ !BOC real,dimension(:,:,:,:),pointer :: rm,rxm,rym,rzm, rxxm, rxym, rxzm, ryym, ryzm, rzzm real,dimension(:,:,:) ,pointer :: m integer i,ie,ii,is,j,l,lrg,n,redfact,lmr real hi,mass,mass_coord,rmm,slope,m_old character(len=5) distr_mode ! ----- start lmr=lm(region) m => m_dat(region)%data rm => mass_dat(region)%rm rxm => mass_dat(region)%rxm rym => mass_dat(region)%rym rzm => mass_dat(region)%rzm rxxm => mass_dat(region)%rxxm rxym => mass_dat(region)%rxym rxzm => mass_dat(region)%rxzm ryym => mass_dat(region)%ryym ryzm => mass_dat(region)%ryzm rzzm => mass_dat(region)%rzzm distr_mode = 'unfrm' ! 'slope' or 'unfrm' do lrg=1,nred(region) j = jred(lrg,region) redfact=clustsize(j,region) do l=1,lmr do i = imredj(j,region),1,-1 ! the i cell will be distributed within the is:ie array section is = (i-1)*redfact + 1 ie = i*redfact !m_uni is the mass-distribution in the non-reduced grid/divided by !the reduced_grid mass. This is used as distribution function!.... mass=m(i,j,l); m(is:ie,j,l)= mass/(ie-is+1) if (distr_mode=='unfrm') then ! mixing ratio and x-slope will be UNiFoRMly distributed ! within is:ie !#ifdef MPI ! do n=1,ntracetloc !#else do n=1,ntracet !#endif rmm = rm(i,j,l,n) rm(is:ie,j,l,n)= m(is:ie,j,l)/mass* rmm rmm = rxm(i,j,l,n) rxm(is:ie,j,l,n)= m(is:ie,j,l)/mass * rmm ! rym and rzm are always distributed uniformly: rmm = rym(i,j,l,n) rym(is:ie,j,l,n)= m(is:ie,j,l)/mass * rmm rmm = rzm(i,j,l,n) rzm(is:ie,j,l,n)= m(is:ie,j,l)/mass * rmm rmm = rxxm(i,j,l,n) rxxm(is:ie,j,l,n)= m(is:ie,j,l)/mass * rmm rmm = rxym(i,j,l,n) rxym(is:ie,j,l,n)= m(is:ie,j,l)/mass * rmm rmm = rxzm(i,j,l,n) rxzm(is:ie,j,l,n)= m(is:ie,j,l)/mass * rmm rmm = ryym(i,j,l,n) ryym(is:ie,j,l,n)= m(is:ie,j,l)/mass * rmm rmm = ryzm(i,j,l,n) ryzm(is:ie,j,l,n)= m(is:ie,j,l)/mass * rmm rmm = rzzm(i,j,l,n) rzzm(is:ie,j,l,n)= m(is:ie,j,l)/mass * rmm end do end if !cmkelseif(distr_mode=='slope') then end do ! update cell(0,...) according to the periodic bc's: rm(0,j,l,1:ntracet) = rm(im(region),j,l,1:ntracet) rxm(0,j,l,:) = rxm(im(region),j,l,:) rym(0,j,l,:) = rym(im(region),j,l,:) rzm(0,j,l,:) = rzm(im(region),j,l,:) m(0,j,l) = m(im(region),j,l) end do end do nullify(m) nullify(rm) nullify(rxm) nullify(rym) nullify(rzm) nullify(rxxm) nullify(rxym) nullify(rxzm) nullify(ryym) nullify(ryzm) nullify(rzzm) end subroutine red2uni_em_2nd !EOC #endif end module redgridZoom