GFDL - Geophysical Fluid Dynamics Laboratory

module barotropic_dynamics_mod

Contact: Isaac Held

Reviewers: Peter Phillipps

OVERVIEW

The dynamical core of the spectral transform model for two-dimensional, non-divergent flow on the surface of the sphere.

DESCRIPTION

Integrates the barotropic vorticity equation for nondivergent flow on the sphere using the spectral transform technique. Also allows for the inclusion of a passive tracer advected by the same spectral advection algorithm as the vorticity, and a gridpoint tracer advected with a finite volume algorithm on the transform grid. The default initial condition provided as an example is a zonal flow resembling that in the Northern winter, plus a sinusoidal disturbance localized in midlatitudes. For a full description of the model and algorithms used, see barotropic.pdf For higher level routines for running this barotropic spectral model, see atmosphere_mod

OTHER MODULES USED

  • fms_mod
  • constants_mod
  • time_manager_mod
  • transforms_mod
  • spectral_damping_mod
  • leapfrog_mod
  • fv_advection_mod

PUBLIC INTERFACE

  use barotropic_dynamics_mod [,only: barotropic_dynamics_init,       
                                      barotropic_dynamics,
                          barotropic_dynamics_end,
                                      dynamics_type,
                      grid_type,
                      spectral_type,
                      tendency_type]

PUBLIC DATA


type grid_type
   real, pointer, dimension(:,:,:) :: u, v, vor, trs, tr, pv
   real, pointer, dimension(:,:)   :: stream
end type

   allocated space for grid fields

   (:,:,:) => (lon, lat, time_level)
   (:,:)   => (lon, lat)
      (lon, lat) on local computational domain
      time_level stores the two time levels needed for the
          leapfrog step
   
   u      -- eastward velocity (m/s)
   v      -- northward velocity (m/s)
   vor    -- vorticity (1/s)
   trs    -- tracer advected spectrally
   tr     -- tracer advected on grid
   pv     -- absolute vorticity, f + vor, where f = 2*omega*sin(lat) (1/s)
   stream -- streamfunction (m^2/s) at current time

type spectral_type
   complex, pointer, dimension(:,:,:) :: vor, trs
end type

   allocated space for spectral fields
   
      (:,:,:) => (zonal, meridional, time_level)

   vor -- spectral vorticity
   trs -- spectral tracer

type tendency_type
    real, pointer, dimension(:,:) :: u, v, trs, tr
end type

   allocated space for accumulating tendencies, d/dt, in grid space, 
           for prognostic variables
   
      (:,:,:) => (lon, lat)

type dynamics_type
   type(grid_type)     :: grid
   type(spectral_type) :: spec
   type(tendency_type) :: tend
   integer             :: num_lon, num_lat  ! size of global domain
   logical             :: grid_tracer, spec_tracer 
end type

   grid_tracer = .true. => tracer with gridpoint advection is being integrated
   similarly for spec_tracer

PUBLIC ROUTINES

subroutine  barotropic_dynamics_init

subroutine  barotropic _dynamics

subroutine  barotropic_dynamics_end

type (grid_type)

type (spectral_type)

type (tendency_type)

type (dynamics_type)

 subroutine barotropic_dynamics_init(Dyn,  Time, Time_init)
 
   type(dynamics_type), intent(inout)  :: Dyn
         type containing all dynamical fields and related information
     (see type (dynamics_type))

type(time_type)    , intent(in)     :: Time, Time_init

         current time and time at which integeration began

     time_type defined by time_manager_mod

   Initializes the module;

   Reads restart from 'INPUT/barotropic_dynamics.res' if Time = Time_init;

     otherwise uses default initial conditions


  subroutine barotropic_dynamics &
     (Time, Time_init, Dyn, previous, current, future, delta_t)
 
     type(time_type)    , intent(inout)  :: Time, Time_init
     type(dynamics_type), intent(inout)  :: Dyn
     integer            , intent(in   )  :: previous, current, future
     real               , intent(in   )  :: delta_t
      
     previous, current and future = 1 or 2
       these integers refer to the third dimension of the 
         three-dimensional fields in Dyn
       the fields at time t - delta_t are assumed to be in (:,:,previous)
       the fields at time t           are assumed to be in (:,:,current)
       the fields at time t + delta_t are placed        in (:,:,future)
          overwriting whatever is already there

     delta_t = time step in seconds
     
     updates dynamical fields by one time step

   subroutine barotropic_dynamics_end(Dyn, previous, current)
   
     type(dynamics_type), intent(inout)  :: Dyn
     integer, intent(in) :: previous, current
   
      
    Terminates module; 
     writes restart file to 'RESTART/barotropic_dynamics.res'

NAMELIST

&barotropic_dynamics_nml

  integer :: num_lat            = 128
        number of latitudes in global grid

  integer :: num_lon            = 256
        number of longitudes in global grid
        should equal 2*num_lat for Triangular truncation

  integer :: num_fourier        = 85
        the retained fourier wavenumber are n*fourier_inc, where
        n ranges from 0 to num_fourier

  integer :: num_spherical      = 86
        the maximum number of meridional modes for any zonal wavenumber
        for triangular truncation, set num_spherical = num_fourier +1

  integer :: fourier_inc        = 1
        creates a "sector" model if fourier_inc > 1; integration domain is
    (360 degrees longitude)/fourier_inc

  (the default values for num_lat, num_lon, num_fourier and num_spherical define a standard T85 resolution)

  logical :: check_fourier_imag = .false.
        if true, checks to see if fields to be transformed to grid
    domain have zero imaginary part to their zonally symmetric
    modes; useful for debugging

  logical :: south_to_north     = .true.
        true => grid runs from south to north
    false => grid runs from north to south

  logical :: triangular_trunc   = .true.
        true  => shape of truncation is triangular
    false => shape of truncation is rhomboidal

  real    :: robert_coeff       = 0.04
        x(current) => (1-2r)*x(current) + r*(x(future)+x(previous))
    where r = robert_coeff (non-dimensional)

  real    :: longitude_origin   = 0.0
        longitude of first longitude, in degrees
    (if you want the western boundary of first grid box to be at 
         0.0, set longitude_origin = 0.5*360./float(num_lon))

  character :: damping_option   = 'resolution_dependent'
  integer :: damping_order      = 4
  real    :: damping_coeff      = 1.e-04

        damping = nu*(del^2)^n where n = damping order
    damping_option = 'resolution_dependent' or 'resolution_independent'

      = 'resolution_dependent' => nu is set so that the damping rate for the 
            mode (m=0,n=num_spherical-1) equals damping_coeff (in 1/s)
            For triangular truncation, damping_coeff is then the 
            rate of damping of the highest retained mode

      = 'resolution_independent' => nu = damping_coef

  real     :: zeta_0     = 8.e-05 (1/sec)
  integer  :: m_0        = 4
  real     :: eddy_width = 15.0 (degrees longitude)
  real     :: eddy_lat   = 45.0 (degrees latitude)

         eddy component of the initial condition is sinusoidal with
     wavenumber m_0 and with a gaussian distribution of
     vorticity in latitude, centered at eddy_lat with half-width

  logical :: spec_tracer      = .true.
  logical :: grid_tracer      = .true.

       spec_tracer = true => a passive tracer is carried that is advected
                     spectrally, with the same algorithm as the vorticity

       grid_tracer = true => a passive tracer is carried that is advected
                 on the spectral transform grid by a finite-volume algorithm

Both tracers can be carried simultaeneously.
The vorticity and the tracers are initialized within subroutine barotropic_dynamics_init

  real, dimension(2) :: valid_range_v = -1000., 1000.

        A valid range for meridional wind. Model terminates if meridional wind
        goes outside the valid range. Allows model to terminate gracefully when,
        for example, the model becomes numerically unstable.

  character :: initial_zonal_wind = 'two_jets'

       initial_zonal_wind = 'two_jets' => A jet in each hemisphere centered near 30 deg latitude
       initial_zonal_wind = 'zero'     => Zero zonal wind


ERROR MESSAGES

   "Dynamics has not been initialized"
      -- barotropic_dynamics_init must be called before any other
         routines in the module are called
     
   "restart does not exist" 
      -- Time is not equal to Time_init at initalization, but the file
          'INPUT/barotropic_dynamics.res' does not exit