LCOV - code coverage report
Current view: | top level - physics/cam - gw_front.F90 (source / functions) |  |  | Hit | Total | Coverage
Test: | coverage.info |  | Lines: | 33 | 42 | 78.6 %
Date: | 2024-12-17 17:57:11 |  | Functions: | 2 | 5 | 40.0 %

              Line data    Source code
           1             : module gw_front
           2             :
           3             : !
           4             : ! This module handles gravity waves from frontal sources, and was extracted
           5             : ! from gw_drag in May 2013.
           6             : !
           7             :
           8             : use gw_utils, only: r8
           9             : use gw_common, only: GWBand, pi, pver
          10             :
          11             : implicit none
          12             : private
          13             : save
          14             :
          15             : public :: CMSourceDesc
          16             : public :: flat_cm_desc
          17             : public :: gaussian_cm_desc
          18             : public :: gw_cm_src
          19             :
          20             : ! Tuneable settings.
          21             : type CMSourceDesc
          22             :    ! Source level.
          23             :    integer :: ksrc
          24             :    ! Level at which to check whether the frontogenesis function is above
          25             :    ! the critical threshold.
          26             :    integer :: kfront
          27             :    ! Frontogenesis function critical threshold.
          28             :    real(r8) :: frontgfc = huge(1._r8)
          29             :    ! The stress spectrum to produce at the source level.
          30             :    real(r8), allocatable :: src_tau(:)
          31             : end type CMSourceDesc
          32             :
          33             : contains
          34             :
          35             : ! Create a flat profile to be launched (all wavenumbers have the same
          36             : ! source strength, except that l=0 is excluded).
          37           0 : function flat_cm_desc(band, ksrc, kfront, frontgfc, taubgnd) result(desc)
          38             :   ! Wavelengths triggered by frontogenesis.
          39             :   type(GWBand), intent(in) :: band
          40             :   ! The following are used to set the corresponding object components.
          41             :   integer, intent(in) :: ksrc
          42             :   integer, intent(in) :: kfront
          43             :   real(r8), intent(in) :: frontgfc
          44             :   ! Amount of stress to launch from each wavelength.
          45             :   real(r8), intent(in) :: taubgnd
          46             :
          47             :   type(CMSourceDesc) :: desc
          48             :
          49           0 :   desc%ksrc = ksrc
          50           0 :   desc%kfront = kfront
          51           0 :   desc%frontgfc = frontgfc
          52             :
          53           0 :   allocate(desc%src_tau(-band%ngwv:band%ngwv))
          54           0 :   desc%src_tau = taubgnd
          55             :
          56             :   ! Prohibit wavenumber 0.
          57           0 :   desc%src_tau(0) = 0._r8
          58             :
          59           0 : end function flat_cm_desc
          60             :
          61             : ! Create a source tau profile that is a gaussian over wavenumbers (l=0 is
          62             : ! excluded).
          63        1536 : function gaussian_cm_desc(band, ksrc, kfront, frontgfc, height, width) &
          64             :      result(desc)
          65             :
          66             :   use shr_spfn_mod, only: erfc => shr_spfn_erfc
          67             :
          68             :   ! Wavelengths triggered by frontogenesis.
          69             :   type(GWBand), intent(in) :: band
          70             :   ! The following are used to set the corresponding object components.
          71             :   integer, intent(in) :: ksrc
          72             :   integer, intent(in) :: kfront
          73             :   real(r8), intent(in) :: frontgfc
          74             :   ! Parameters of gaussian.
          75             :   real(r8), intent(in) :: height
          76             :   real(r8), intent(in) :: width
          77             :
          78             :   type(CMSourceDesc) :: desc
          79             :
          80             :   ! Bounds used to average bins of the gaussian.
          81        3072 :   real(r8) :: gaussian_bounds(2*band%ngwv+2)
          82             :
          83             :   ! Wavenumber index.
          84             :   integer :: l
          85             :
          86        1536 :   desc%ksrc = ksrc
          87        1536 :   desc%kfront = kfront
          88        1536 :   desc%frontgfc = frontgfc
          89             :
          90        4608 :   allocate(desc%src_tau(-band%ngwv:band%ngwv))
          91             :
          92             :   ! Find the boundaries of each bin.
          93      101376 :   gaussian_bounds(:2*band%ngwv+1) = band%cref - 0.5_r8*band%dc
          94        1536 :   gaussian_bounds(2*band%ngwv+2) = band%cref(band%ngwv) + 0.5_r8*band%dc
          95             :
          96             :   ! Integral of the gaussian at bin interfaces (from the point to
          97             :   ! positive infinity).
          98             :   gaussian_bounds = &
          99      101376 :        [( erfc(gaussian_bounds(l)/width)*height*width*sqrt(pi)/2._r8, &
         100      204288 :        l = 1, 2*band%ngwv+2 )]
         101             :
         102             :   ! Get average in each bin using integral from right to left side.
         103             :   desc%src_tau = &
         104      102912 :        (gaussian_bounds(:2*band%ngwv+1) - gaussian_bounds(2:)) / band%dc
         105             :
         106             :   ! Prohibit wavenumber 0.
         107        1536 :   desc%src_tau(0) = 0._r8
         108             :
         109        1536 : end function gaussian_cm_desc
         110             :
         111             : !==========================================================================
         112     2978352 : subroutine gw_cm_src(ncol, band, desc, u, v, frontgf, &
         113     1489176 :      src_level, tend_level, tau, ubm, ubi, xv, yv, c)
         114             :   use gw_utils, only: get_unit_vector, dot_2d, midpoint_interp
         115             :   !-----------------------------------------------------------------------
         116             :   ! Driver for multiple gravity wave drag parameterization.
         117             :   !
         118             :   ! The parameterization is assumed to operate only where water vapor
         119             :   ! concentrations are negligible in determining the density.
         120             :   !-----------------------------------------------------------------------
         121             :
         122             :   !------------------------------Arguments--------------------------------
         123             :   ! Column dimension.
         124             :   integer, intent(in) :: ncol
         125             :
         126             :   ! Wavelengths triggered by frontogenesis.
         127             :   type(GWBand), intent(in) :: band
         128             :
         129             :   ! Bandification of how to produce the gravity wave bandtrum.
         130             :   type(CMSourceDesc), intent(in) :: desc
         131             :
         132             :   ! Midpoint zonal/meridional winds.
         133             :   real(r8), intent(in) :: u(ncol,pver), v(ncol,pver)
         134             :   ! Frontogenesis function.
         135             :   real(r8), intent(in) :: frontgf(:,:)
         136             :
         137             :   ! Indices of top gravity wave source level and lowest level where wind
         138             :   ! tendencies are allowed.
         139             :   integer, intent(out) :: src_level(ncol)
         140             :   integer, intent(out) :: tend_level(ncol)
         141             :
         142             :   ! Wave Reynolds stress.
         143             :   real(r8), intent(out) :: tau(ncol,-band%ngwv:band%ngwv,pver+1)
         144             :   ! Projection of wind at midpoints and interfaces.
         145             :   real(r8), intent(out) :: ubm(ncol,pver), ubi(ncol,pver+1)
         146             :   ! Unit vectors of source wind (zonal and meridional components).
         147             :   real(r8), intent(out) :: xv(ncol), yv(ncol)
         148             :   ! Phase speeds.
         149             :   real(r8), intent(out) :: c(ncol,-band%ngwv:band%ngwv)
         150             :
         151             :   !---------------------------Local Storage-------------------------------
         152             :   ! Column and wavenumber indices.
         153             :   integer :: k, l
         154             :
         155             :   ! Whether or not to launch waves in this column.
         156     2978352 :   logical :: launch_wave(ncol)
         157             :
         158             :   ! Zonal/meridional wind averaged over source region.
         159     2978352 :   real(r8) :: usrc(ncol), vsrc(ncol)
         160             :
         161             :   !------------------------------------------------------------------------
         162             :   ! Determine the source layer wind and unit vectors, then project winds.
         163             :   !------------------------------------------------------------------------
         164             :
         165             :   ! Just use the source level interface values for the source wind speed
         166             :   ! and direction (unit vector).
         167    24865776 :   src_level = desc%ksrc
         168    24865776 :   tend_level = desc%ksrc
         169    24865776 :   usrc = 0.5_r8*(u(:,desc%ksrc+1)+u(:,desc%ksrc))
         170    24865776 :   vsrc = 0.5_r8*(v(:,desc%ksrc+1)+v(:,desc%ksrc))
         171             :
         172             :   ! Get the unit vector components and magnitude at the surface.
         173     1489176 :   call get_unit_vector(usrc, vsrc, xv, yv, ubi(:,desc%ksrc+1))
         174             :
         175             :   ! Project the local wind at midpoints onto the source wind.
         176   102753144 :   do k = 1, desc%ksrc
         177   102753144 :      ubm(:,k) = dot_2d(u(:,k), v(:,k), xv, yv)
         178             :   end do
         179             :
         180             :   ! Compute the interface wind projection by averaging the midpoint winds.
         181             :   ! Use the top level wind at the top interface.
         182    26354952 :   ubi(:,1) = ubm(:,1)
         183             :
         184     1489176 :   ubi(:,2:desc%ksrc) = midpoint_interp(ubm(:,1:desc%ksrc))
         185             :
         186             :   !-----------------------------------------------------------------------
         187             :   ! Gravity wave sources
         188             :   !-----------------------------------------------------------------------
         189             :
         190 >15207*10^7 :   tau = 0._r8
         191             :
         192             :   ! GW generation depends on frontogenesis at specified level (may be below
         193             :   ! actual source level).
         194    24865776 :   launch_wave = (frontgf(:,desc%kfront) > desc%frontgfc)
         195             :
         196    98285616 :   do l = -band%ngwv, band%ngwv
         197  1617764616 :      where (launch_wave)
         198    96796440 :         tau(:,l,desc%ksrc+1) = desc%src_tau(l)
         199             :      end where
         200             :   end do
         201             :
         202             :   ! Set phase speeds as reference speeds plus the wind speed at the source
         203             :   ! level.
         204             :   c = spread(band%cref, 1, ncol) + &
         205  1617764616 :        spread(ubi(:,desc%ksrc+1), 2, 2*band%ngwv+1)
         206             :
         207     1489176 : end subroutine gw_cm_src
         208             :
         209           0 : end module gw_front

Generated by: LCOV version 1.14