LCOV - code coverage report
Current view: | top level - physics/cam - gw_oro.F90 (source / functions) |  |  | Hit | Total | Coverage
Test: | coverage.info |  | Lines: | 42 | 49 | 85.7 %
Date: | 2025-01-13 21:54:50 |  | Functions: | 1 | 1 | 100.0 %

              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

Generated by: LCOV version 1.14