Line data    Source code

       1             : module gw_oro
       2             : 
       3             : !
       4             : ! This module handles gravity waves from orographic sources, and was
       5             : ! extracted from gw_drag in May 2013.
       6             : !
       7             : use gw_utils, only: r8
       8             : use coords_1d, only: Coords1D
       9             : 
      10             : implicit none
      11             : private
      12             : save
      13             : 
      14             : ! Public interface
      15             : public :: gw_oro_src
      16             : 
      17             : contains
      18             : 
      19             : !==========================================================================
      20             : 
      21     1489176 : subroutine gw_oro_src(ncol, band, p, &
      22     1489176 :      u, v, t, sgh, zm, nm, &
      23     1489176 :      src_level, tend_level, tau, ubm, ubi, xv, yv, c)
      24             :   use gw_common, only: GWBand, pver, rair
      25             :   use gw_utils, only: get_unit_vector, dot_2d, midpoint_interp
      26             :   !-----------------------------------------------------------------------
      27             :   ! Orographic source for multiple gravity wave drag parameterization.
      28             :   !
      29             :   ! The stress is returned for a single wave with c=0, over orography.
      30             :   ! For points where the orographic variance is small (including ocean),
      31             :   ! the returned stress is zero.
      32             :   !------------------------------Arguments--------------------------------
      33             :   ! Column dimension.
      34             :   integer, intent(in) :: ncol
      35             :   ! Band to emit orographic waves in.
      36             :   ! Regardless, we will only ever emit into l = 0.
      37             :   type(GWBand), intent(in) :: band
      38             :   ! Pressure coordinates.
      39             :   type(Coords1D), intent(in) :: p
      40             : 
      41             :   ! Midpoint zonal/meridional winds.
      42             :   real(r8), intent(in) :: u(ncol,pver), v(ncol,pver)
      43             :   ! Midpoint temperatures.
      44             :   real(r8), intent(in) :: t(ncol,pver)
      45             :   ! Standard deviation of orography.
      46             :   real(r8), intent(in) :: sgh(ncol)
      47             :   ! Midpoint altitudes.
      48             :   real(r8), intent(in) :: zm(ncol,pver)
      49             :   ! Midpoint Brunt-Vaisalla frequencies.
      50             :   real(r8), intent(in) :: nm(ncol,pver)
      51             : 
      52             :   ! Indices of top gravity wave source level and lowest level where wind
      53             :   ! tendencies are allowed.
      54             :   integer, intent(out) :: src_level(ncol)
      55             :   integer, intent(out) :: tend_level(ncol)
      56             : 
      57             :   ! Wave Reynolds stress.
      58             :   real(r8), intent(out) :: tau(ncol,-band%ngwv:band%ngwv,pver+1)
      59             :   ! Projection of wind at midpoints and interfaces.
      60             :   real(r8), intent(out) :: ubm(ncol,pver), ubi(ncol,pver+1)
      61             :   ! Unit vectors of source wind (zonal and meridional components).
      62             :   real(r8), intent(out) :: xv(ncol), yv(ncol)
      63             :   ! Phase speeds.
      64             :   real(r8), intent(out) :: c(ncol,-band%ngwv:band%ngwv)
      65             : 
      66             :   !---------------------------Local Storage-------------------------------
      67             :   ! Column and level indices.
      68             :   integer :: i, k
      69             : 
      70             :   ! Surface streamline displacement height (2*sgh).
      71     2978352 :   real(r8) :: hdsp(ncol)
      72             :   ! Max orographic standard deviation to use.
      73             :   real(r8) :: sghmax
      74             :   ! c=0 stress from orography.
      75     2978352 :   real(r8) :: tauoro(ncol)
      76             :   ! Averages over source region.
      77     2978352 :   real(r8) :: nsrc(ncol) ! B-V frequency.
      78     2978352 :   real(r8) :: rsrc(ncol) ! Density.
      79     2978352 :   real(r8) :: usrc(ncol) ! Zonal wind.
      80     2978352 :   real(r8) :: vsrc(ncol) ! Meridional wind.
      81             : 
      82             :   ! Difference in interface pressure across source region.
      83     1489176 :   real(r8) :: dpsrc(ncol)
      84             : 
      85             :   ! Limiters (min/max values)
      86             :   ! min surface displacement height for orographic waves
      87             :   real(r8), parameter :: orohmin = 10._r8
      88             :   ! min wind speed for orographic waves
      89             :   real(r8), parameter :: orovmin = 2._r8
      90             : 
      91             : !--------------------------------------------------------------------------
      92             : ! Average the basic state variables for the wave source over the depth of
      93             : ! the orographic standard deviation. Here we assume that the appropiate
      94             : ! values of wind, stability, etc. for determining the wave source are
      95             : ! averages over the depth of the atmosphere penterated by the typical
      96             : ! mountain.
      97             : ! Reduces to the bottom midpoint values when sgh=0, such as over ocean.
      98             : !--------------------------------------------------------------------------
      99             : 
     100    24865776 :   hdsp = 2.0_r8 * sgh
     101             : 
     102     1489176 :   k = pver
     103    24865776 :   src_level = k-1
     104    24865776 :   rsrc = p%mid(:,k)/(rair*t(:,k)) * p%del(:,k)
     105    24865776 :   usrc = u(:,k) * p%del(:,k)
     106    24865776 :   vsrc = v(:,k) * p%del(:,k)
     107    24865776 :   nsrc = nm(:,k)* p%del(:,k)
     108             : 
     109     1489176 :   do k = pver-1, 1, -1
     110    24865776 :      do i = 1, ncol
     111    24865776 :         if (hdsp(i) > sqrt(zm(i,k)*zm(i,k+1))) then
     112           0 :            src_level(i) = k-1
     113             :            rsrc(i) = rsrc(i) + &
     114           0 :                 p%mid(i,k) / (rair*t(i,k)) * p%del(i,k)
     115           0 :            usrc(i) = usrc(i) + u(i,k) * p%del(i,k)
     116           0 :            vsrc(i) = vsrc(i) + v(i,k) * p%del(i,k)
     117           0 :            nsrc(i) = nsrc(i) + nm(i,k)* p%del(i,k)
     118             :         end if
     119             :      end do
     120             :      ! Break the loop when all source levels found.
     121    24865776 :      if (all(src_level >= k)) exit
     122             :   end do
     123             : 
     124    24865776 :   do i = 1, ncol
     125    24865776 :      dpsrc(i) = p%ifc(i,pver+1) - p%ifc(i,src_level(i)+1)
     126             :   end do
     127             : 
     128    24865776 :   rsrc = rsrc / dpsrc
     129    24865776 :   usrc = usrc / dpsrc
     130    24865776 :   vsrc = vsrc / dpsrc
     131    24865776 :   nsrc = nsrc / dpsrc
     132             : 
     133             :   ! Get the unit vector components and magnitude at the surface.
     134     1489176 :   call get_unit_vector(usrc, vsrc, xv, yv, ubi(:,pver+1))
     135             : 
     136             :   ! Project the local wind at midpoints onto the source wind.
     137    40207752 :   do k = 1, pver
     138    40207752 :      ubm(:,k) = dot_2d(u(:,k), v(:,k), xv, yv)
     139             :   end do
     140             : 
     141             :   ! Compute the interface wind projection by averaging the midpoint winds.
     142             :   ! Use the top level wind at the top interface.
     143    24865776 :   ubi(:,1) = ubm(:,1)
     144             : 
     145     1489176 :   ubi(:,2:pver) = midpoint_interp(ubm)
     146             : 
     147             :   ! Determine the orographic c=0 source term following McFarlane (1987).
     148             :   ! Set the source top interface index to pver, if the orographic term is
     149             :   ! zero.
     150    24865776 :   do i = 1, ncol
     151    24865776 :      if ((ubi(i,pver+1) > orovmin) .and. (hdsp(i) > orohmin)) then
     152           0 :         sghmax = band%fcrit2 * (ubi(i,pver+1) / nsrc(i))**2
     153             :         tauoro(i) = 0.5_r8 * band%kwv * min(hdsp(i)**2, sghmax) * &
     154           0 :              rsrc(i) * nsrc(i) * ubi(i,pver+1)
     155             :      else
     156    23376600 :         tauoro(i) = 0._r8
     157    23376600 :         src_level(i) = pver
     158             :      end if
     159             :   end do
     160             : 
     161             :   ! Set the phase speeds and wave numbers in the direction of the source
     162             :   ! wind. Set the source stress magnitude (positive only, note that the
     163             :   ! sign of the stress is the same as (c-u).
     164   713072880 :   tau = 0._r8
     165    26354952 :   do k = pver, minval(src_level), -1
     166    26354952 :      where (src_level <= k) tau(:,0,k+1) = tauoro
     167             :   end do
     168             : 
     169             :   ! Allow wind tendencies all the way to the model bottom.
     170    24865776 :   tend_level = pver
     171             : 
     172             :   ! No spectrum; phase speed is just 0.
     173    26354952 :   c = 0._r8
     174             : 
     175     1489176 : end subroutine gw_oro_src
     176             : 
     177             : end module gw_oro