phys/module_fr_fire_model.F

   1 !
   2 #define DEBUG_OUT
   3
   4 module module_fr_fire_model
   5
   6 use module_fr_fire_core
   7 use module_fr_fire_util
   8 use module_fr_fire_phys
   9
  10 contains
  11
  12 subroutine fire_model (                    &
  13     id,                                     & ! unique number for prints and debug
  14     ifun,                                   & ! what to do see below
  15     restart,                                &
  16     need_lfn_update,                          & ! if lfn needs to be synced between tiles
  17     run_fuel_moisture,                      & ! if need update fuel moisture in pass 4
  18     num_ignitions,                          & ! number of ignitions before advancing
  19     ifuelread,nfuel_cat0,                   & ! initialize fuel categories
  20     ifds,ifde,jfds,jfde,                    & ! fire domain dims - the whole domain
  21     ifms,ifme,jfms,jfme,                    & ! fire memory dims - how declared
  22     ifps,ifpe,jfps,jfpe,                    & ! patch - nodes owned by this process
  23     ifts,ifte,jfts,jfte,                    & ! fire tile dims  - this thread
  24     time_start,dt,                          & ! time and increment
  25     fdx,fdy,                                & ! fire mesh spacing,
  26     ignition_line,                          & ! small array of ignition line descriptions
  27     ignitions_done,ignited_tile,            &
  28     coord_xf,coord_yf,unit_xf,unit_yf,      & ! fire mesh coordinates
  29     lfn,                                & ! state: level function
  30     lfn_hist,                           & ! PAJ: to init obs fire perimeter.
  31     fire_is_real_perim,                 & ! PAJ: to init obs fire perimeter.
  32     lfn_0,lfn_1,lfn_2,lfn_s0,lfn_s1,lfn_s2,lfn_s3,flame_length,ros_front, & ! state
  33     lfn_out,tign,fuel_frac,fire_area,   & ! state: level function, ign time, fuel left, area burning
  34     burnt_area_dt,                          &
  35     grnhfx,grnqfx,                          & ! output: heat fluxes
  36     ros,                                    & ! output: rate of spread
  37     nfuel_cat,                              & ! fuel data per point
  38     fuel_time,                              & ! save derived internal data
  39     fp, &
  40     grid, &
  41     ids_hu,ide_hu,jds_hu,jde_hu,kds_hu,kde_hu, &
  42     ims_hu,ime_hu,jms_hu,jme_hu,kms_hu,kme_hu, &
  43     ips_hu,ipe_hu,jps_hu,jpe_hu,kps_hu,kpe_hu  &
  44 )
  45
  46 ! This subroutine implements the fire spread model.
  47 ! All quantities are on the fire grid. It inputs
  48 ! winds given on the nodes of the fire grid
  49 ! and outputs the heat fluxes on the cells of the fire grid.
  50 ! This subroutine has no knowledge of any atmospheric model.
  51 ! This code was written to conform with the WRF parallelism model, however it
  52 ! does not depend on it. It can be called with domain equal to tile.
  53 ! Wind and height must be given on 1 more node beyond the domain bounds.
  54 ! The subroutine changes only array entries of the arguments in the tile.
  55 ! Upon exit with ifun=2 (time step), lfn_out is to be copied into lfn by the caller.
  56 ! When this subroutine is used on separate tiles that make a domain the value, the
  57 ! it uses lfn on a strip of width 2 from neighboring tiles.
  58 !
  59 ! All computation is done on one tile.
  60 !
  61 ! This subroutine is intended to be called in a loop like
  62 !
  63 !
  64 ! do ifun=1,6 (if initizalize run, otherwise 3,6)
  65 !   start parallel loop over tiles
  66 !       if ifun=1, set z and fuel data
  67 !       if ifun=3, set the wind arrays
  68 !       call fire_model(....)
  69 !   end parallel loop over tiles
  70 !
  71 !   if need_lfn_update, halo exchange on lfn width 2
  72 !
  73 !   if ifun=0
  74 !       halo exchange on z width 2
  75 !       halo exchange on fuel data width 1
  76 !   endif
  77 !
  78 !   if ifun=3, halo exchange on winds width 2
  79 !
  80 ! enddo
  81
  82 USE module_domain , only: domain
  83 #ifdef DM_PARALLEL
  84     USE module_dm        , ONLY : ntasks_x,ntasks_y,local_communicator,mytask,ntasks
  85     USE module_comm_dm , ONLY : halo_fire_lfn_sub
  86 #endif
  87 USE module_configure, only: grid_config_rec_type
  88
  89 implicit none
  90
  91 !*** arguments
  92
  93 ! DME added for halo update
  94     type(domain) , target :: grid
  95     integer, intent(in):: ids_hu,ide_hu,jds_hu,jde_hu,kds_hu,kde_hu, &
  96                           ims_hu,ime_hu,jms_hu,jme_hu,kms_hu,kme_hu, &
  97                           ips_hu,ipe_hu,jps_hu,jpe_hu,kps_hu,kpe_hu
  98 !
  99
 100 ! control switches
 101 integer, intent(in) :: id
 102 integer, intent(in) :: ifun                 ! 1 = initialize run pass 1
 103                                             ! 2 = initialize run pass 2
 104                                             ! 3 = initialize timestep
 105                                             ! 4 = do one timestep
 106                                             ! 5 = copy timestep output to input
 107                                             ! 6 = compute output fluxes
 108 logical, intent(in):: restart               ! if true, use existing state
 109 logical, intent(out)::need_lfn_update       ! if true, halo update on lfn afterwards
 110 logical, intent(in)::run_fuel_moisture      !
 111 ! scalar data
 112 integer, intent(in) :: num_ignitions        ! number of ignition lines
 113 integer, intent(in) :: ifuelread,nfuel_cat0 ! for set_fire_params
 114 integer, intent(in) :: ifds,ifde,jfds,jfde,&  ! fire domain bounds
 115         ifps,ifpe,jfps,jfpe                ! patch - nodes owned by this process
 116 integer, intent(in) :: ifts,ifte,jfts,jfte  ! fire tile bounds
 117 integer, intent(in) :: ifms,ifme,jfms,jfme  ! fire memory array bounds
 118 REAL,INTENT(in) :: time_start,dt            ! starting time, time step
 119 REAL,INTENT(in) :: fdx,fdy                  ! spacing of the fire mesh
 120 ! array data
 121 type(ignition_line_type), dimension (num_ignitions), intent(in):: ignition_line ! descriptions of ignition lines
 122 integer, intent(out):: ignited_tile(num_ignitions),ignitions_done
 123 real, dimension(ifms:ifme, jfms:jfme), intent(in):: &
 124     coord_xf,coord_yf                       !  node coordinates
 125 real, intent(in):: unit_xf,unit_yf          !  coordinate units in m
 126
 127 ! state
 128     ! PAJ: to init obs fire perimeter
 129   real, intent(in), dimension(ifms:ifme,jfms:jfme):: lfn_hist
 130   logical, intent(in) :: fire_is_real_perim
 131
 132 REAL, INTENT(inout), dimension(ifms:ifme,jfms:jfme):: &
 133     lfn   , &                               ! level function: fire is where lfn<0 (node)
 134     tign  , &                               ! absolute time of ignition (node)
 135     fuel_frac                               ! fuel fraction (node), currently redundant
 136
 137 REAL, INTENT(inout), dimension(ifms:ifme,jfms:jfme):: &
 138     lfn_0,lfn_1,lfn_2,lfn_s0,lfn_s1,lfn_s2,lfn_s3,flame_length,ros_front           ! level function stages
 139
 140 REAL, INTENT(out), dimension(ifms:ifme,jfms:jfme):: &
 141     fire_area                               ! fraction of each cell burning
 142
 143 REAL, INTENT(out), dimension(ifms:ifme,jfms:jfme):: &
 144     burnt_area_dt
 145
 146 ! output
 147 REAL, INTENT(out), dimension(ifms:ifme,jfms:jfme):: &
 148     lfn_out, &                              !
 149     grnhfx,grnqfx, &                        ! heat fluxes J/m^2/s  (cell)
 150     ros                                     ! output: rate of spread
 151
 152 ! constant arrays - set at initialization
 153 real, intent(inout), dimension(ifms:ifme, jfms:jfme)::nfuel_cat ! cell based, data, constant
 154 real,intent(inout),dimension(ifms:ifme,jfms:jfme):: fuel_time
 155 type(fire_params),intent(inout)::fp
 156
 157 !*** local
 158
 159 integer :: xifms,xifme,xjfms,xjfme  ! memory bounds for pass-through arguments to normal spread
 160 real, dimension(ifts:ifte,jfts:jfte)::fuel_frac_burnt,fuel_frac_end
 161 integer::ignited,ig,i,j,itso,iteo,jtso,jteo
 162 real::tbound,err,erri,errj,maxgrad,grad,tfa,thf,mhf,tqf,mqf,aw,mw
 163 character(len=128)::msg
 164 logical:: freeze_fire
 165 integer:: stat_lev=1
 166
 167     ! PAJ:
 168   real :: start_time_ig, end_time_ig
 169   real, parameter :: EPSILON = 0.00001
 170
 171 !*** executable
 172
 173 call check_mesh_2dim(ifts-1,ifte+1,jfts-1,jfte+1,ifms,ifme,jfms,jfme)
 174
 175 xifms=ifms  ! dimensions for the include file
 176 xifme=ifme
 177 xjfms=jfms
 178 xjfme=jfme
 179
 180
 181 ! init flags
 182 need_lfn_update=.false.
 183 ignitions_done=0
 184 freeze_fire = fire_const_time > 0. .and. time_start < fire_const_time
 185
 186 if(ifun.eq.1)then       ! do nothing, init pass 1 is outside only
 187 elseif(ifun.eq.2)then
 188         ! initialize all arrays that the model will not change later
 189
 190         ! assuming halo on zsf done
 191         ! extrapolate on 1 row of cells beyond the domain boundary
 192         ! including on the halo regions
 193
 194         call continue_at_boundary(1,1,0., & ! do x direction or y direction
 195             ifms,ifme,jfms,jfme,           &                ! memory dims
 196             ifds,ifde,jfds,jfde, &                     ! domain dims
 197             ifps,ifpe,jfps,jfpe, &            ! patch dims - winds defined up to +1
 198             ifts,ifte,jfts,jfte, &                ! tile dims
 199             itso,iteo,jtso,jteo, &              ! where set now
 200             fp%zsf)                               ! array
 201
 202 !       compute the gradients once for all
 203         err=0.
 204         maxgrad=0.
 205         do j=jfts,jfte
 206             do i=ifts,ifte
 207                 erri = fp%dzdxf(i,j) - (fp%zsf(i+1,j)-fp%zsf(i-1,j))/(2.*fdx)
 208                 errj = fp%dzdyf(i,j) - (fp%zsf(i,j+1)-fp%zsf(i,j-1))/(2.*fdy)
 209                 err=max(err,abs(erri),abs(errj))
 210                 grad=sqrt(fp%dzdxf(i,j)**2+fp%dzdyf(i,j)**2)
 211                 maxgrad=max(maxgrad,grad)
 212             enddo
 213         enddo
 214 !$OMP CRITICAL(FIRE_MODEL_CRIT)
 215         write(msg,*)'max gradient ',maxgrad,' max error against zsf',err
 216 !$OMP END CRITICAL(FIRE_MODEL_CRIT)
 217         call message(msg)
 218
 219         if(.not.restart)call set_nfuel_cat( &
 220             ifms,ifme,jfms,jfme, &
 221             ifts,ifte,jfts,jfte, &
 222             ifuelread,nfuel_cat0,&
 223             fp%zsf,nfuel_cat)            ! better not use the extrapolated zsf!!
 224
 225         ! uses nfuel_cat to set the other fuel data arrays
 226         ! needs zsf on halo width 1 to compute the terrain gradient
 227         if(.not.restart)call set_fire_params(   &
 228             ifds,ifde,jfds,jfde, &
 229             ifms,ifme,jfms,jfme, &
 230             ifts,ifte,jfts,jfte, &
 231             fdx,fdy,nfuel_cat0,  &
 232             nfuel_cat,fuel_time, &
 233             fp  &
 234 )
 235
 236         ! initialize model state to no fire
 237         if(.not.restart)then
 238             call init_no_fire  ( &
 239             ifds,ifde,jfds,jfde, &
 240             ifms,ifme,jfms,jfme, &
 241             ifts,ifte,jfts,jfte, &
 242             fdx,fdy,time_start,  &
 243             fuel_frac,fire_area,lfn,tign)
 244
 245             need_lfn_update=.true. ! because we have set lfn
 246
 247         endif
 248
 249 elseif(ifun.eq.3)then   ! ignition if so specified
 250
 251
 252 elseif (ifun.eq.4) then  ! do the timestep
 253
 254     if(run_fuel_moisture)then
 255         ! fuel moisture may have changed, reset the precomputed ros parameters
 256         ! uses nfuel_cat to set the other fuel data arrays
 257         ! needs zsf on halo width 1 to compute the terrain gradient
 258         call set_fire_params(    &       ! also on restart
 259             ifds,ifde,jfds,jfde, &
 260             ifms,ifme,jfms,jfme, &
 261             ifts,ifte,jfts,jfte, &
 262             fdx,fdy,nfuel_cat0,  &
 263             nfuel_cat,fuel_time, &
 264             fp)
 265     endif
 266
 267     if(fire_print_msg.ge.stat_lev)then
 268       aw=fun_real(RNRM_SUM,  &
 269         ifms,ifme,1,1,jfms,jfme, &                ! memory dims
 270         ifds,ifde,1,1,jfds,jfde, &                ! domain dims
 271         ifts,ifte,1,1,jfts,jfte, &                ! patch or tile dims
 272         0,0,0,       &                            ! staggering
 273         fp%vx,fp%vy)/((ifde-ifds+1)*(jfde-jfds+1))
 274       mw=fun_real(RNRM_MAX,  &
 275         ifms,ifme,1,1,jfms,jfme, &                ! memory dims
 276         ifds,ifde,1,1,jfds,jfde, &                ! domain dims
 277         ifts,ifte,1,1,jfts,jfte, &                ! patch or tile dims
 278         0,0,0,       &                            ! staggering
 279         fp%vx,fp%vy)
 280 !$OMP MASTER
 281       write(msg,91)time_start,'Average wind        ',aw,'m/s'
 282       call message(msg,stat_lev)
 283       write(msg,91)time_start,'Maximum wind        ',mw,'m/s'
 284       call message(msg,stat_lev)
 285 !$OMP END MASTER
 286     endif
 287
 288 !   compute fuel fraction at start
 289 !    call fuel_left( &
 290 !        ifms,ifme,jfms,jfme, &
 291 !        ifts,ifte,jfts,jfte, &
 292 !        ifms,ifme,jfms,jfme, &
 293 !        lfn,tign,fuel_time,time_start,fuel_frac,fire_area) ! fuel frac is shared
 294
 295     call print_2d_stats(ifts,ifte,jfts,jfte, &
 296                    ifms,ifme,jfms,jfme, &
 297                    fuel_frac,'model: fuel_frac start')
 298
 299     ! advance the model from time_start to time_start+dt
 300     ! return the fuel fraction burnt this call in each fire cell
 301     ! will call module_fr_fire_speed::normal_spread for propagation speed
 302     ! We cannot simply compute the spread rate here because that will change with the
 303     ! angle of the wind and the direction of propagation, thus it is done in subroutine
 304     ! normal_spread at each fire time step. Instead, we pass arguments that
 305     ! the speed function may use as fp.
 306
 307 !   propagate level set function in time
 308 !   set lfn_out tign
 309 !   lfn does not change, tign has no halos
 310
 311     if(.not. freeze_fire)then
 312
 313          call prop_ls_rk3(id,                           &
 314              ifds,ifde,jfds,jfde,                       &
 315              ifms,ifme,jfms,jfme,                       &
 316              ifps,ifpe,jfps,jfpe,                       &
 317              ifts,ifte,jfts,jfte,                       &
 318              time_start,dt,fdx,fdy,tbound,              &
 319              lfn,                                       &
 320              lfn_0,lfn_1,lfn_2,                         &
 321              lfn_out,tign,ros, fp,                      &
 322              grid,                                      &
 323              ids_hu,ide_hu,jds_hu,jde_hu,kds_hu,kde_hu, &
 324              ims_hu,ime_hu,jms_hu,jme_hu,kms_hu,kme_hu, &
 325              ips_hu,ipe_hu,jps_hu,jpe_hu,kps_hu,kpe_hu  &
 326          )
 327
 328        call tign_update(ifts,ifte,jfts,jfte, &
 329                         ifms,ifme,jfms,jfme, &
 330                         ifds,jfds,ifde,jfde, &
 331                         time_start,dt,       &
 332                         lfn,lfn_out,tign     &
 333                        )
 334
 335        call calc_flame_length(ifts,ifte,jfts,jfte, &
 336                               ifms,ifme,jfms,jfme, &
 337                               ros,fp%iboros,flame_length,ros_front,fire_area)
 338
 339        if (fire_lsm_reinit) then ! DME added call to reinitialize level-set function
 340
 341            call reinit_ls_rk3(id,                                    &
 342                           ifts,ifte,jfts,jfte,                       &
 343                           ifms,ifme,jfms,jfme,                       &
 344                           ifds,ifde,jfds,jfde,                       &
 345                           ifps,ifpe,jfps,jfpe,                       &
 346                           time_start,dt,fdx,fdy,                     &
 347                           lfn,                                       &
 348                           lfn_2,lfn_s0,lfn_s1,lfn_s2,lfn_s3,         &
 349                           lfn_out,tign,                              &
 350                           grid,                                      &
 351                           ids_hu,ide_hu,jds_hu,jde_hu,kds_hu,kde_hu, &
 352                           ims_hu,ime_hu,jms_hu,jme_hu,kms_hu,kme_hu, &
 353                           ips_hu,ipe_hu,jps_hu,jpe_hu,kps_hu,kpe_hu  &
 354                          )
 355        endif ! fire_lsm_reinit
 356
 357     else
 358         call message('fire_model: EXPERIMENTAL: skipping fireline propagation')
 359
 360     endif
 361
 362 elseif (ifun.eq.5) then ! copy the result of timestep back to input
 363     ! this cannot be done in the time step itself because of race condition
 364     ! some thread may still be using lfn as input in their tile halo
 365
 366     if (.not. freeze_fire) then
 367
 368        do j=jfts,jfte
 369            do i=ifts,ifte
 370                lfn(i,j)=lfn_out(i,j)
 371                ! if want to try timestep again treat tign the same way here
 372                ! even if tign does not need a halo
 373            enddo
 374        enddo
 375
 376     endif
 377
 378     ! check for ignitions
 379 !paj
 380   ig = 1
 381   start_time_ig = ignition_line(ig)%start_time
 382   end_time_ig  = ignition_line(ig)%end_time
 383
 384   if ( fire_is_real_perim .and. time_start >= start_time_ig .and. time_start < start_time_ig + dt) then
 385     ignited = 0
 386     do j = jfts, jfte
 387       do i = ifts, ifte
 388         lfn(i, j) = lfn_hist(i, j)
 389         if (abs(lfn(i, j)) < EPSILON) then
 390           tign(i, j)= time_start
 391           ignited = ignited + 1
 392         end if
 393       enddo
 394     enddo
 395   elseif (.not. fire_is_real_perim) then
 396     do ig = 1,num_ignitions
 397
 398 !  for now, check for ignition every time step...
 399 !        if(ignition_line(ig)%end_time>=time_start.and.ignition_line(ig)%start_time<time_start+dt)then
 400             call ignite_fire(                             &
 401                 ifds,ifde,jfds,jfde,                      & ! fire domain dims - the whole domain
 402                 ifms,ifme,jfms,jfme,                      &
 403                 ifts,ifte,jfts,jfte,                      &
 404                 ignition_line(ig),                        &
 405                 time_start,time_start+dt,                 &
 406                 coord_xf,coord_yf,unit_xf,unit_yf,        &
 407                 lfn,tign,ignited)
 408
 409             ignitions_done=ignitions_done+1
 410             ignited_tile(ignitions_done)=ignited
 411
 412 !            need_lfn_update=.true. ! if ignition, lfn changed
 413 #ifdef DEBUG_OUT
 414             call write_array_m(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,lfn,'lfn_ig',id)
 415             call write_array_m(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,coord_xf,'coord_xf_ig',id)
 416             call write_array_m(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,coord_yf,'coord_yf_ig',id)
 417 #endif
 418 !        endif
 419
 420     enddo
 421  end if
 422
 423
 424     call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme, &
 425                    lfn,'fire_model: lfn out')
 426
 427
 428     need_lfn_update=.true. ! duh
 429
 430 elseif (ifun.eq.6) then ! timestep postprocessing
 431
 432   if(.not. freeze_fire)then
 433
 434     ! compute the heat fluxes from the fuel burned
 435     ! needs lfn and tign from neighbors so halo must be updated before
 436     call fuel_left(&
 437         ifms,ifme,jfms,jfme, &
 438         ifts,ifte,jfts,jfte, &
 439         ifts,ifte,jfts,jfte, &
 440         lfn,tign,fuel_time,time_start+dt,fuel_frac_end,fire_area) !fuel_frac_end is private and tile based
 441
 442     call print_2d_stats(ifts,ifte,jfts,jfte, &
 443                    ifts,ifte,jfts,jfte, &
 444                    fuel_frac_end,'model: fuel_frac end')
 445
 446     do j=jfts,jfte
 447         do i=ifts,ifte
 448             fuel_frac_burnt(i,j)=fuel_frac(i,j)-fuel_frac_end(i,j) ! fuel lost this timestep
 449             burnt_area_dt(i,j)=fuel_frac_burnt(i,j)
 450             fuel_frac(i,j)=fuel_frac_end(i,j) ! copy new value to state array
 451         enddo
 452     enddo
 453
 454     call print_2d_stats(ifts,ifte,jfts,jfte, &
 455                    ifts,ifte,jfts,jfte, &
 456                    fuel_frac_burnt,'model: fuel_frac burned')
 457
 458     call heat_fluxes(dt,fp,                       &
 459         ifms,ifme,jfms,jfme,                      &
 460         ifts,ifte,jfts,jfte,                      &
 461         ifts,ifte,jfts,jfte,                      &  ! fuel_frac_burned is tile dimensioned
 462         fp%fgip,                                     &
 463         fuel_frac_burnt,                          & !
 464         grnhfx,grnqfx)                              !out
 465
 466     if(fire_print_msg.ge.stat_lev)then
 467       tfa=fun_real(REAL_SUM,  &
 468         ifms,ifme,1,1,jfms,jfme, &                ! memory dims
 469         ifds,ifde,1,1,jfds,jfde, &                ! domain dims
 470         ifts,ifte,1,1,jfts,jfte, &                ! patch or tile dims
 471         0,0,0,       &                            ! staggering
 472         fire_area,fire_area) * fdx * fdy
 473       thf=fun_real(REAL_SUM,  &
 474         ifms,ifme,1,1,jfms,jfme, &                ! memory dims
 475         ifds,ifde,1,1,jfds,jfde, &                ! domain dims
 476         ifts,ifte,1,1,jfts,jfte, &                ! patch or tile dims
 477         0,0,0,       &                            ! staggering
 478         grnhfx,grnhfx) * fdx * fdy
 479       mhf=fun_real(REAL_MAX,  &
 480         ifms,ifme,1,1,jfms,jfme, &                ! memory dims
 481         ifds,ifde,1,1,jfds,jfde, &                ! domain dims
 482         ifts,ifte,1,1,jfts,jfte, &                ! patch or tile dims
 483         0,0,0,       &                            ! staggering
 484         grnhfx,grnhfx)
 485       tqf=fun_real(REAL_SUM,  &
 486         ifms,ifme,1,1,jfms,jfme, &                ! memory dims
 487         ifds,ifde,1,1,jfds,jfde, &                ! domain dims
 488         ifts,ifte,1,1,jfts,jfte, &                ! patch or tile dims
 489         0,0,0,       &                            ! staggering
 490         grnqfx,grnqfx) * fdx * fdy
 491       mqf=fun_real(REAL_MAX,  &
 492         ifms,ifme,1,1,jfms,jfme, &                ! memory dims
 493         ifds,ifde,1,1,jfds,jfde, &                ! domain dims
 494         ifts,ifte,1,1,jfts,jfte, &                ! patch or tile dims
 495         0,0,0,       &                            ! staggering
 496         grnqfx,grnqfx)
 497 !$OMP MASTER
 498       write(msg,91)time_start,'Fire area           ',tfa,'m^2'
 499       call message(msg,stat_lev)
 500       write(msg,91)time_start,'Heat output         ',thf,'W'
 501       call message(msg,stat_lev)
 502       write(msg,91)time_start,'Max heat flux       ',mhf,'W/m^2'
 503       call message(msg,stat_lev)
 504       write(msg,91)time_start,'Latent heat output  ',tqf,'W'
 505       call message(msg,stat_lev)
 506       write(msg,91)time_start,'Max latent heat flux',mqf,'W/m^2'
 507       call message(msg,stat_lev)
 508 !$OMP END MASTER
 509 91  format('Time ',f11.3,' s ',a,e12.3,1x,a)
 510     endif
 511
 512
 513   else
 514      call message('fire_model: EXPERIMENTAL: skipping fuel burnt computation')
 515
 516      if (fire_const_grnhfx >= 0. .and. fire_const_grnqfx >= 0.) then
 517
 518 !$OMP CRITICAL(FIRE_MODEL_CRIT)
 519         write(msg,'(a,2e12.3,a)')'fire_model: EXPERIMENTAL output constant heat flux', &
 520            fire_const_grnhfx, fire_const_grnqfx, ' W/s'
 521 !$OMP END CRITICAL(FIRE_MODEL_CRIT)
 522         call message(msg)
 523
 524         do j=jfts,jfte
 525             do i=ifts,ifte
 526                 grnhfx(i,j)=fire_const_grnhfx
 527                 grnqfx(i,j)=fire_const_grnqfx
 528             enddo
 529         enddo
 530
 531       endif
 532
 533    endif
 534
 535     call print_2d_stats(ifts,ifte,jfts,jfte, &
 536                    ifms,ifme,jfms,jfme, &
 537                    grnhfx,'model: heat flux(J/m^2/s)')
 538
 539 else
 540 !$OMP CRITICAL(FIRE_MODEL_CRIT)
 541     write(msg,*)'fire_model: bad ifun=',ifun
 542 !$OMP END CRITICAL(FIRE_MODEL_CRIT)
 543     call crash(msg)
 544 endif
 545
 546 end subroutine fire_model
 547
 548 !
 549 !*****************
 550 !
 551
 552 end module module_fr_fire_model