rtmmod.f90

来自「CCSM Research Tools: Community Atmospher」· F90 代码 · 共 910 行 · 第 1/3 页

#include <misc.h>
#include <preproc.h>

module RtmMod

#if (defined RTM) 

!----------------------------------------------------------------------- 
! 
! Purpose: 
! River Routing Model
! Contains routines: Rtmgridini, Rtmlandini, Rtmriverflux, Rtm
! 
! Method: 
! (U. of Texas River Transport Model) 
!
! Author: Sam Levis
! 
!-----------------------------------------------------------------------
! $Id: RtmMod.F90,v 1.10.6.5 2002/04/27 15:38:55 erik Exp $
!-----------------------------------------------------------------------

  use precision
  use clm_varpar, only : lsmlon, lsmlat, rtmlon, rtmlat 
  implicit none

! RTM grid info

  integer , private :: numlon_r(rtmlat)               !number of lon points at each lat
  real(r8), private, dimension(4) :: rtmedge = (/ 90., 180., -90., -180. /)  !N,E,S,W edges of rtm grid

  real(r8), public, allocatable :: latixy_r(:,:)      !rtm latitudes  of grid cells (degrees)       
  real(r8), public, allocatable :: longxy_r(:,:)      !rtm longitudes of grid cells (degrees)       
  real(r8), public :: area_r(rtmlon,rtmlat)           !rtm gridcell area (km^2)          
  integer , public :: mask_r(rtmlon,rtmlat)           !rtm landmask (land=1,ocean=0)

! land model to RTM mapping. for each rtm grid cell:

  integer , private :: mxovr_s2r                      !max number of overlapping cells
  integer , private :: novr_s2r(rtmlon,rtmlat)        !number    of overlapping cells
  integer , private, allocatable :: iovr_s2r(:,:,:)   !lon index of overlapping cells
  integer , private, allocatable :: jovr_s2r(:,:,:)   !lat index of overlapping cells
  real(r8), private, allocatable :: wovr_s2r(:,:,:)   !weight    of overlapping cells

! RTM runoff for coupled communication

  integer , public, allocatable :: ocnrof_iindx(:)    !rtm longitude index of ocean runoff point
  integer , public, allocatable :: ocnrof_jindx(:)    !rtm latitude index of ocean runoff point
  real(r8), public, allocatable :: ocnrof_vec(:)      !rtm runoff vector (1/2 deg grid, kg/m^2/s)

! RTM history file variables

  real(r8), public, allocatable :: qchan2(:)          !river (channel) flow (m**3 H2O /s)
  real(r8), public, allocatable :: qchocn2(:)         !river (channel) flow into ocean (m**3/s)

! time averaging for rtm calculatino

  real(r8), public, allocatable :: totrunin_ave(:)    !time averaged vector of input fluxes
  real(r8), public, allocatable :: prec_ave(:)        !time averaged vector of precipitation
  real(r8), public, allocatable :: evap_ave(:)        !time averaged vector of evaporation
  real(r8), public :: delt_rtm                        !rtm time step
  integer , public :: ncount_rtm                      !number of time samples to average over

! fluxes

  integer , private :: rdirc(0:rtmlon+1,0:rtmlat+1)   !rtm river flow direction (0-8)
  real(r8), private :: fluxout(0:rtmlon+1,0:rtmlat+1) !water flux out of cell (m^3/s)
  real(r8), private :: ddist(rtmlon,rtmlat)           !downstream distance (m)
  real(r8), private :: rivarea(rtmlon,rtmlat)         !cell area (m^2)
  real(r8), public  :: volr(rtmlon,rtmlat)            !water volume in cell (m^3)

  real(r8), private, allocatable :: latsh(:)          !southern edge of cells at rtm grid     
  real(r8), private, allocatable :: lonwh(:,:)        !western  edge of cells at rtm grid     

! inputs to RTM at 1/2 degree resolution

  real(r8), private :: totrunin_r(rtmlon,rtmlat)      !surface runoff (mm/s)

! outputs returned from RTM at 1/2 degree resolution

  real(r8), private :: flxlnd_r(rtmlon,rtmlat)        !river flux (m**3/s)
  real(r8), private :: flxocn_r(rtmlon,rtmlat)        !river flux to the ocean (m**3/s)
  real(r8), private :: dvolrdt_r(rtmlon,rtmlat)       !change in storage (mm/s)
  real(r8), private :: volrtm(rtmlon,rtmlat)          !change in storage (m**3/s)
  real(r8), private :: runrtm(rtmlon,rtmlat)          !input runoff on rtm grid (m**3/s)

! RTM water flux into cell

  real(r8), private :: sfluxin(rtmlon,rtmlat)         !water flux into cell (m3/s)

! global averaging

  character(len=*),parameter :: F40="('(diag) ',a17,'    date  ', &
  &     '   prec        evap        runoff(lnd)   runoff(rtm) dvoldt(rtm) runoff-ocn(rtm)  (m^3/sec)')"
  character(len=*),parameter :: F41="('(diag) ',a17,'   nstep  ', &
  &     '   prec        evap        runoff(lnd)   runoff(rtm) dvoldt(rtm) runoff-ocn(rtm)  (m^3/sec)')"
  character(len=*),parameter :: F21="('(diag) ',a17,' ----------------------', &
  &     7('----------'))"
  character(len=*),parameter :: F22="('(diag) ',a17,i8,6(d13.4))"

  real(r8) prec_global            !total precipitation (m^3/sec) 
  real(r8) evap_global            !total evaporation (m^3/sec)
  real(r8) runlnd_global          !total input runoff on land grid (m^3/sec)
  real(r8) runrtm_global          !total input runoff on rtm grid (m^3/sec)
  real(r8) ocnrtm_global          !total ocean runoff on rtm grid (m^3/sec)
  real(r8) volrtm_global          !total change in storage on rtm (m^3/sec)
  integer  ncount_global          !global counter 
  integer  yrold                  !old year

  SAVE

!=======================================================================
CONTAINS
!=======================================================================

subroutine Rtmgridini

!----------------------------------------------------------------------- 
! 
! Purpose: 
! Initialize RTM grid and land mask (U. of Texas River Transport Model)
! 
! Method: 
! 
! Author: Sam Levis
! 
!-----------------------------------------------------------------------

  use spmdMod   , only : masterproc
  use areaMod   , only : celledge, cellarea   
  use clm_varctl, only : frivinp_rtm
  use clm_varcon, only : re
  use shr_const_mod, only: SHR_CONST_PI

! ------------------------ local variables ---------------------------
  integer  :: ioff(0:8) = (/0,0,1,1,1,0,-1,-1,-1/) !calc dist as in hydra
  integer  :: joff(0:8) = (/0,1,1,0,-1,-1,-1,0,1/) !of grid cell down stream
  integer  :: i,j,k,n                       !loop indices
  integer  :: i2,j2                         !downstream i and j
  real(r8) :: deg2rad                       !pi/180
  real(r8) :: dx                            !lon dist. between grid cells (m)
  real(r8) :: dy                            !lat dist. between grid cells (m)
  real(r8) :: dist(rtmlon,rtmlat)           !dist. of the grid cell down stream (m) 
  real(r8) :: tempg(rtmlon,rtmlat)          !temporary buffer
  integer  :: tempgp(0:rtmlon+1,0:rtmlat+1) !temporary buffer  
! --------------------------------------------------------------------
  
  if (masterproc) then

! --------------------------------------------------------------------
! Useful constants and initial values
! --------------------------------------------------------------------

     write(6,*)'Columns in RTM = ',rtmlon
     write(6,*)'Rows in RTM    = ',rtmlat

     allocate(latixy_r(rtmlon,rtmlat))
     allocate(longxy_r(rtmlon,rtmlat))
     allocate(latsh(rtmlat+1))        !southern edge of cells at rtm grid     
     allocate(lonwh(rtmlon+1,rtmlat)) !western  edge of cells at rtm grid     
     
     deg2rad = SHR_CONST_PI / 180.
     volr = 0.

! --------------------------------------------------------------------
! Open and read input data (river direction file)
! rtm operates from south to north and from the dateline
! --------------------------------------------------------------------

     open (1,file=frivinp_rtm)
     write(6,*)'opened river direction data'
     do j = 1,rtmlat
        numlon_r(j) = 0
        do i = 1,rtmlon
           read(1,*) latixy_r(i,j),longxy_r(i,j),tempg(i,j)
           if (longxy_r(i,j) /= 1.e36) numlon_r(j) = numlon_r(j) + 1
           tempgp(i,j) = nint(tempg(i,j))
        enddo
     enddo
     close(1)
     write(6,*)'closed river direction data'
     write(6,*)
     
! --------------------------------------------------------------------
! Determine RTM celledges, areas and interpolation masks
! --------------------------------------------------------------------

     call celledge (rtmlat    , rtmlon    , numlon_r  , longxy_r  , &
                    latixy_r  , rtmedge(1), rtmedge(2), rtmedge(3), &
                    rtmedge(4), latsh     , lonwh     )

     call cellarea (rtmlat    , rtmlon    , numlon_r  , latsh     , lonwh , &
                    rtmedge(1), rtmedge(2), rtmedge(3), rtmedge(4), area_r) 

! --------------------------------------------------------------------
! Determine rtm mask, downstream distance and area
! --------------------------------------------------------------------

! determine rtm ocn/land mask

     do i=1,rtmlon
        tempgp(i,0)        = tempgp(mod(i+rtmlon/2-1,rtmlon)+1,1)
        tempgp(i,rtmlat+1) = tempgp(mod(i+rtmlon/2-1,rtmlon)+1,rtmlat)
        if (tempgp(i,0)        /= 0) tempgp(i,0)        = mod(tempgp(i,0)       +4-1,8)+1
        if (tempgp(i,rtmlat+1) /= 0) tempgp(i,rtmlat+1) = mod(tempgp(i,rtmlat+1)+4-1,8)+1
     enddo
     do j=0,rtmlat+1
        tempgp(0,j) =tempgp(rtmlon,j)
        tempgp(rtmlon+1,j)=tempgp(1,j)
     enddo
     
     do j=0,rtmlat+1
        do i=0,rtmlon+1
           rdirc(i,j)=tempgp(i,j)
        enddo
     enddo
     
     do j=1,rtmlat
        do i=1,rtmlon
           if (rdirc(i,j) == 0) then
              mask_r(i,j) = 0
           else
              mask_r(i,j) = 1
           end if
        enddo
     enddo

! determine downstream distance - instead of reading a distance file 
! calculate the downstream distance as in hydra

     do j=1,rtmlat
        do i=1,rtmlon
           i2 = i + ioff(tempgp(i,j))
           j2 = j + joff(tempgp(i,j))
           if (i2 == 0) i2 = 2                 !avoids i2 out of bounds in the following
           if (i2 == rtmlon+1) i2 = rtmlon-1   !avoids i2 out of bounds in the following  
           dy = deg2rad * abs(latixy_r(i,j)-latixy_r(i2,j2)) * re*1000.
           dx = deg2rad * abs(longxy_r(i,j)-longxy_r(i2,j2)) * re*1000. &
                *0.5*(cos(latixy_r(i,j)*deg2rad)+cos(latixy_r(i2,j2)*deg2rad))
           dist(i,j) = sqrt(dx*dx + dy*dy)
           ddist(i,j) = dist(i,j)
           rivarea(i,j)=1.e6 * area_r(i,j)     !convert into m**2
        enddo
     enddo
  
  endif  ! end of if-masterproc block

end subroutine Rtmgridini

!=======================================================================

subroutine Rtmlandini

!----------------------------------------------------------------------- 
! 
! Purpose: 
! Initialize RTM-land interpolation weights (U. of Texas River Transport Model)
! and variables related to runoff time averaging
! 
! Method: 
! 
! Author: Mariana Vertenstein
! 
!-----------------------------------------------------------------------

  use spmdMod     , only : masterproc
  use areaMod     , only : areaini            
  use clm_varsur  , only : numlon, area, lats, lonw, landmask 
  use clm_varmap  , only : numpatch
  use time_manager, only : get_curr_date

! ------------------------ local variables ---------------------------
  integer  :: i,j,k,n                    !loop indices
  integer  :: is,js                      !land model grid indices
  real(r8) :: maskone_s(lsmlon,lsmlat)   !dummy field: see below                 
  real(r8) :: maskone_r(rtmlon,rtmlat)   !dummy field: see below                 
  integer  :: ocnrof_mask(rtmlon,rtmlat) !rtm mask for ocean points with possible nonzero runoff
  integer  :: ocnrof_num                 !number of valid ocean points with possible nonzero runoff
  integer  :: yrnew                      !year (0, ...)
  integer  :: mon                        !month (1, ..., 12) 
  integer  :: day                        !day of month (1, ..., 31)
  integer  :: ncsec                      !seconds of current date
! --------------------------------------------------------------------

  if (masterproc) then

! --------------------------------------------------------------------
! The following section allows RTM and land model to coexist at different
! horizontal resolutions
! --------------------------------------------------------------------

     write(6,*)
     write(6,*) 'Initializing area-averaging interpolation for RTM.....'

! To find fraction of each land model grid cell that is land based on rtm grid.
! For this purpose, want all rtm grid cells to contribute to grid cell 
! average on land model grid, i.e., all cells used regardless of whether land 
! or ocean. Do this by setting [maskone_s] = 1 

! [maskone_s] = 1 means all grid cells on land model grid, regardless of whether
! land or ocean, will contribute to rtm grid.

     do j = 1, lsmlat
        do i = 1, numlon(j)