!### 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 advect use GO, only : gol, goPr, goErr implicit none ! ------------ interface --------------------- private public :: Advect_Init, Advect_Done public :: dynam0 ! --- const -------------------------------------- character(len=*), parameter :: mname = 'advect' contains ! ================================================================ subroutine Advect_Init( status ) use dims , only : nregions use Meteo, only : Set use Meteo, only : sp1_dat, sp2_dat, sp_dat #ifdef with_prism use Meteo, only : tsp_dat #endif use Meteo, only : mfu_dat, mfv_dat, mfw_dat use Meteo, only : pu_dat, pv_dat, pw_dat ! --- in/out -------------------------------- integer, intent(out) :: status ! --- const ------------------------------ character(len=*), parameter :: rname = mname//'/Advect_Init' ! --- local -------------------------------- integer :: n ! --- begin -------------------------------- ! loop over all (zoom) regions: do n = 1, nregions ! surface pressures and mass fluxes for advection call Set( sp_dat(n), status, used=.true. ) call Set( sp1_dat(n), status, used=.true. ) call Set( sp2_dat(n), status, used=.true. ) #ifdef with_prism call Set( tsp_dat(n), status, used=.true. ) #endif call Set( mfu_dat(n), status, used=.true. ) call Set( mfv_dat(n), status, used=.true. ) call Set( mfw_dat(n), status, used=.true. ) call Set( pu_dat(n), status, used=.true. ) call Set( pv_dat(n), status, used=.true. ) call Set( pw_dat(n), status, used=.true. ) end do ! regions ! ok status = 0 end subroutine Advect_Init ! *** subroutine Advect_Done( status ) ! --- in/out -------------------------------- integer, intent(out) :: status ! --- const ------------------------------ character(len=*), parameter :: rname = mname//'/Advect_Done' ! --- begin -------------------------------- ! nothing to be done ! ok status = 0 end subroutine Advect_Done ! ================================================================ !----------------------------------------------------------------------- ! !**** dynam0 - calculates vertical massfluxes v 9.1 ! ! programmed by mh mpi HH 1-oct-1991 ! ! purpose ! ------- ! This sr is CALLed each time after new massfluxes are read in. ! Calculate the new air-masses in each grid box from the ! surface pressure. ! Calculate the amount of air moved in each advection substep, ! also in the vertical direction. ! ! interface ! --------- ! CALL dynam0 ! ! method ! ------ ! integrate air conservation equation in the vertical direction ! in order to obtain the vertical massfluxes. ! ! externals ! --------- ! none ! ! reference ! --------- ! see manual ! ! modified by mk (1999) for zoom version. !----------------------------------------------------------------------- subroutine dynam0(region) use dims, only : im, jm, lm, ndyn, tref, zero, bt, xcyc use global_data, only : wind_dat use meteo , only : pu_dat, pv_dat #ifdef MPI use mpi_const, only : ntracetloc use mpi_comm, only : barrier_t #endif implicit none ! input integer,intent(in) :: region ! local real,dimension(:,:,:),pointer :: pu,pv,am,bm,cm real,dimension(im(region),jm(region)) :: pit real,dimension(im(region),jm(region),lm(region)-1) :: sd real,dimension(im(region),jm(region),lm(region)) :: conv_adv real :: dtu, dtv, dtw integer :: i, j, l, lmr ! start #ifdef MPI if ( ntracetloc == 0 ) return !WP! only processors with nonzero ntracet #endif ! leave if fluxes are not in use: if ( .not. pu_dat(region)%used ) return if ( .not. pv_dat(region)%used ) return ! set pointers: pu => pu_dat(region)%data pv => pv_dat(region)%data am => wind_dat(region)%am_t bm => wind_dat(region)%bm_t cm => wind_dat(region)%cm_t !---- length of advection substeps (in seconds) in the three directions dtu=ndyn/(2.*tref(region)) ! again removed the revert! cmk dtv=ndyn/(2.*tref(region)) dtw=ndyn/(2.*tref(region)) ! number of levels in this region: lmr = lm(region) ! init to zero: am = zero bm = zero cm = zero ! ! compute conv_adv, the horizontal mass convergence ! the arrays pu and pv contain the horizontal air mass fluxes crossing ! the grid box boundaries. pu(i,j,l) is the eastward mass flux in kg/sec ! at the eastern edge of box i,j,l; pv(i,j,l) is the northward ! mass flux at the southern edge of box i,j,l ! do l=1,lmr do j=1,jm(region) do i=1,im(region) conv_adv(i,j,l)=pu(i-1,j,l)-pu(i,j,l)+pv(i,j,l)-pv(i,j+1,l) end do end do end do ! ! compute pit, the vertiCALLy integrated conv_adv ! do j=1,jm(region) do i=1,im(region) pit(i,j)=conv_adv(i,j,1) do l=2,lmr pit(i,j)=pit(i,j)+conv_adv(i,j,l) end do ! ! compute the vertical massflux on the box boundaries (in kg/s) !mk in the hybrid system, the tendency in ps is given by bt. Bt is !mk already normalized. Now distribute pit such that the new mass !mk distribution after advection exactly matches the !mk exactly coincides with the distribution that is obtained !mk from the new surface pressure. ! sd(i,j,lmr-1)=conv_adv(i,j,lmr) - (bt(lmr)-bt(lmr+1))*pit(i,j) do l=lmr-2,1,-1 sd(i,j,l)=sd(i,j,l+1)+conv_adv(i,j,l+1)-(bt(l+1)-bt(l+2))*pit(i,j) end do end do end do ! ! compute amount of air moved each advection substep, am, bm or cm (kg) ! ! am(i,j,l) is the amount of air moved eastward at the eastern ! boundary of grid box i,j,l ! compute cm ! cm(i,j,l) is the amount of air moved in upward direction at the ! upper boundary of grid box i,j,l ! compute bm ! bm(i,j,l) is the amount of air moved northward at the southern ! boundary of grid box i,j,l ! do l=1,lmr do j=1,jm(region) do i=0,im(region) am(i,j,l)=dtu*pu(i,j,l) end do end do end do do l=1,lmr do j=1,jm(region)+1 do i=1,im(region) bm(i,j,l)=dtv*pv(i,j,l) end do end do end do do l=1,lmr-1 do j=1,jm(region) do i=1,im(region) cm(i,j,l)=-dtw*sd(i,j,l) end do end do end do if ( xcyc(region) == 1 ) then am(0,:,:) = am(im(region),:,:) am(im(region)+1,:,:) = am(1,:,:) end if nullify(pu) nullify(pv) nullify(am) nullify(bm) nullify(cm) #ifdef MPI call barrier_t #endif end subroutine dynam0 end module advect