wrfv2_fire/phys/module_fr_sfire_model.F

   1 !
   2 !*** Jan Mandel October 2007 email: jmandel@ucar.edu or Jan.Mandel@gmail.com
   3 !
   4
   5 #define DEBUG_OUT
   6
   7 module module_fr_sfire_model
   8
   9 use module_fr_sfire_core
  10 use module_fr_sfire_util
  11 use module_fr_sfire_phys
  12
  13 contains
  14
  15 subroutine sfire_model (                    &
  16     id,                                     & ! unique number for prints and debug
  17     ifun,                                   & ! what to do see below
  18     need_lfn_update,                          & ! if lfn needs to be synced between tiles
  19     num_ignitions,                          & ! number of ignitions before advancing
  20     ifuelread,nfuel_cat0,                   & ! initialize fuel categories
  21     ifds,ifde,jfds,jfde,                    & ! fire domain dims - the whole domain
  22     ifms,ifme,jfms,jfme,                    & ! fire memory dims - how declared
  23     ifps,ifpe,jfps,jfpe,                    & ! patch - nodes owned by this process
  24     ifts,ifte,jfts,jfte,                    & ! fire tile dims  - this thread
  25     time_start,dt,                          & ! time and increment
  26     fdx,fdy,                                & ! fire mesh spacing,
  27     ignition_start_x,ignition_start_y,      & ! ignition - small arrays
  28     ignition_end_x,ignition_end_y,          &
  29     ignition_radius,                        &
  30     ignition_time,                          &
  31     coord_xf,coord_yf,unit_xf,unit_yf,      & ! fire mesh coordinates
  32     zsf,                                    & ! terrain height (for gradient)
  33     vx,vy,                                  & ! input: wind
  34     lfn,lfn_out,tign,fuel_frac,fire_area,   & ! state: level function, ign time, fuel left, area burning
  35     grnhfx,grnqfx,                          & ! output: heat fluxes
  36     nfuel_cat,                              & ! fuel data per point
  37     fuel_time,                              & ! save derived internal data
  38     bbb,betafl,phiwc,r_0,fgip,ischap &
  39 )
  40
  41 ! This subroutine implements the fire spread model.
  42 ! All quantities are on the fire grid. It inputs
  43 ! winds given on the nodes of the fire grid
  44 ! and outputs the heat fluxes on the cells of the fire grid.
  45 ! This subroutine has no knowledge of any atmospheric model.
  46 ! This code was written to conform with the WRF parallelism model, however it
  47 ! does not depend on it. It can be called with domain equal to tile.
  48 ! Wind and height must be given on 1 more node beyond the domain bounds.
  49 ! The subroutine changes only array entries of the arguments in the tile.
  50 ! Upon exit with ifun=2 (time step), lfn_out is to be copied into lfn by the caller.
  51 ! When this subroutine is used on separate tiles that make a domain the value, the
  52 ! it uses lfn on a strip of width 2 from neighboring tiles.
  53 !
  54 ! All computation is done on one tile.
  55 !
  56 ! This subroutine is intended to be called in a loop like
  57 !
  58 !
  59 ! do ifun=1,6 (if initizalize run, otherwise 3,6)
  60 !   start parallel loop over tiles
  61 !       if ifun=1, set z and fuel data
  62 !       if ifun=3, set the wind arrays
  63 !       call sfire_model(....)
  64 !   end parallel loop over tiles
  65 !
  66 !   if need_lfn_update, halo exchange on lfn width 2
  67 !
  68 !   if ifun=0
  69 !       halo exchange on z width 2
  70 !       halo exchange on fuel data width 1
  71 !   endif
  72 !
  73 !   if ifun=3, halo exchange on winds width 2
  74 !
  75 ! enddo
  76
  77 implicit none
  78
  79 !*** arguments
  80
  81 ! control switches
  82 integer, intent(in) :: id
  83 integer, intent(in) :: ifun                 ! 1 = initialize run pass 1
  84                                             ! 2 = initialize run pass 2
  85                                             ! 3 = initialize timestep
  86                                             ! 4 = do one timestep
  87                                             ! 5 = copy timestep output to input
  88                                             ! 6 = compute output fluxes
  89 logical, intent(out)::need_lfn_update       ! if true, halo update on lfn afterwards
  90 ! scalar data
  91 integer, intent(in) :: num_ignitions        ! number of ignition locations/times
  92 integer, intent(in) :: ifuelread,nfuel_cat0 ! for set_fire_params
  93 integer, intent(in) :: ifds,ifde,jfds,jfde,&  ! fire domain bounds
  94         ifps,ifpe,jfps,jfpe                ! patch - nodes owned by this process
  95 integer, intent(in) :: ifts,ifte,jfts,jfte  ! fire tile bounds
  96 integer, intent(in) :: ifms,ifme,jfms,jfme  ! fire memory array bounds
  97 REAL,INTENT(in) :: time_start,dt            ! starting time, time step
  98 REAL,INTENT(in) :: fdx,fdy                  ! spacing of the fire mesh
  99 ! array data
 100 real, dimension(num_ignitions), intent(in):: &
 101     ignition_start_x,ignition_start_y, &
 102     ignition_end_x,ignition_end_y,ignition_radius, & ! start, end, radius, time
 103     ignition_time                           !  of ignition lines
 104 real, dimension(ifms:ifme, jfms:jfme), intent(in):: &
 105     coord_xf,coord_yf                       !  node coordinates
 106 real, intent(in):: unit_xf,unit_yf          !  coordinate units in m
 107 REAL, INTENT(in), dimension(ifms:ifme,jfms:jfme):: &
 108     vx,vy                                   ! wind m/s (node based), data, variable
 109
 110 ! state
 111 REAL, INTENT(inout), dimension(ifms:ifme,jfms:jfme):: &
 112     zsf,   &                                ! terrain height, node based, data, constant after extr
 113     lfn   , &                               ! level function: fire is where lfn<0 (node)
 114     tign  , &                               ! absolute time of ignition (node)
 115     fuel_frac                               ! fuel fraction (node), currently redundant
 116
 117 REAL, INTENT(out), dimension(ifms:ifme,jfms:jfme):: &
 118     fire_area                               ! fraction of each cell burning
 119
 120 ! output
 121 REAL, INTENT(out), dimension(ifms:ifme,jfms:jfme):: &
 122     lfn_out, &                              !
 123     grnhfx,grnqfx                           ! heat fluxes J/m^2/s  (cell)
 124
 125 ! constant arrays - set at initialization
 126 integer, intent(inout), dimension(ifms:ifme, jfms:jfme)::nfuel_cat ! cell based, data, constant
 127 real,intent(inout),dimension(ifms:ifme,jfms:jfme):: fuel_time
 128 real,intent(inout),dimension(ifms:ifme,jfms:jfme):: bbb,betafl,phiwc,r_0 ! (node) spread formula coefficients
 129 real,intent(inout),dimension(ifms:ifme,jfms:jfme):: fgip                 ! (cell) init mass of surface fuel (kg/m^2)
 130 integer,intent(inout),dimension(ifms:ifme,jfms:jfme):: ischap            ! (node) .ne.0 if chapparal
 131
 132 !*** local
 133
 134 integer :: xifms,xifme,xjfms,xjfme  ! memory bounds for pass-through arguments to normal spread
 135 real, dimension(ifts:ifte,jfts:jfte)::fuel_frac_burnt,fuel_frac_end
 136 integer::ignited,ig,i,j
 137 real::tbound
 138 character(len=128)::msg
 139
 140 !*** executable
 141
 142 call check_mesh_2dim(ifts-1,ifte+1,jfts-1,jfte+1,ifms,ifme,jfms,jfme)
 143
 144
 145 xifms=ifms  ! dimensions for the include file
 146 xifme=ifme
 147 xjfms=jfms
 148 xjfme=jfme
 149
 150
 151 need_lfn_update=.false.
 152
 153 if(ifun.eq.1)then       ! do nothing, init pass 1 is outside only
 154 elseif(ifun.eq.2)then
 155         ! initialize all arrays that the model will not change later
 156
 157         ! assuming halo on zsf done
 158         ! extrapolate on 1 row of cells beyond the domain boundary
 159         ! including on the halo regions
 160
 161         call continue_at_boundary(1,1,0., & ! do x direction or y direction
 162             ifms,ifme,jfms,jfme,           &                ! memory dims
 163             ifds,ifde,jfds,jfde, &                     ! domain dims
 164             ifps,ifpe,jfps,jfpe, &            ! patch dims - winds defined up to +1
 165             ifts,ifte,jfts,jfte, &                ! tile dims
 166             zsf)                               ! array
 167
 168         call set_nfuel_cat( &
 169             ifms,ifme,jfms,jfme, &
 170             ifts,ifte,jfts,jfte, &
 171             ifuelread,nfuel_cat0,&
 172             zsf,nfuel_cat)            ! better not use the extrapolated zsf!!
 173
 174         ! uses nfuel_cat to set the other fuel data arrays
 175         ! needs zsf on halo width 1 to compute the terrain gradient
 176         call set_fire_params(   &
 177             ifds,ifde,jfds,jfde, &
 178             ifms,ifme,jfms,jfme, &
 179             ifts,ifte,jfts,jfte, &
 180             fdx,fdy,             &
 181             nfuel_cat,fuel_time &
 182 #           include "fr_sfire_params_args.h"
 183 )
 184
 185         call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,zsf,'model: terrain height')
 186
 187         ! initialize model state to no fire
 188         call init_no_fire   ( &
 189             ifds,ifde,jfds,jfde, &
 190             ifms,ifme,jfms,jfme, &
 191             ifts,ifte,jfts,jfte, &
 192             fdx,fdy,time_start,  &
 193             fuel_frac,fire_area,lfn,tign)
 194
 195         need_lfn_update=.true. ! because we have set lfn
 196
 197 elseif(ifun.eq.3)then   ! ignition if so specified
 198
 199     ! check for ignitions
 200     do ig = 1,num_ignitions
 201
 202         if(ignition_time(ig)>=time_start.and.ignition_time(ig)<time_start+dt)then
 203             call ignite_fire(                             &
 204                 ifds,ifde,jfds,jfde,                      & ! fire domain dims - the whole domain
 205                 ifms,ifme,jfms,jfme,                      &
 206                 ifts,ifte,jfts,jfte,                      &
 207                 ignition_start_x(ig),ignition_start_y(ig),&
 208                 ignition_end_x(ig),ignition_end_y(ig),    &
 209                 ignition_radius(ig),                      &
 210                 ignition_time(ig),                        &
 211                 coord_xf,coord_yf,unit_xf,unit_yf,        &
 212                 lfn,tign,ignited)
 213
 214             need_lfn_update=.true. ! if ignition, lfn changed
 215 #ifdef DEBUG_OUT
 216             call write_array_m(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,lfn,'lfn_ig',id)
 217             call write_array_m(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,coord_xf,'coord_xf_ig',id)
 218             call write_array_m(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,coord_yf,'coord_yf_ig',id)
 219 #endif
 220         endif
 221
 222     enddo
 223
 224 elseif (ifun.eq.4) then  ! do the timestep
 225
 226 !   compute fuel fraction at start
 227 !    call fuel_left( &
 228 !        ifms,ifme,jfms,jfme, &
 229 !        ifts,ifte,jfts,jfte, &
 230 !        ifms,ifme,jfms,jfme, &
 231 !        lfn,tign,fuel_time,time_start,fuel_frac,fire_area) ! fuel frac is shared
 232
 233     call print_2d_stats(ifts,ifte,jfts,jfte, &
 234                    ifms,ifme,jfms,jfme, &
 235                    fuel_frac,'model: fuel_frac start')
 236
 237     ! advance the model from time_start to time_start+dt
 238     ! return the fuel fraction burnt this call in each fire cell
 239     ! will call module_fr_sfire_speed::normal_spread for propagation speed
 240     ! We cannot simply compute the spread rate here because that will change with the
 241     ! angle of the wind and the direction of propagation, thus it is done in subroutine
 242     ! normal_spread at each fire time step. Instead, we pass arguments that
 243     ! the speed function may use. The include is to guarantee this is done consistently
 244     ! over the call chain.
 245
 246 !   propagate level set function in time
 247 !   set lfn_out tign
 248 !   lfn does not change, tign has no halos
 249
 250     call prop_ls(id,     &
 251         ifds,ifde,jfds,jfde,                      & ! fire domain dims - the whole domain
 252         ifms,ifme,jfms,jfme,                      &
 253         ifps,ifpe,jfps,jfpe, &                ! patch - nodes owned by this process
 254         ifts,ifte,jfts,jfte,                      &
 255         time_start,dt,fdx,fdy,tbound,  &
 256         lfn,lfn_out,tign &
 257 #       include "fr_sfire_params_args.h"
 258     )
 259
 260 elseif (ifun.eq.5) then ! copy the result of timestep back to input
 261     ! this cannot be done in the time step itself because of race condition
 262     ! some thread may still be using lfn as input in their tile halo
 263
 264     do j=jfts,jfte
 265         do i=ifts,ifte
 266             lfn(i,j)=lfn_out(i,j)
 267             ! if want to try timestep again treat tign the same way here
 268             ! even if tign does not need a halo
 269         enddo
 270     enddo
 271
 272     need_lfn_update=.true. ! duh
 273
 274 elseif (ifun.eq.6) then ! timestep postprocessing
 275
 276     ! compute the heat fluxes from the fuel burned
 277     ! needs lfn and tign from neighbors so halo must be updated before
 278     call fuel_left(&
 279         ifms,ifme,jfms,jfme, &
 280         ifts,ifte,jfts,jfte, &
 281         ifts,ifte,jfts,jfte, &
 282         lfn,tign,fuel_time,time_start+dt,fuel_frac_end,fire_area) !fuel_frac_end is private and tile based
 283
 284     call print_2d_stats(ifts,ifte,jfts,jfte, &
 285                    ifts,ifte,jfts,jfte, &
 286                    fuel_frac_end,'model: fuel_frac end')
 287
 288     do j=jfts,jfte
 289         do i=ifts,ifte
 290             fuel_frac_burnt(i,j)=fuel_frac(i,j)-fuel_frac_end(i,j) ! fuel lost this timestep
 291             fuel_frac(i,j)=fuel_frac_end(i,j) ! copy new value to state array
 292         enddo
 293     enddo
 294
 295     call print_2d_stats(ifts,ifte,jfts,jfte, &
 296                    ifts,ifte,jfts,jfte, &
 297                    fuel_frac_burnt,'model: fuel_frac burned')
 298
 299     call heat_fluxes(dt,                          &
 300         ifms,ifme,jfms,jfme,                      &
 301         ifts,ifte,jfts,jfte,                      &
 302         ifts,ifte,jfts,jfte,                      &  ! fuel_frac_burned is tile dimensioned
 303         fgip,                                     &
 304         fuel_frac_burnt,                          & !
 305         grnhfx,grnqfx)                              !out
 306
 307
 308     call print_2d_stats(ifts,ifte,jfts,jfte, &
 309                    ifms,ifme,jfms,jfme, &
 310                    grnhfx,'model: heat flux(J/m^2/s)')
 311
 312 else
 313     write(msg,*)'sfire_model: bad ifun=',ifun
 314     call crash(msg)
 315 endif
 316
 317 end subroutine sfire_model
 318
 319 !
 320 !*****************
 321 !
 322
 323 end module module_fr_sfire_model