LCOV - code coverage report
Current view: | top level - physics/cam - gw_diffusion.F90 (source / functions) |  |  | Hit | Total | Coverage
Test: | coverage.info |  | Lines: | 27 | 31 | 87.1 %
Date: | 2024-12-17 17:57:11 |  | Functions: | 2 | 2 | 100.0 %

              Line data    Source code
           1             : module gw_diffusion
           2             :
           3             : !
           4             : ! This module contains code computing the effective diffusion of
           5             : ! constituents and dry static energy due to gravity wave breaking.
           6             : !
           7             :
           8             : use gw_utils, only: r8
           9             : use linear_1d_operators, only: TriDiagDecomp
          10             :
          11             : implicit none
          12             : private
          13             : save
          14             :
          15             : public :: gw_ediff
          16             : public :: gw_diff_tend
          17             :
          18             : contains
          19             :
          20             : !==========================================================================
          21             :
          22    19359288 : subroutine gw_ediff(ncol, pver, ngwv, kbot, ktop, tend_level, &
          23    38718576 :      gwut, ubm, nm, rho, dt, prndl, gravit, p, c, vramp, &
          24           0 :      egwdffi, decomp, ro_adjust)
          25             : !
          26             : ! Calculate effective diffusivity associated with GW forcing.
          27             : !
          28             : ! Author: F. Sassi, Jan 31, 2001
          29             : !
          30             :   use gw_utils, only: midpoint_interp
          31             :   use coords_1d, only: Coords1D
          32             :   use vdiff_lu_solver, only: fin_vol_lu_decomp
          33             :
          34             : !-------------------------------Input Arguments----------------------------
          35             :
          36             :   ! Column, level, and gravity wave spectrum dimensions.
          37             :   integer, intent(in) :: ncol, pver, ngwv
          38             :   ! Bottom and top levels to operate on.
          39             :   integer, intent(in) :: kbot, ktop
          40             :   ! Per-column bottom index where tendencies are applied.
          41             :   integer, intent(in) :: tend_level(ncol)
          42             :   ! GW zonal wind tendencies at midpoint.
          43             :   real(r8), intent(in) :: gwut(ncol,pver,-ngwv:ngwv)
          44             :   ! Projection of wind at midpoints.
          45             :   real(r8), intent(in) :: ubm(ncol,pver)
          46             :   ! Brunt-Vaisalla frequency.
          47             :   real(r8), intent(in) :: nm(ncol,pver)
          48             :
          49             :   ! Density at interfaces.
          50             :   real(r8), intent(in) :: rho(ncol,pver+1)
          51             :   ! Time step.
          52             :   real(r8), intent(in) :: dt
          53             :   ! Inverse Prandtl number.
          54             :   real(r8), intent(in) :: prndl
          55             :   ! Acceleration due to gravity.
          56             :   real(r8), intent(in) :: gravit
          57             :   ! Pressure coordinates.
          58             :   type(Coords1D), intent(in) :: p
          59             :   ! Wave phase speeds for each column.
          60             :   real(r8), intent(in) :: c(ncol,-ngwv:ngwv)
          61             :
          62             :   ! Coefficient to ramp down diffusion coeff.
          63             :   real(r8), pointer, intent(in) :: vramp(:)
          64             :
          65             :   ! Adjustment parameter for IGWs.
          66             :   real(r8), intent(in), optional :: &
          67             :        ro_adjust(ncol,-ngwv:ngwv,pver+1)
          68             :
          69             : !-----------------------------Output Arguments-----------------------------
          70             :   ! Effective gw diffusivity at interfaces.
          71             :   real(r8), intent(out) :: egwdffi(ncol,pver+1)
          72             :   ! LU decomposition.
          73             :   type(TriDiagDecomp), intent(out) :: decomp
          74             :
          75             : !-----------------------------Local Workspace------------------------------
          76             :
          77             :   ! Effective gw diffusivity at midpoints.
          78    38718576 :   real(r8) :: egwdffm(ncol,pver)
          79             :   ! Temporary used to hold gw_diffusivity for one level and wavenumber.
          80    38718576 :   real(r8) :: egwdff_lev(ncol)
          81             :   ! (dp/dz)^2 == (gravit*rho)^2
          82    19359288 :   real(r8) :: dpidz_sq(ncol,pver+1)
          83             :   ! Level and wave indices.
          84             :   integer :: k, l
          85             :
          86             :   ! Density scale height.
          87             :   real(r8), parameter :: dscale=7000._r8
          88             :
          89             : !--------------------------------------------------------------------------
          90             :
          91 30405337560 :   egwdffi = 0._r8
          92 30082082472 :   egwdffm = 0._r8
          93             :
          94             :   ! Calculate effective diffusivity at midpoints.
          95   229333104 :   do l = -ngwv, ngwv
          96 13430050491 :      do k = ktop, kbot
          97             :
          98             :         egwdff_lev = &
          99 >22048*10^7 :              prndl * 0.5_r8 * gwut(:,k,l) * (c(:,l)-ubm(:,k)) / nm(:,k)**2
         100             :
         101             :         ! IGWs need ro_adjust factor.
         102 13200717387 :         if (present(ro_adjust)) then
         103           0 :            egwdff_lev = egwdff_lev * ro_adjust(:,l,k)**2
         104             :         end if
         105             :
         106 >22069*10^7 :         egwdffm(:,k) = egwdffm(:,k) + egwdff_lev
         107             :
         108             :      end do
         109             :   end do
         110             :
         111    19359288 :   if (associated(vramp)) then
         112           0 :      do k = ktop,kbot
         113           0 :         egwdffm(:,k) = egwdffm(:,k) * vramp(k)
         114             :      end do
         115             :   endif
         116             :
         117             :   ! Interpolate effective diffusivity to interfaces.
         118             :   ! Assume zero at top and bottom interfaces.
         119    19359288 :   egwdffi(:,ktop+1:kbot) = midpoint_interp(egwdffm(:,ktop:kbot))
         120             :
         121             :   ! Do not calculate diffusivities below level where tendencies are
         122             :   ! actually allowed.
         123  1686948747 :   do k = ktop+1, kbot
         124 27865189158 :      where (k > tend_level) egwdffi(:,k) = 0.0_r8
         125             :   enddo
         126             :
         127             :   ! Calculate (dp/dz)^2.
         128 30405337560 :   dpidz_sq = rho*gravit
         129 30405337560 :   dpidz_sq = dpidz_sq*dpidz_sq
         130             :
         131             :   ! Decompose the diffusion matrix.
         132             :   decomp = fin_vol_lu_decomp(dt, p%section([1,ncol],[ktop,kbot]), &
         133 28589136486 :        coef_q_diff=egwdffi(:,ktop:kbot+1)*dpidz_sq(:,ktop:kbot+1))
         134             :
         135    19359288 : end subroutine gw_ediff
         136             :
         137             : !==========================================================================
         138             :
         139   232311456 : subroutine gw_diff_tend(ncol, pver, kbot, ktop, q, dt, decomp, dq)
         140             :
         141             : !
         142             : ! Calculates tendencies from effective diffusion due to gravity wave
         143             : ! breaking.
         144             : !
         145             : ! Method:
         146             : ! A constituent flux on interfaces is given by:
         147             : !
         148             : !              rho * (w'q') = rho * Deff qz
         149             : !
         150             : ! where (all evaluated on interfaces):
         151             : !
         152             : !        rho   = density
         153             : !        qz    = constituent vertical gradient
         154             : !        Deff  = effective diffusivity
         155             : !
         156             : ! An effective diffusivity is calculated by adding up the diffusivities
         157             : ! from all waves (see gw_ediff). The tendency is calculated by invoking LU
         158             : ! decomposition and solving as for a regular diffusion equation.
         159             : !
         160             : ! Author: Sassi - Jan 2001
         161             : !--------------------------------------------------------------------------
         162             :
         163             : !---------------------------Input Arguments--------------------------------
         164             :
         165             :   ! Column and level dimensions.
         166             :   integer, intent(in) :: ncol, pver
         167             :   ! Bottom and top levels to operate on.
         168             :   integer, intent(in) :: kbot, ktop
         169             :
         170             :   ! Constituent to diffuse.
         171             :   real(r8), intent(in) :: q(ncol,pver)
         172             :   ! Time step.
         173             :   real(r8), intent(in) :: dt
         174             :
         175             :   ! LU decomposition.
         176             :   type(TriDiagDecomp), intent(in) :: decomp
         177             :
         178             : !--------------------------Output Arguments--------------------------------
         179             :
         180             :   ! Constituent tendencies.
         181             :   real(r8), intent(out) :: dq(ncol,pver)
         182             :
         183             : !--------------------------Local Workspace---------------------------------
         184             :
         185             :   ! Temporary storage for constituent.
         186   232311456 :   real(r8) :: qnew(ncol,pver)
         187             :
         188             : !--------------------------------------------------------------------------
         189             :
         190 >36098*10^7 :   dq   = 0.0_r8
         191 >36098*10^7 :   qnew = q
         192             :
         193   232311456 :   call decomp%left_div(qnew(:,ktop:kbot))
         194             :
         195             :   ! Evaluate tendency to be reported back.
         196 >36098*10^7 :   dq = (qnew-q) / dt
         197             :
         198   232311456 : end subroutine gw_diff_tend
         199             :
         200             : end module gw_diffusion

Generated by: LCOV version 1.14