inidat.f90

CCSM Research Tools: Community Atmosphere Model (CAM)

#include <misc.h>
#include <params.h>

module inidat
!BOP
!
! !MODULE: inidat --- dynamics-physics coupling module
!
! !USES:
   use precision
   use comspe
   use chemistry, only: chem_init_mix

! !PUBLIC MEMBER FUNCTIONS:
   public read_inidat, copy_inidat

! !PUBLIC DATA MEMBERS:
   real(r8), allocatable :: ps_tmp(:,:)
   real(r8), allocatable :: u3s_tmp(:,:,:)
   real(r8), allocatable :: v3s_tmp(:,:,:)
   real(r8), allocatable :: uv_local(:,:,:)
   real(r8), allocatable :: t3_tmp(:,:,:)
   real(r8), allocatable :: q3_tmp(:,:,:,:)
   real(r8), allocatable :: q3_local(:,:,:,:)
   real(r8), allocatable :: qcwat_tmp(:,:,:) 
   real(r8), allocatable :: lcwat_tmp(:,:,:) 
   real(r8), allocatable :: phis_tmp(:,:)
   real(r8), allocatable :: landfrac_tmp(:,:)                  
   real(r8), allocatable :: landm_tmp(:,:)
   real(r8), allocatable :: sgh_tmp(:,:)
   real(r8), allocatable :: ts_tmp(:,:)
   real(r8), allocatable :: tsice_tmp(:,:)
   real(r8), allocatable :: tssub_tmp(:,:,:)
   real(r8), allocatable :: sicthk_tmp(:,:)
   real(r8), allocatable :: snowhice_tmp(:,:)                

   real(r8) zgsint_tmp

!
! !DESCRIPTION:
!
!      This module provides 
!
!      \begin{tabular}{|l|l|} \hline \hline
!        read\_inidat    &   \\ \hline
!        copy\_inidat    &   \\ \hline 
!                                \hline
!      \end{tabular}
!
! !REVISION HISTORY:
!   YY.MM.DD   ?????      Creation
!   00.06.01   Grant      First attempt at modifying for LRDC
!   01.10.01   Lin        Various revisions
!   01.01.15   Sawyer     Bug fixes for SPMD mode
!   01.03.26   Sawyer     Added ProTeX documentation
!
!EOP
!-----------------------------------------------------------------------

contains

!-----------------------------------------------------------------------
!BOP
! !IROUTINE: read_inidat --- read initial dataset
!
! !INTERFACE: 
   subroutine read_inidat

! !USES:
      use precision
      use pmgrid
      use pspect
      use rgrid
      use comsrf,    only: plevmx,srfflx_state
      use commap
      use physconst, only: gravit
      use history, only: fillvalue
      use constituents, only: pcnst, pnats, cnst_name, qmin
      use tracers, only: nusr_adv, nusr_nad, ixuadv, ixunad, ixcldw

      implicit none

      include 'netcdf.inc'

!------------------------------Commons----------------------------------

#include <comctl.h>
#include <comqfl.h>
#include <comlun.h>
#include <perturb.h>

! !DESCRIPTION:
!
!   Read initial dataset and spectrally truncate as appropriate.
!
! !REVISION HISTORY:
!
!   00.06.01   Grant      First attempt at modifying for LRDC
!   00.10.01   Lin        Various revisions
!   01.01.15   Sawyer     Bug fixes for SPMD mode
!   01.03.09   Eaton      Modifications
!   01.03.26   Sawyer     Added ProTeX documentation
!
!EOP
!-----------------------------------------------------------------------
!BOC
!
! !LOCAL VARIABLES:

      integer i,j,k,m,lat       ! grid and constituent indices
      integer ihem              ! hemisphere index
      real(r8) pdelb(plond,plev)! pressure diff between interfaces
      real(r8) pertval          ! perturbation value
      real(r8) zgssum           ! partial sums of phis
      integer ii, ic
!
! Netcdf related variables
!
      integer lonsiz, latsiz, levsiz ! Dimension sizes
      integer londimid, levdimid, latdimid ! Dimension ID's
      integer tid, qid  ! Variable ID's
      integer tracid(pcnst+pnats) ! Variable ID's
      integer phisid, sghid, psid ! Variable ID's
      integer landmid
#if ( ! defined COUP_CSM )
      integer tsid, ts1id, ts2id, ts3id, ts4id,tsiceid ! Variable ID's
#endif
      integer sicid,  snowhiceid           ! Variable ID's
      integer landfracid        ! Variable ID's
      integer usid, vsid
      integer strt2d(3)         ! start lon, lat, time indices for netcdf 2-d
      integer strt3d(4)         ! start lon, lev, lat, time for netcdf 3-d
      data strt2d/3*1/          ! Only index 2 will ever change
      data strt3d/4*1/          ! Only indices 2,3 will ever change

      integer cnt2d(3)          ! lon, lat, time counts for netcdf 2-d
      integer cnt3d(4)          ! lon, lat, lev, time counts for netcdf 2-d
      data cnt2d/plon,1,1/      ! 2-d arrs: Always grab only a "plon" slice
      data cnt3d/plon,plev,plat,1/ ! 3-d arrs: Always grab a full time slice

      integer ndims2d           ! number of dimensions
      integer dims2d(NF_MAX_VAR_DIMS) ! variable shape
      integer ndims3d           ! number of dimensions
      integer dims3d(NF_MAX_VAR_DIMS) ! variable shape
      integer tmptype
      integer natt, ret, attlen ! netcdf return values
      logical phis_hires        ! true => PHIS came from hi res topo
      real(r8) arr3d(plon,plev,plat)

      character*(NF_MAX_NAME) tmpname
      character*256 text
      character*80 trunits      ! tracer untis

      real(r8) splon_arr3d(plon,plev,plat)
!     real(r8) splat_arr3d(plon,plev,splat)   ! Glenn Grant's original code
      real(r8) splat_arr3d(plon,plev,plat-1)  ! temporary patch until splat = plat-1

      integer slatid, slatdimid, slatsiz
      integer slonid, slondimid, slonsiz
      integer cnt3dus(4)        ! index counts for netcdf U staggered grid
      integer cnt3dvs(4)        ! index counts for netcdf V staggered grid

!      data cnt3dus/plon,plev,splat,1/ ! 3-d arrs: Always grab a full time slice
! SJL
      integer platm1
      parameter (platm1=plat-1)
      data cnt3dus/plon,plev,platm1,1/ ! temporary patch
      data cnt3dvs/plon,plev,plat,1/ ! 3-d arrs: Always grab a full time slice

!
!-----------------------------------------------------------------------
! Allocate memory for temporary arrays
!-----------------------------------------------------------------------
!
! Note if not masterproc still might need to allocate array for spmd case
! since each processor calls MPI_scatter 
!
      allocate ( ps_tmp(plond,plat) )
      allocate ( u3s_tmp(plon,plat,plev) )
      allocate ( v3s_tmp(plon,plat,plev) )
      allocate ( t3_tmp(plond,plev,plat) )             
      allocate ( q3_tmp(plond,plev,pcnst+pnats,plat) ) 
      allocate ( qcwat_tmp(plond,plev,plat) ) 
      allocate ( lcwat_tmp(plond,plev,plat) ) 
      allocate ( phis_tmp(plond,plat) )        
      allocate ( landm_tmp(plond,plat) )                  
      allocate ( sgh_tmp(plond,plat) )                 
      allocate ( ts_tmp(plond,plat) )                   
      allocate ( tsice_tmp(plond,plat) )                   
      allocate ( tssub_tmp(plond,plevmx,plat) )         
      allocate ( sicthk_tmp(plond,plat) )               
      allocate ( snowhice_tmp(plond,plat) )                
      allocate ( landfrac_tmp(plond,plat) )                
!
!-----------------------------------------------------------------------
! Read in input variables
!-----------------------------------------------------------------------

      if (masterproc) then
!
! Get dimension IDs and lengths 
!
         call wrap_inq_dimid  (ncid_ini, 'lat', latdimid)
         call wrap_inq_dimlen (ncid_ini, latdimid, latsiz)
         call wrap_inq_dimid  (ncid_ini, 'lev', levdimid)
         call wrap_inq_dimlen (ncid_ini, levdimid, levsiz)
         call wrap_inq_dimid  (ncid_ini, 'lon', londimid)
         call wrap_inq_dimlen (ncid_ini, londimid, lonsiz)
         call wrap_inq_dimid  (ncid_ini, 'slat', slatdimid)
         call wrap_inq_dimlen (ncid_ini, slatdimid, slatsiz)
         call wrap_inq_dimid  (ncid_ini, 'slon', slondimid)
         call wrap_inq_dimlen (ncid_ini, slondimid, slonsiz)

!
! Get variable id's 
! Check that all tracer units are in mass mixing ratios
!
!         call wrap_inq_varid (ncid_ini, 'U'   , uid)
!         call wrap_inq_varid (ncid_ini, 'V'   , vid)

         call wrap_inq_varid (ncid_ini, 'slat', slatid)
         call wrap_inq_varid (ncid_ini, 'slon', slonid)
         call wrap_inq_varid (ncid_ini, 'US'  , usid)
         call wrap_inq_varid (ncid_ini, 'VS'  , vsid)

         call wrap_inq_varid (ncid_ini, 'T'   , tid)
         call wrap_inq_varid (ncid_ini, 'Q'   , qid)
         call wrap_inq_varid (ncid_ini, 'PS'  , psid)
         call wrap_inq_varid (ncid_ini, 'PHIS', phisid)
         call wrap_inq_varid (ncid_ini, 'SGH' , sghid)
         call wrap_inq_varid (ncid_ini, 'LANDM', landmid)
#if ( ! defined COUP_CSM )
!
! For land-fraction check if the variable name LANDFRAC is on the dataset if not assume FLAND
!
         if ( nf_inq_varid(ncid_ini, 'LANDFRAC', landfracid ) == NF_NOERR ) then
            call wrap_inq_varid (ncid_ini, 'LANDFRAC', landfracid)
         else
            call wrap_inq_varid (ncid_ini, 'FLAND', landfracid)
         end if
         call wrap_inq_varid (ncid_ini, 'TS', tsid)
         call wrap_inq_varid (ncid_ini, 'TSICE', tsiceid)
         call wrap_inq_varid (ncid_ini, 'TS1', ts1id)
         call wrap_inq_varid (ncid_ini, 'TS2', ts2id)
         call wrap_inq_varid (ncid_ini, 'TS3', ts3id)
         call wrap_inq_varid (ncid_ini, 'TS4', ts4id)
         call wrap_inq_varid (ncid_ini, 'SNOWHICE', snowhiceid)
#if ( defined COUP_SOM )
         call wrap_inq_varid (ncid_ini, 'SICTHK', sicid)
#endif
#endif
         if (readtrace) then
            do m=2,pcnst+pnats
               call wrap_inq_varid (NCID_INI,cnst_name(m), tracid(m))
               call wrap_get_att_text (NCID_INI,tracid(m),'units', trunits)
               if (trunits(1:5) .ne. 'KG/KG' .and. trunits(1:5) .ne. 'kg/kg') then
                  write(6,*)'INIDAT: tracer units for tracer = ', &
                            cnst_name(m),' must be in KG/KG'
                  call endrun
               endif
            end do
         end if
!
! Check dimension ordering for one 2-d and one 3-d field.
! Assume other arrays of like rank will have dimensions ordered the same.
!
         call wrap_inq_var (ncid_ini, psid, tmpname, tmptype, &
                            ndims2d, dims2d, natt)
         if (dims2d(1).ne.londimid .or. dims2d(2).ne.latdimid .or. &
             ndims2d.gt.3) then
            write(6,*)'INIDAT: Bad number of dims or ordering on 2d fld'
            call endrun
         end if
!
! Check for presence of 'from_hires' attribute to decide whether to filter
!
         ret = nf_inq_attlen (ncid_ini, phisid, 'from_hires', attlen)
         if (ret.eq.NF_NOERR .and. attlen.gt.256) then
            write(6,*)'INIDAT: from_hires attribute length is too long'
            call endrun
         end if
         ret = nf_get_att_text (ncid_ini, phisid, 'from_hires', text)
         if (ret.eq.NF_NOERR .and. text(1:4).eq.'true') then
            phis_hires = .true.
!            write(6,*)'INIDAT: Will filter input PHIS: attribute ', &
!                      'from_hires is true'
         else
            phis_hires = .false.
!            write(6,*)'INIDAT: Will not filter input PHIS: attribute ', &
!                      'from_hires is either false or not present'
         end if
!
! Read in 2d fields.  
! For stand alone run: get surface temp and 4 (sub)surface temp fields
! For stand alone run with slab-ocean: get sea-ice thickness and snow cover
!
         do j=1,plat
            strt2d(2) = j
            if (ideal_phys .or. aqua_planet) then
               do i=1,nlon(j)
                  phis_tmp(i,j) = 0.
                  sgh_tmp (i,j) = 0.
               end do
            else
               call wrap_get_vara_realx (ncid_ini, phisid, strt2d, cnt2d, &
                                         phis_tmp(1,j))
               call wrap_get_vara_realx (ncid_ini, sghid , strt2d, cnt2d, &
                                         sgh_tmp(1,j))
            endif
            call wrap_get_vara_realx (ncid_ini, landmid, strt2d, cnt2d, &
                                      landm_tmp(1,j))
            call wrap_get_vara_realx (ncid_ini, psid, strt2d, cnt2d, &
                                      ps_tmp(1,j))
#if ( ! defined COUP_CSM )
            if (aqua_planet) then
               do i=1,nlon(j)
                  landfrac_tmp(i,j) = 0.
               end do
            else
               call wrap_get_vara_realx (ncid_ini, landfracid, strt2d, cnt2d, &
                                         landfrac_tmp(1,j))
            endif
            call wrap_get_vara_realx (ncid_ini, tsid, strt2d, cnt2d, &
                                           ts_tmp(1,j))
            call wrap_get_vara_realx (ncid_ini, tsiceid, strt2d, cnt2d, &
                                           tsice_tmp(1,j))
            call wrap_get_vara_realx (ncid_ini, ts1id, strt2d, cnt2d, &
                                           tssub_tmp(1,1,j))
            call wrap_get_vara_realx (ncid_ini, ts2id, strt2d, cnt2d, &
                                           tssub_tmp(1,2,j))
            call wrap_get_vara_realx (ncid_ini, ts3id, strt2d, cnt2d, &
                                           tssub_tmp(1,3,j))
            call wrap_get_vara_realx (ncid_ini, ts4id, strt2d, cnt2d, &
                                           tssub_tmp(1,4,j))
!
! Set sea-ice thickness and snow cover:
!
#if ( defined COUP_SOM )
            call wrap_get_vara_realx(ncid_ini, sicid, strt2d, cnt2d, sicthk_tmp(1,j))
#endif
            call wrap_get_vara_realx(ncid_ini, snowhiceid, strt2d, cnt2d, snowhice_tmp(1,j))
#endif
         end do
!
! Read in 3d fields.  

! Staggered grid variables and transpose

         call wrap_get_vara_realx(ncid_ini,usid, strt3d, cnt3dus, splat_arr3d)

         do k = 1, plev
!
! SJL: initialize j=1 because later on u3s_tmp will be copied to u3s using f90 array syntax
               do i = 1, plon
                  u3s_tmp(i,1,k) = fillvalue
               enddo
!
            do j = 1, plat-1
               do i = 1, plon
                  u3s_tmp(i,j+1,k) = splat_arr3d(i,k,j)
               enddo
            enddo
         enddo

         call wrap_get_vara_realx(ncid_ini,vsid, strt3d, cnt3dvs, splon_arr3d)