wrfv2_fire/phys/module_fr_sfire_atm.F

   1 !WRF:MEDIATION_LAYER:FIRE_MODEL
   2 ! Routines dealing with the atmosphere
   3
   4 module module_fr_sfire_atm
   5
   6 use module_model_constants, only: cp,xlv
   7 use module_fr_sfire_util
   8
   9 contains
  10
  11 SUBROUTINE fire_tendency( &
  12     ids,ide, kds,kde, jds,jde,   & ! dimensions
  13     ims,ime, kms,kme, jms,jme,   &
  14     its,ite, kts,kte, jts,jte,   &
  15     grnhfx,grnqfx,canhfx,canqfx, & ! heat fluxes summed up to  atm grid
  16     alfg,alfc,z1can,             & ! coeffients, properties, geometry
  17     zs,z_at_w,dz8w,mu,rho,       &
  18     rthfrten,rqvfrten)             ! theta and Qv tendencies
  19
  20 ! This routine is atmospheric physics
  21 ! it does NOT go into module_fr_sfire_phys because it is not related to fire physical processes
  22
  23 ! --- this routine takes fire generated heat and moisture fluxes and
  24 !     calculates their influence on the theta and water vapor
  25 ! --- note that these tendencies are valid at the Arakawa-A location
  26
  27    IMPLICIT NONE
  28
  29 ! --- incoming variables
  30
  31    INTEGER , INTENT(in) :: ids,ide, kds,kde, jds,jde, &
  32                            ims,ime, kms,kme, jms,jme, &
  33                            its,ite, kts,kte, jts,jte
  34
  35    REAL, INTENT(in), DIMENSION( ims:ime,jms:jme ) :: grnhfx,grnqfx  ! W/m^2
  36    REAL, INTENT(in), DIMENSION( ims:ime,jms:jme ) :: canhfx,canqfx  ! W/m^2
  37    REAL, INTENT(in), DIMENSION( ims:ime,jms:jme ) :: zs  ! topography (m abv sealvl)
  38    REAL, INTENT(in), DIMENSION( ims:ime,jms:jme ) :: mu  ! dry air mass (Pa)
  39
  40    REAL, INTENT(in), DIMENSION( ims:ime,kms:kme,jms:jme ) :: z_at_w ! m abv sealvl
  41    REAL, INTENT(in), DIMENSION( ims:ime,kms:kme,jms:jme ) :: dz8w   ! dz across w-lvl
  42    REAL, INTENT(in), DIMENSION( ims:ime,kms:kme,jms:jme ) :: rho    ! density
  43
  44    REAL, INTENT(in) :: alfg ! extinction depth surface fire heat (m)
  45    REAL, INTENT(in) :: alfc ! extinction depth crown  fire heat (m)
  46    REAL, INTENT(in) :: z1can    ! height of crown fire heat release (m)
  47
  48 ! --- outgoing variables
  49
  50    REAL, INTENT(out), DIMENSION( ims:ime,kms:kme,jms:jme ) ::   &
  51        rthfrten, & ! theta tendency from fire (in mass units)
  52        rqvfrten    ! Qv tendency from fire (in mass units)
  53 ! --- local variables
  54
  55    INTEGER :: i,j,k
  56    INTEGER :: i_st,i_en, j_st,j_en, k_st,k_en
  57
  58    REAL :: cp_i
  59    REAL :: rho_i
  60    REAL :: xlv_i
  61    REAL :: z_w
  62    REAL :: fact_g, fact_c
  63    REAL :: alfg_i, alfc_i
  64
  65    REAL, DIMENSION( its:ite,kts:kte,jts:jte ) :: hfx,qfx
  66
  67 !!   character(len=128)::msg
  68
  69         do j=jts,jte
  70             do k=kts,min(kte+1,kde)
  71                do i=its,ite
  72                    rthfrten(i,k,j)=0.
  73                    rqvfrten(i,k,j)=0.
  74                enddo
  75             enddo
  76         enddo
  77
  78
  79 ! --- set some local constants
  80
  81
  82    cp_i = 1./cp     ! inverse of specific heat
  83    xlv_i = 1./xlv   ! inverse of latent heat
  84    alfg_i = 1./alfg
  85    alfc_i = 1./alfc
  86
  87 !!write(msg,'(8e11.3)')cp,cp_i,xlv,xlv_i,alfg,alfc,z1can
  88 !!call message(msg)
  89
  90    call print_2d_stats(its,ite,jts,jte,ims,ime,jms,jme,grnhfx,'fire_tendency:grnhfx')
  91    call print_2d_stats(its,ite,jts,jte,ims,ime,jms,jme,grnqfx,'fire_tendency:grnqfx')
  92
  93 ! --- set loop indicies : note that
  94
  95    i_st = MAX(its,ids+1)
  96    i_en = MIN(ite,ide-1)
  97    k_st = kts
  98    k_en = MIN(kte,kde-1)
  99    j_st = MAX(jts,jds+1)
 100    j_en = MIN(jte,jde-1)
 101
 102 ! --- distribute fluxes
 103
 104    DO j = j_st,j_en
 105       DO k = k_st,k_en
 106          DO i = i_st,i_en
 107
 108             ! --- set z (in meters above ground)
 109
 110             z_w = z_at_w(i,k,j) - zs(i,j) ! should be zero when k=k_st
 111
 112             ! --- heat flux
 113
 114             fact_g = cp_i * EXP( - alfg_i * z_w )
 115             IF ( z_w < z1can ) THEN
 116                fact_c = cp_i
 117             ELSE
 118                fact_c = cp_i * EXP( - alfc_i * (z_w - z1can) )
 119             END IF
 120             hfx(i,k,j) = fact_g * grnhfx(i,j) + fact_c * canhfx(i,j)
 121
 122 !!            write(msg,2)i,k,j,z_w,grnhfx(i,j),hfx(i,k,j)
 123 !!2           format('hfx:',3i4,6e11.3)
 124 !!            call message(msg)
 125
 126             ! --- vapor flux
 127
 128             fact_g = xlv_i * EXP( - alfg_i * z_w )
 129             IF (z_w < z1can) THEN
 130                fact_c = xlv_i
 131             ELSE
 132                fact_c = xlv_i * EXP( - alfc_i * (z_w - z1can) )
 133             END IF
 134             qfx(i,k,j) = fact_g * grnqfx(i,j) + fact_c * canqfx(i,j)
 135
 136 !!            if(hfx(i,k,j).ne.0. .or. qfx(i,k,j) .ne. 0.)then
 137 !!                write(msg,1)i,k,j,hfx(i,k,j),qfx(i,k,j)
 138 !!1               format('tend:',3i6,2e11.3)
 139 !!                call message(msg)
 140 !            endif
 141
 142          END DO
 143       END DO
 144    END DO
 145
 146 ! --- add flux divergence to tendencies
 147 !
 148 !   multiply by dry air mass (mu) to eliminate the need to
 149 !   call sr. calculate_phy_tend (in dyn_em/module_em.F)
 150
 151    DO j = j_st,j_en
 152       DO k = k_st,k_en-1
 153          DO i = i_st,i_en
 154
 155             rho_i = 1./rho(i,k,j)
 156
 157             rthfrten(i,k,j) = - mu(i,j) * rho_i * (hfx(i,k+1,j)-hfx(i,k,j)) / dz8w(i,k,j)
 158             rqvfrten(i,k,j) = - mu(i,j) * rho_i * (qfx(i,k+1,j)-qfx(i,k,j)) / dz8w(i,k,j)
 159
 160          END DO
 161       END DO
 162    END DO
 163
 164    call print_3d_stats(its,ite,kts,kte,jts,jte,ims,ime,kms,kme,jms,jme,rthfrten,'fire_tendency:rthfrten')
 165    call print_3d_stats(its,ite,kts,kte,jts,jte,ims,ime,kms,kme,jms,jme,rqvfrten,'fire_tendency:rqvfrten')
 166
 167    RETURN
 168
 169 END SUBROUTINE fire_tendency
 170
 171 !
 172 !***
 173 !
 174
 175 end module module_fr_sfire_atm