phys/module_fr_sfire_driver.F

   1 ! SFIRE - Spread fire model in WRF-Fire
   2 !
   3 !*** Jan Mandel August 2007 - March 2011
   4 !*** email: Jan.Mandel@gmail.com
   5
   6 ! For support please subscribe to the wrf-fire mailing list at NCAR at
   7 ! http://mailman.ucar.edu/mailman/listinfo/wrf-fire
   8 ! or go to http://www.openwfm.org/wiki/WRF-Fire_user_support
   9
  10 ! Current drafts of the technical documentation and
  11 ! user's guide can be found at
  12
  13 ! http://www.openwfm.org/wiki/WRF-Fire_documentation
  14 ! http://www.openwfm.org/wiki/WRF-Fire_publications
  15
  16 ! This module is the only entry point from WRF-ARW to the wildland
  17 ! fire model. The call to sfire_driver advances the fire model by
  18 ! one timestep. The fire model inputs the wind and outputs
  19 ! temperature and humidity tendencies. The fire model also inputs a
  20 ! number of constant arrays (fuel data, topography). Additional
  21 ! arguments are model state (for data assimilation) and constant arrays
  22 ! the model gives to WRF for safekeeping because it is not allowed
  23 ! to save anything.
  24
  25 ! This code as of mid-2011 is described in [1]. If you use this code,
  26 ! please acknowledge our work by citing [1].
  27 ! Thank you.
  28
  29 ! Acknowledgements
  30 !
  31 ! The fire physics code is adapted from the CAWFE code [2].
  32 ! The coupling with WRF is adapted from a code by Ned Patton,
  33 ! coupling a Fortran 90 port of the CAWFE fire module to WRF [3].
  34 ! Support of refined fire grids in WRF was provided by John Michalakes.
  35 ! Jonathan D. Beezley has set up and maintained the WRF build and
  36 ! execution environment, provided software engineering infrastructure
  37 ! including synchronization with the WRF repository, and was responsibe
  38 ! for all aspects of WRF modification. UCD students Minjeong Kim and
  39 ! Volodymyr Kondratenko have contributed to the implementation of the
  40 ! fire propagation by the level set method.
  41
  42 ! Refefences
  43 !
  44 ! [1] Jan Mandel, Jonathan D. Beezley, and Adam K. Kochanski, "Coupled
  45 ! atmosphere-wildland fire modeling with WRF 3.3 and SFIRE 2011,
  46 ! Geoscientific Model Development (GMD) 4, 591-610, 2011.
  47 ! doi:10.5194/gmd-4-591-2011
  48 !
  49 ! [2] T. L. Clark, J. Coen, and D. Latham, Description of a coupled
  50 ! atmosphere-fire model, Intl. J. Wildland Fire, vol. 13, pp. 49-64,
  51 ! 2004
  52 !
  53 ! [3] Edward G. Patton and Janice L. Coen, WRF-Fire: A Coupled
  54 ! Atmosphere-Fire Module for WRF, Preprints of Joint MM5/Weather
  55 ! Research and Forecasting Model Users' Workshop, Boulder, CO,
  56 ! June 22-25, 2004, pp. 221-223, NCAR
  57 !
  58 ! ---------------------------------------------
  59 !
  60 ! CURRENT ACTIVITY
  61 !
  62 ! For current activity and development trends please check out
  63 ! http://ccm.ucdenver.edu/wiki/User:Jmandel/blog
  64 ! http://www.openwfm.org/wiki/WRF-Fire_development_notes
  65 !
  66
  67 module module_fr_sfire_driver
  68 ! use this module for standalone call, you only need to provide some mock-up wrf modules
  69
  70 use module_fr_sfire_model, only: sfire_model
  71 use module_fr_sfire_phys, only: fire_params, init_fuel_cats, fuel_moisture, &
  72    advance_moisture, moisture_classes, &
  73    fire_rate_of_spread, ibeh, balbi_msglevel
  74 use module_fr_sfire_atm, only: apply_windrf,interpolate_wind2fire_height,interpolate_atm2fire, &
  75    interpolate_z2fire,setup_wind_log_interpolation,find_trees_fmesh,massman_fwh, rho2fire
  76 use module_fr_sfire_util
  77 !use module_fr_sfire_util, only: lines_type,fire_max_lines
  78 ! Driver layer modules
  79 #ifdef DM_PARALLEL
  80     USE module_dm        , ONLY : ntasks_x,ntasks_y,local_communicator,mytask,ntasks,wrf_dm_sum_reals
  81     USE module_comm_dm , ONLY : halo_fire_fuel1_sub, halo_fire_tign_sub, halo_fire_wind_f_sub, &
  82 HALO_SFIRE_WIND_A_sub, halo_fire_ph_sub, halo_fire_zsf_sub, halo_fire_longlat_sub, &
  83 halo_fire_phb_sub, halo_fire_z0_sub, halo_fire_lfn1_sub, HALO_FIRE_LFN_OUT_sub, &
  84 HALO_FIRE_MAG_sub, HALO_FIRE_MFG_sub, halo_fire_ndwi_sub, &
  85 HALO_SFIRE_RHO2FIRE_ATM_sub, HALO_SFIRE_RHO2FIRE_SUB_sub
  86 #endif
  87 use module_fr_sfire_atm, only: read_emissions_table, add_fire_emissions
  88
  89
  90 ! WRF dependencies
  91 USE module_domain, only: domain
  92 USE module_configure, only: grid_config_rec_type
  93 use module_model_constants, only: reradius,    & ! 1/earth radiusw
  94                                   pi2            ! 2*pi
  95
  96 implicit none
  97
  98
  99 private
 100 public sfire_driver_em,use_atm_vars,set_flags, &
 101        set_fp_from_grid, fire_ignition_convert
 102 public ifun_beg, ifun_step, ifun_end
 103
 104 logical:: use_atm_vars=.true.   !  interpolate wind from atm mesh and average output fluxes back
 105 logical:: interpolate_long_lat=.true. ! get fxlong fxlat by interpolation
 106
 107 logical:: fmoist_run, fmoist_interp, fire_run  ! which kind of model to run overall
 108
 109 integer, parameter:: ifun_beg=1, ifun_step=3, ifun_end=6
 110
 111 TYPE(lines_type):: ignition, hfx
 112
 113 contains
 114
 115 ! to write debugging information
 116 #define DEBUG_OUT
 117
 118 subroutine sfire_driver_em ( grid , config_flags                    &
 119             ,time_step_start,dt                                     &
 120             ,fire_ifun_start,fire_ifun_end,tsteps                   &
 121             ,ids,ide, kds,kde, jds,jde                              &
 122             ,ims,ime, kms,kme, jms,jme                              &
 123             ,ips,ipe, kps,kpe, jps,jpe                              &
 124             ,ifds,ifde, jfds,jfde                                   &
 125             ,ifms,ifme, jfms,jfme                                   &
 126             ,ifps,ifpe, jfps,jfpe                                   &
 127             ,rho,z_at_w,dz8w                                        &
 128 )
 129
 130 !*** purpose: driver from grid structure
 131
 132
 133
 134     implicit none
 135 !*** arguments
 136     TYPE(domain) , TARGET :: grid                             ! state
 137     TYPE (grid_config_rec_type) , INTENT(IN)  :: config_flags ! namelist
 138     real, intent(in):: time_step_start, dt
 139     integer, intent(in)::     fire_ifun_start,fire_ifun_end,tsteps ! driver cycle controls
 140     integer, intent(in):: &
 141              ids,ide, kds,kde, jds,jde,                              &
 142              ims,ime, kms,kme, jms,jme,                              &
 143              ips,ipe, kps,kpe, jps,jpe,                              &
 144              ifds,ifde, jfds,jfde,                                   &
 145              ifms,ifme, jfms,jfme,                                   &
 146              ifps,ifpe, jfps,jfpe
 147     real,dimension(ims:ime, kms:kme, jms:jme),intent(in), optional::rho,z_at_w,dz8w
 148
 149 !*** local
 150     integer::fire_ifun,ir,jr,istep,itimestep,i,ipe1,kpe1,jpe1,j
 151     logical::restart,replay
 152     real:: corner_ll,corner_ul,corner_ur,corner_lr, max_u, max_v, max_w, max_rho, min_rho
 153     character(len=128)msg, msg2
 154     type(fire_params)::fp
 155     real:: moisture_time
 156
 157     logical:: run_advance_moisture,run_fuel_moisture, moisture_initializing
 158     real::    dt_moisture
 159
 160
 161
 162 !*** executable
 163
 164     call sfire_debug_hook(config_flags%fire_debug_hook_sec)
 165     call time_start
 166
 167 !   if(fire_ifun_start.le.1)call print_id  ! print id only once, during initialization
 168
 169
 170 ! **** THE FOLLOWING REALLY SHOULD BE DONE ONCE NOT EVERY TIMESTEP
 171
 172 ! populate our structures from wrf
 173
 174     call set_fp_from_grid(grid,fp)
 175
 176 ! copy configuration flags to sfire internal structures
 177     call set_flags(config_flags)
 178
 179     write(msg,'(a,i3)')' domain', grid%grid_id
 180
 181     call message(msg,level=1)
 182
 183 ! prevent memory errors because of variables missing in registry packages
 184     call check_grid_alloc(grid,config_flags)
 185
 186 #ifndef SFIRE_STANDALONE
 187
 188 !   see what we got from wrf
 189 !! need to replace ipe by min(ide-1,ipe) and similarly jpe
 190
 191     if(fire_print_msg.ge.2 .and. fire_ifun_start .gt. 1)then
 192
 193       ipe1=min(ide-1,ipe)
 194       jpe1=min(jde-1,jpe)
 195       kpe1=kpe-1
 196
 197       max_u=fun_real(REAL_AMAX,  &
 198         ims,ime,kms,kme,jms,jme, &                ! memory dims
 199         ids,ide,kds,kde,jds,jde, &                ! domain dims
 200         ips,ipe1,kps,kpe1,jps,jpe1, &                ! patch or tile dims
 201         1,0,0,       &                            ! staggering
 202         grid%u_2,grid%u_2)
 203
 204       max_v=fun_real(REAL_AMAX,  &
 205         ims,ime,kms,kme,jms,jme, &                ! memory dims
 206         ids,ide,kds,kde,jds,jde, &                ! domain dims
 207         ips,ipe1,kps,kpe1,jps,jpe1, &                ! patch or tile dims
 208         0,0,1,       &                            ! staggering
 209         grid%v_2,grid%v_2)
 210
 211       max_w=fun_real(REAL_AMAX,  &
 212         ims,ime,kms,kme,jms,jme, &                ! memory dims
 213         ids,ide,kds,kde,jds,jde, &                ! domain dims
 214         ips,ipe1,kps,kpe1,jps,jpe1, &                ! patch or tile dims
 215         0,1,0,       &                            ! staggering
 216         grid%w_2,grid%w_2)
 217
 218       !write(msg,93)time_step_start,'Maximal u v w  wind',max_u,max_v,max_w,'m/s'
 219       !call message(msg,0)
 220       !write(msg,92)time_step_start,'Min and max rho    ',min_rho,max_rho,'kg/m^3'
 221       !call message(msg,0)
 222
 223       write(msg,91)time_step_start,'Maximal u wind      ',max_u,'m/s'
 224       call message(msg,0)
 225       write(msg,91)time_step_start,'Maximal v wind      ',max_v,'m/s'
 226       call message(msg,0)
 227       write(msg,91)time_step_start,'Maximal w wind      ',max_w,'m/s'
 228       call message(msg,0)
 229
 230       if (present(rho)) then
 231
 232         max_rho=fun_real(REAL_MAX,  &
 233           ims,ime,kms,kme,jms,jme, &                ! memory dims
 234           ids,ide,kds,kde,jds,jde, &                ! domain dims
 235           ips,ipe1,kps,kpe1,jps,jpe1, &                ! patch or tile dims
 236           0,0,0,       &                            ! staggering
 237           rho,rho)
 238
 239         min_rho=fun_real(REAL_MIN,  &
 240           ims,ime,kms,kme,jms,jme, &                ! memory dims
 241           ids,ide,kds,kde,jds,jde, &                ! domain dims
 242           ips,ipe1,kps,kpe1,jps,jpe1, &                ! patch or tile dims
 243           0,0,0,       &                            ! staggering
 244           rho,rho)
 245
 246
 247         write(msg,91)time_step_start,'Minimal rho         ',min_rho,'kg/m^3'
 248         call message(msg,0)
 249         write(msg,91)time_step_start,'Maximal rho         ',max_rho,'kg/m^3'
 250         call message(msg,0)
 251
 252  !       write(msg,*)'max diff rho and grid%rho ', &
 253  !           maxval(abs(rho(ips:ipe,kps:kpe1,jps:jpe)-grid%rho(ips:ipe,kps:kpe1,jps:jpe)))
 254  !       call message(msg,0)
 255
 256       endif
 257
 258
 259 93    format('Time ',f11.3,' s ',a,3e12.3,1x,a)
 260 92    format('Time ',f11.3,' s ',a,2e12.3,1x,a)
 261 91    format('Time ',f11.3,' s ',a,e12.3,1x,a)
 262
 263
 264     endif
 265 #endif
 266
 267
 268     ! refinement r
 269     ir=grid%sr_x ! refinement ratio
 270     jr=grid%sr_y
 271     write(msg,'(a,2i4)')'fire mesh refinement ratios', ir,jr
 272     call message(msg)
 273     if(ir.le.0.or.jr.le.0)then
 274         call crash('fire mesh refinement ratio must be positive')
 275     endif
 276
 277 call print_2d_stats(ifps,min(ifpe,ifde-ir),jfps,min(jfpe,jfde-jr),ifms,ifme,jfms,jfme, grid%tign_g,'sfire_driver_em: grid%tign_g')
 278 call print_2d_stats(ifps,min(ifpe,ifde-ir),jfps,min(jfpe,jfde-jr),ifms,ifme,jfms,jfme, grid%nfuel_cat,'sfire_driver_em: grid%nfuel_cat')
 279
 280     itimestep=grid%itimestep
 281     restart=config_flags%restart .or. config_flags%cycling .or. config_flags%fire_restart ! skip state initialization
 282     replay= time_step_start+dt .le.  config_flags%fire_perimeter_time
 283 94  format('Time step',i11,' from',f11.3,' to',f11.3,' perimeter_time',f11.3,' setting replay ',l1)
 284     write(msg,94)itimestep,time_step_start,time_step_start+dt,config_flags%fire_perimeter_time,replay
 285     call message(msg)
 286 95  format('namelist.input restart ',l1,' cycling ',l1,' fire_restart ',l1,' setting restart ',l1)
 287     write(msg,95)config_flags%restart,config_flags%cycling,config_flags%fire_restart,restart
 288     call message(msg)
 289
 290
 291     ! **** moisture model
 292
 293     ! decide what to run - moisture, interpolation, or fire model itself
 294     fmoist_run    = config_flags%fmoist_run
 295     fmoist_interp = config_flags%fmoist_interp
 296     if(fire_fmc_read.ne.0.and.fmoist_run)call crash('fmoist_run=T requires fire_fmc_read=0')
 297     fire_run = .not. config_flags%fmoist_only
 298
 299     !decide what to run
 300     moisture_time = time_step_start
 301     run_advance_moisture = .false. ! default
 302     run_fuel_moisture = .false. ! default
 303     moisture_initializing = fire_ifun_start < 3
 304
 305
 306
 307     if(fmoist_run)then
 308         if(moisture_initializing)then
 309             if(fire_ifun_end>2)call crash('initialization must be run separately')
 310             grid%fmoist_lasttime=moisture_time ! initialize the last time the model has run to start of run
 311             grid%fmoist_nexttime=moisture_time
 312             call message('moisture initialization')
 313             run_advance_moisture = .true.
 314         else ! regular timestep
 315             if(config_flags%fmoist_freq > 0)then  ! regular timestep. go by multiples?
 316                 if(mod(grid%itimestep,config_flags%fmoist_freq) .eq. 0)then
 317                     write(msg,'(a,i10,a,i10)')'moisture model runs because timestep ',grid%itimestep,' is a multiple of ',config_flags%fmoist_freq
 318                     call message(msg)
 319                     run_advance_moisture = .true.
 320                 endif
 321             else
 322                 if(.not. moisture_time  < grid%fmoist_nexttime) then ! no, by time interval
 323                     write(msg,'(a,f12.2,a)')'moisture model runs because time ',grid%fmoist_nexttime,'s has arrived'
 324                     call message(msg)
 325                     run_advance_moisture = .true.
 326                 endif
 327             endif
 328             if(run_advance_moisture)then ! decide on timing
 329                 dt_moisture  = moisture_time - grid%fmoist_lasttime  ! Time since moisture model run the last time. Should be long.
 330                 grid%fmoist_lasttime = moisture_time
 331                 if(config_flags%fmoist_freq > 0)then
 332                     write(msg,'(a,f12.2,a,i10,a)')'moisture time step is ',dt_moisture,'s running every ',config_flags%fmoist_freq,' steps'
 333                     call message(msg)
 334                 else
 335                     grid%fmoist_nexttime = moisture_time + config_flags%fmoist_dt
 336                     write(msg,'(a,f12.2,a,f12.2,a)')'moisture time step is ',dt_moisture,'s next run at ',grid%fmoist_nexttime,'s'
 337                     call message(msg)
 338                 endif
 339                 if(fmoist_interp)then
 340                     call message('moisture interpolation to fuels will run because moisture model does')
 341                     run_fuel_moisture=.true.
 342                 endif
 343             endif
 344         endif
 345     elseif(itimestep.eq.1.and.fmoist_interp)then
 346             call message('initializing, moisture interpolation to fuels will run from input data')
 347             run_fuel_moisture=.true.
 348     endif
 349
 350 !$OMP CRITICAL(SFIRE_DRIVER_CRIT)
 351     write(msg,'(a,i1,a,i1,a,l1,a,l1)') &
 352        'sfire_driver_em: ifun from ',fire_ifun_start,' to ',fire_ifun_end, &
 353            ' restart=',restart,' replay=',replay
 354 !$OMP END CRITICAL(SFIRE_DRIVER_CRIT)
 355     call message(msg)
 356
 357     do istep=0,tsteps ! istep >0 is for testing only, exit after the first call
 358       itimestep = grid%itimestep + istep ! in the first call, do fire_test_steps steps of the fire model
 359
 360       do fire_ifun=fire_ifun_start,fire_ifun_end
 361
 362         write(msg,'(a,i1,a)')'*** stage ',fire_ifun,' ***'
 363         call message(msg)
 364
 365         ! 1 = moisture_initialize run pass 1: interpolate height to zsf=terrain
 366         ! 2 = initialize run pass 2: set fuel data, terrain gradient
 367         ! 3 = initialize timestep: interpolate winds, check for ignition, time step on moisture model
 368         ! 4 = do one timestep
 369         ! 5 = copy timestep output to input
 370         ! 6 = compute output fluxes
 371
 372 #ifdef DM_PARALLEL
 373
 374        if(fire_run)then
 375
 376         if(fire_ifun.eq.1)then
 377
 378 !       halo exchange on topography
 379 #include "HALO_FIRE_LONGLAT.inc"
 380 !!            if(fire_topo_from_atm.eq.1)then
 381 !!#include "HALO_FIRE_HT.inc"
 382 !!            endif
 383 ! base geopotential and roughness
 384 #include "HALO_FIRE_PHB.inc"
 385 #include "HALO_FIRE_Z0.inc"
 386         if(kfmc_ndwi > 0 .and. fndwi_from_ndwi .eq.1)then
 387 #include "HALO_FIRE_NDWI.inc"
 388         endif
 389
 390         elseif(fire_ifun.eq.2)then
 391 !           halo exchange on zsf width 2
 392 #include "HALO_FIRE_ZSF.inc"
 393
 394             if(replay)then
 395                 call message('replay, halo exchange on lfn and tign')
 396 #include "HALO_FIRE_LFN1.inc"
 397 #include "HALO_FIRE_TIGN.inc"
 398             endif
 399
 400         elseif(fire_ifun.eq.3)then
 401 !           halo exchange on atm winds and geopotential for interpolation
 402 #include "HALO_SFIRE_WIND_A.inc"
 403 #include "HALO_FIRE_PH.inc"
 404             if(ibeh.eq.2)then
 405 #include "HALO_SFIRE_RHO2FIRE_ATM.inc"
 406             endif
 407         elseif(fire_ifun.eq.4)then
 408 !           halo exchange on fire winds width 2 for a 2-step RK method
 409 #include "HALO_FIRE_WIND_F.inc"
 410             if(ibeh.eq.2)then
 411 #include "HALO_SFIRE_RHO2FIRE_SUB.inc"
 412             endif
 413             if(run_fuel_moisture)then
 414             ! have interpolated to the fire grid
 415 #include "HALO_FIRE_MFG.inc"
 416             endif
 417
 418         elseif(fire_ifun.eq.5)then
 419 #include "HALO_FIRE_LFN_OUT.inc"
 420
 421         elseif(fire_ifun.eq.6)then
 422 !           computing fuel_left needs ignition time from neighbors
 423             call message('halo exchange on lfn width 2 and tign')
 424 #include "HALO_FIRE_TIGN.inc"
 425 #include "HALO_FIRE_LFN1.inc"
 426
 427         endif ! fire_ifun
 428        endif ! fire_run
 429 #endif
 430
 431         if(fire_ifun.eq.2)then
 432             write(msg,'(a,i4,a,i4)') 'chem_opt=', config_flags%chem_opt, &
 433                ' tracer_opt=',config_flags%tracer_opt
 434             call message(msg)
 435             if(config_flags%chem_opt>0 .or. config_flags%tracer_opt > 0)then
 436                 ! need reading fuel categories first
 437                 call read_emissions_table(config_flags%chem_opt,config_flags%tracer_opt)
 438             endif ! config_flags
 439         endif  ! fire_ifun
 440
 441         ! print *,'dt: ',dt,grid%dt,' diff ', dt-grid%dt
 442         ! need domain by 1 smaller, in last row.col winds are not set properly
 443         call sfire_driver_phys ( &
 444             fire_ifun,                  &
 445             config_flags,                                           &
 446             ids,ide-1, kds,kde, jds,jde-1,                          &
 447             ims,ime, kms,kme, jms,jme,                          &
 448             ips,min(ipe,ide-1), kps,kpe, jps,min(jpe,jde-1),                          &
 449             ifds,ifde-ir, jfds,jfde-jr,                    &
 450             ifms,ifme, jfms,jfme,                    &
 451             ifps,min(ifpe,ifde-ir), jfps,min(jfpe,jfde-jr),      &
 452             ir,jr,                                      & ! atm/fire grid ratio
 453             grid%num_tiles,                             & ! atm grid tiling
 454             grid%i_start,min(grid%i_end,ide-1),                    &
 455             grid%j_start,min(grid%j_end,jde-1),                    &
 456             itimestep,restart,replay,config_flags%fire_fuel_read,config_flags%fire_fuel_cat, &  ! in scalars
 457             time_step_start,dt,grid%dx,grid%dy,                    &
 458             grid%u_frame,grid%v_frame,                  &
 459             config_flags%fire_ext_grnd,config_flags%fire_ext_crwn,config_flags%fire_crwn_hgt, &
 460             ignition,hfx,              &  ! lines
 461             grid%u_2,grid%v_2,           &          ! atm arrays in
 462             grid%ph_2,grid%phb,               & ! geopotential
 463             grid%z0,                        & ! roughness height
 464             grid%ht,                        &                         ! terrain height
 465             grid%xlong,grid%xlat,                         & ! coordinates of atm grid centers, for ignition location
 466             grid%tign_in,                                 &
 467             grid%lfn,grid%tign_g,grid%fuel_frac,          & ! state arrays, fire grid
 468             grid%fire_area,                               & ! redundant, for display, fire grid
 469             grid%fuel_frac_burnt,                         &
 470             grid%lfn_out,                                 & ! work - one timestep
 471             grid%avg_fuel_frac,                           & ! out redundant arrays, atm grid
 472             grid%grnhfx,grid%grnqfx,grid%canhfx,grid%canqfx, & ! out redundant arrays, atm grid
 473             grid%uah,grid%vah,                            &
 474             grid%fgrnhfx,grid%fgrnqfx,grid%fcanhfx,grid%fcanqfx, & ! out redundant arrays, atm grid
 475             grid%ros,grid%flineint,grid%flineint2,         & ! diagnostic variables
 476             grid%f_ros0,grid%f_rosx,grid%f_rosy,grid%f_ros,& ! fire risk spread
 477             grid%f_int,grid%f_lineint,grid%f_lineint2,     & ! fire risk intensities
 478             grid%f_ros11,grid%f_ros12,grid%f_ros13,grid%f_ros21,  & ! fire spread in nodal directions
 479             grid%f_ros23,grid%f_ros31,grid%f_ros32,grid%f_ros33,  & ! fire spread in nodal directions
 480             grid%fxlong,grid%fxlat,                           &
 481             grid%fire_hfx,                                & !
 482             grid%nfuel_cat,                               & ! input, or internal for safekeeping
 483             grid%fuel_time,                      &
 484             grid%wz0,                            &          !wrf roughness length on fire mesh
 485             grid%fz0, grid%fwh,grid%can_top,grid%cuf,grid%cvf, &
 486             grid%tempf, grid%rhof,grid%ffwidth,             & ! variables for Balbi ROS model
 487             fp,                                    & ! structure with pointers passed to spread rate calculation
 488             config_flags%nfmc,         & ! moisture model variables start
 489             run_advance_moisture,run_fuel_moisture,dt_moisture,     &    ! moisture model control
 490             config_flags%fmep_decay_tlag,                               & ! moisture extended model assim. diffs decay time lag
 491             grid%rainc, grid%rainnc,                       & ! accumulated rain from different sources
 492             grid%t2, grid%q2, grid%psfc,               & ! temperature (K), vapor contents (kg/kg), pressure (Pa) at the surface
 493             grid%rain_old,                   & ! previous value of accumulated rain
 494             grid%t2_old, grid%q2_old, grid%psfc_old,   & ! previous values of the atmospheric state at surface
 495             grid%rh_fire,                    & ! relative humidity, diagnostics
 496             grid%fmc_gc,                      & ! fuel moisture fields updated, by class, assumed set to something reasonable
 497             grid%fmep,                      & ! fuel moisture extended model parameters
 498             grid%fmc_equi,                      & ! fuel moisture fields updated, by class, equilibrium diagnostic
 499             grid%fmc_lag,                      & ! fuel moisture fields updated, by class, tendency diagnostic
 500             fp%fmc_g, &                             !  write-only alias. need to exit before using fp again
 501             fp%fmc_gc_f, &                             !  write-only alias. need to exit before using fp again
 502             grid%ndwi, &
 503             grid%fndwi)
 504
 505
 506 #ifdef DM_PARALLEL
 507           if(fire_run)then
 508               if(fire_ifun.eq.2)then
 509 !                 halo exchange on all fuel data width 2
 510 #include "HALO_FIRE_FUEL1.inc"
 511 !                 fire state was initialized
 512                   call message('halo exchange on lfn width 2')
 513 #include "HALO_FIRE_LFN1.inc"
 514               endif
 515               if(run_fuel_moisture)then
 516                   if(fire_ifun.eq.3)then
 517                      ! prepare for interpolation to the fire grid
 518 #include "HALO_FIRE_MAG.inc"
 519                   endif
 520               endif
 521           endif
 522 #endif
 523
 524         if(fire_ifun.eq.2)then
 525 !           store mesh units computed in fire_ignition_convert in stage 2
 526             grid%unit_fxlong = ignition%unit_fxlong
 527             grid%unit_fxlat = ignition%unit_fxlat
 528         endif  ! fire_ifun
 529
 530      if(fire_ifun.eq.6)then
 531          if(config_flags%chem_opt>0 .or. config_flags%tracer_opt>0)then
 532              if(.not.(present(rho).and.present(dz8w)))then
 533                  call crash('sfire_driver_em: must have rho and dz8w to call add_fire_emissions')
 534              endif
 535              call add_fire_emissions( &
 536                  config_flags%chem_opt,config_flags%tracer_opt,dt,grid%dx,grid%dy,      &
 537                  ifms,ifme,jfms,jfme,    &
 538                  ifps,ifpe,jfps,jfpe,    &               ! use patch instead of tile
 539                  ids,ide,kds,kde,jds,jde,          &
 540                  ims,ime,kms,kme,jms,jme,          &
 541                  ips,ipe,kps,kpe,jps,jpe,          &
 542                  rho,dz8w,                         & ! from atmosphere state
 543                  grid%fgip, grid%fuel_frac_burnt, grid%nfuel_cat, & ! from fire state
 544                  grid%chem,grid%tracer)              ! update/output
 545           endif
 546      endif
 547
 548
 549       enddo
 550     enddo
 551
 552     if(tsteps>0)call crash('sfire_driver_em: test run of uncoupled fire model completed')
 553     call time_end('sfire')
 554
 555 end subroutine sfire_driver_em
 556
 557 !
 558 !*******************
 559 !
 560
 561 subroutine sfire_driver_phys (ifun,    &
 562     config_flags,                                 & ! namelist.input
 563     ids,ide, kds,kde, jds,jde,                    & ! atm grid dimensions
 564     ims,ime, kms,kme, jms,jme,                    &
 565     ips,ipe, kps,kpe, jps,jpe,                    &
 566     ifds, ifde, jfds, jfde,                       & ! fire grid dimensions
 567     ifms, ifme, jfms, jfme,                       &
 568     ifps, ifpe, jfps, jfpe,                       & ! fire patch in - will use smaller
 569     ir,jr,                                        & ! atm/fire grid ratio
 570     num_tiles,i_start,i_end,j_start,j_end,        & ! atm grid tiling
 571     itimestep,restart,replay,ifuelread,nfuel_cat0,       & ! in scalars
 572     time_step_start,dt,dx,dy,                     & ! in scalars
 573     u_frame,v_frame,                              &
 574     fire_ext_grnd,fire_ext_crwn,fire_crwn_hgt,    &
 575     ignition,hfx,                                 & ! lines
 576     u,v,                                          & ! in arrays, atm grid
 577     ph,phb,                                       &
 578     z0,zs,                                        &
 579     xlong,xlat,                                   &
 580     tign_in,                                      &
 581     lfn,tign,fuel_frac,                           & ! state arrays, fire grid
 582     fire_area,                                    & ! redundant state, fire grid
 583     fuel_frac_burnt,                              &
 584     lfn_out,                                      & ! out level set function
 585     avg_fuel_frac,                                &
 586     grnhfx,grnqfx,canhfx,canqfx,                  & ! out redundant arrays, atm grid
 587     uah,vah,                                      & ! out atm grid
 588     fgrnhfx,fgrnqfx,fcanhfx,fcanqfx,              & ! out redundant arrays, fire grid
 589     ros,flineint,flineint2,                       & ! diagnostic variables
 590     f_ros0,f_rosx,f_rosy,f_ros,                   & ! fire risk spread
 591     f_int,f_lineint,f_lineint2,                   & ! fire risk intensities
 592     f_ros11,f_ros12,f_ros13,f_ros21,  & ! fire spread in nodal directions
 593     f_ros23,f_ros31,f_ros32,f_ros33,  & ! fire spread in nodal directions
 594     fxlong,fxlat,                                 & !
 595     fire_hfx,                                     & !
 596     nfuel_cat,                                    & ! in array, data, fire grid, or constant internal
 597     fuel_time,                                    & ! save constant internal data, fire grid
 598     wz0,                                          &
 599     fz0,fwh,can_top,cuf,cvf,                      & ! input fuel map/vegetation variables on fire mesh
 600     tempf,rhof,ffwidth,                           & ! variables for Balbi ROS model
 601     fp,                                           & ! fire params
 602     nfmc,                                     & ! number of fuel moisture classes
 603     run_advance_moisture,run_fuel_moisture,dt_moisture,& ! moisture model control
 604     fmep_decay_tlag,                              & ! moist. extended model assim. diffs time lag
 605     rainc,rainnc,               & ! accumulated rain from different sources
 606     t2, q2, psfc,               & ! temperature (K), vapor contents (kg/kg), pressure (Pa) at the surface
 607     rain_old,                   & ! previous value of accumulated rain
 608     t2_old, q2_old, psfc_old,   & ! previous values of the atmospheric state at surface
 609     rh_fire,                    & ! relative humidity, diagnostics
 610     fmc_gc,                     &  ! fuel moisture fields updated, by class, assumed set to something reasonable
 611     fmep,                       &  ! fuel moisture extended model parameters
 612     fmc_equi,                   &  ! fuel moisture fields updated, by class equilibrium diagnostic
 613     fmc_lag,                    &  ! fuel moisture fields updated, by class tendency diagnostic
 614     fmc_g,                      &  ! fuel moisture, alias of fp%fmc_g
 615     fmc_gc_f,                   &  ! fuel moisture, alias of fp%fmc_gc_f
 616     ndwi,                       &
 617     fndwi)                         ! ndwi on fire grid
 618
 619
 620 implicit none
 621
 622 !*** arguments
 623
 624 TYPE (grid_config_rec_type) , INTENT(IN)  :: config_flags ! namelist, needed for ignion
 625 integer, intent(in)::ifun,                        &
 626     ids,ide, kds,kde, jds,jde,                    & ! atm domain bounds
 627     ims,ime, kms,kme, jms,jme,                    & ! atm memory bounds
 628     ips,ipe, kps,kpe, jps,jpe,                    & ! atm patch bounds
 629     ifds, ifde, jfds, jfde,                       & ! fire domain bounds
 630     ifms, ifme, jfms, jfme,                       & ! fire memory bounds
 631     ifps, ifpe, jfps, jfpe,                       & ! fire patch bounds
 632     ir,jr,                                        & ! atm/fire grid refinement ratio
 633     nfmc,                                     & ! number of fuel moisture classes
 634     itimestep,                                    & ! number of this timestep
 635     ifuelread,                                    & ! how to initialize nfuel_cat:
 636                                                        ! -1=not at all, done outside
 637                                                        ! 0=from nfuel_cat0
 638                                                        ! 1=from altitude
 639                                                        ! 2=from file
 640     nfuel_cat0,                                   & ! fuel category to initialize everything to
 641     num_tiles                                       ! number of tiles
 642
 643 logical, intent(in)::restart,replay
 644
 645 integer,dimension(num_tiles),intent(in) :: i_start,i_end,j_start,j_end  ! atm grid tiling
 646
 647 real, intent(in):: &
 648     time_step_start,                              & ! time step start
 649     dt,                                           & ! time step length
 650     dx,dy,                                        & ! atm grid step
 651     u_frame,v_frame,                              & ! velocity offset
 652     fire_crwn_hgt,                                & ! lowest height crown fire heat is released (m)
 653     fire_ext_grnd,                                & ! extinction depth of ground fire heat (m)
 654     fire_ext_crwn                                   ! and for the canopy (m)
 655
 656
 657 TYPE (lines_type), intent(inout):: ignition,hfx
 658
 659 real,intent(in),dimension(ims:ime,kms:kme,jms:jme)::u,v, & ! wind velocity (m/s) (staggered atm grid)
 660                               ph, phb                      ! geopotential (w-points atm grid)
 661 real,intent(in),dimension(ims:ime, jms:jme)::   z0, &    ! roughness height
 662                                                 zs       ! terrain height
 663 real,intent(out),dimension(ims:ime,jms:jme)::&
 664     uah,                                           & ! atm wind at fire_wind_height, diagnostics
 665     vah                                              ! atm wind at fire_wind_height, diagnostics
 666
 667 real, dimension(ims:ime, jms:jme), intent(inout)::xlong, xlat, ndwi ! inout because of extension at bdry
 668
 669 real, intent(inout), dimension(ifms:ifme,jfms:jfme):: & ! fuel data; can be also set inside (cell based, fire grid)
 670     fz0,fwh,can_top,cuf,cvf,wz0,                            &
 671     nfuel_cat,fndwi
 672
 673 real, intent(inout), dimension(ifms:ifme, jfms:jfme)::     &
 674     tign_in,                                   &
 675     lfn,tign,fuel_frac,                        &     ! state: level function, ign time, fuel left
 676     lfn_out                                    ! fire wind velocities
 677
 678 real, intent(out), dimension(ifms:ifme, jfms:jfme)::  &
 679     fire_area, &                                 ! fraction of area of each cell burning
 680     fuel_frac_burnt,                           & ! fraction of mass of fuel burnt
 681     tempf,rhof,ffwidth                           ! variables for Balbi ROS model
 682
 683 real, intent(out), dimension(ims:ime, jms:jme):: &  ! redundant arrays, for display purposes only (atm grid)
 684     avg_fuel_frac,                               &  ! average fuel fraction
 685     grnhfx,                                      &  ! heat flux from ground fire (W/m^2)
 686     grnqfx,                                      &  ! moisture flux from ground fire (W/m^2)
 687     canhfx,                                      &  ! heat flux from crown fire (W/m^2)
 688     canqfx                                         ! moisture flux from crown fire (W/m^2)
 689
 690 real, intent(out), dimension(ifms:ifme, jfms:jfme):: &  ! redundant arrays, for display only, fire grid
 691     fgrnhfx,                                      &  ! heat flux from ground fire (W/m^2)
 692     fgrnqfx,                                      &  ! moisture flux from ground fire (W/m^2)
 693     fcanhfx,                                      &  ! heat flux from crown fire (W/m^2)
 694     fcanqfx,                                      &  ! moisture flux from crown fire (W/m^2)
 695     ros,flineint,flineint2,                       & ! diagnostic variables
 696     f_ros0,f_rosx,f_rosy,f_ros,                   & ! fire risk spread
 697     f_int,f_lineint,f_lineint2,                   &  ! fire risk intensities
 698     f_ros11,f_ros12,f_ros13,f_ros21,  & ! fire spread in nodal directions
 699     f_ros23,f_ros31,f_ros32,f_ros33     ! fire spread in nodal directions
 700
 701
 702 ! moisture model arguments
 703 logical, intent(in)::run_advance_moisture,run_fuel_moisture
 704 real, intent(in)::dt_moisture
 705 real, intent(in)::fmep_decay_tlag
 706 real, intent(in), dimension(ims:ime,jms:jme):: t2, q2, psfc, rainc, rainnc
 707 real, intent(inout), dimension(ims:ime,jms:jme):: t2_old, q2_old, psfc_old, rain_old
 708 real, intent(out),dimension(ims:ime,jms:jme):: rh_fire
 709 real, intent(inout), dimension(ims:ime,nfmc,jms:jme):: fmc_gc
 710 real, intent(inout), dimension(ims:ime,2,jms:jme):: fmep
 711 real, intent(out), dimension(ims:ime,nfmc,jms:jme):: fmc_equi,fmc_lag
 712 real, intent(inout), dimension(ifms:ifme,jfms:jfme):: fmc_g
 713 real, intent(inout), dimension(ifms:ifme,nfmc,jfms:jfme):: fmc_gc_f
 714
 715
 716
 717 !  ***** data (constant in time) *****
 718
 719 real, dimension(ifms:ifme, jfms:jfme), intent(inout)::fxlong,fxlat, & ! fire mesh coordinates
 720     fire_hfx
 721 real, intent(out), dimension(ifms:ifme, jfms:jfme)::fuel_time   ! fire params arrays
 722
 723 type(fire_params),intent(inout)::fp
 724
 725 !*** local
 726 real :: dxf,dyf,time_start,latm, s
 727 integer :: its,ite,jts,jte,kts,kte, &            ! tile
 728     ij,i,j,k,id,pid,ipe1,jpe1,ite1,jte1, &
 729     ii,jj,                                 &
 730     ifts,ifte,jfts,jfte                          ! fire tile
 731 character(len=128)::msg
 732 character(len=3)::kk
 733 real, parameter::zero=0.
 734
 735 !*** executable
 736
 737 ! time - assume dt does not change
 738 ! time_start = (itimestep-1) * dt     ! timestep 1 starts at 0
 739 ! print *,'time_start: ',time_start,time_step_start,' diff ', time_start-time_step_start
 740 time_start = time_step_start ! use the time passed from wrf
 741
 742 ! fire mesh step
 743 dxf=dx/ir
 744 dyf=dy/jr
 745
 746
 747 !$OMP CRITICAL(SFIRE_DRIVER_CRIT)
 748 write(msg,'(a,i5)')'sfire_driver_phys stage ',ifun
 749 call message(msg)
 750 write(msg,'(a,2f15.6)')'atmosphere mesh step:',dx,dy
 751 call message(msg)
 752 write(msg,'(a,2f15.6)')'fire mesh step:      ',dxf,dyf
 753 call message(msg)
 754 write(msg,7001)'atm domain      ','ids',ids,ide,jds,jde
 755 call message(msg)
 756 write(msg,7001)'atm memory      ','ims',ims,ime,jms,jme
 757 call message(msg)
 758 write(msg,7001)'atm patch       ','ips',ips,ipe,jps,jpe
 759 call message(msg)
 760 write(msg,7001)'fire domain     ','ifds',ifds,ifde,jfds,jfde
 761 call message(msg)
 762 write(msg,7001)'fire memory     ','ifms',ifms,ifme,jfms,jfme
 763 call message(msg)
 764 write(msg,7001)'fire patch      ','ifps',ifps,ifpe,jfps,jfpe
 765 call message(msg)
 766 !$OMP END CRITICAL(SFIRE_DRIVER_CRIT)
 767
 768 ! check mesh dimensions
 769 call check_fmesh(ids,ide,ifds,ifde,ir,'id')           ! check if atm and fire grids line up
 770 call check_fmesh(jds,jde,jfds,jfde,jr,'jd')
 771 call check_fmesh(ips,ipe,ifps,ifpe,ir,'ip')
 772 call check_fmesh(jps,jpe,jfps,jfpe,jr,'jp')
 773 call check_mesh_2dim(ips,ipe,jps,jpe,ims,ime,jms,jme)        ! check if atm patch fits in atm array
 774 call check_mesh_2dim(ifps,ifpe,jfps,jfpe,ifms,ifme,jfms,jfme) ! check if fire patch fits in fire array
 775 call check_mesh_2dim(ips,ipe,jps,jpe,ids,ide,jds,jde)        ! check if atm patch fits in atm domain
 776 call check_mesh_2dim(ifps,ifpe,jfps,jfpe,ifds,ifde,jfds,jfde) ! check if fire patch fits in fire domain
 777
 778 pid=0
 779 if(fire_print_file.gt.0)then
 780     if(itimestep.le.fire_print_file.or.mod(itimestep,fire_print_file).eq.0)pid=itimestep ! print 1-fire_print_file then every fire_print_file-th
 781 endif
 782
 783
 784 if(ifun.eq.1)then
 785          call init_fuel_cats(fmoist_run .or. fmoist_interp) ! properties of fuel categories and moisture classes from namelist.fire
 786 endif
 787
 788 if(ifun.eq.2)then
 789  call print_chsum(itimestep,ims,ime,kms,kme,jms,jme,ids,ide,kds,kde,jds,jde,ips,ipe,kps,kpe,jps,jpe,1,0,0,u,'u')
 790  call print_chsum(itimestep,ims,ime,kms,kme,jms,jme,ids,ide,kds,kde,jds,jde,ips,ipe,kps,kpe,jps,jpe,0,0,1,v,'v')
 791  call print_chsum(itimestep,ims,ime,kms,kme,jms,jme,ids,ide,kds,kde,jds,jde,ips,ipe,kps,kpe,jps,jpe,0,1,0,ph,'ph')
 792 endif ! ifun.eq.3
 793
 794 call print_chsum(itimestep,ifms,ifme,1,1,jfms,jfme,ifds,ifde,1,1,jfds,jfde,ifps,ifpe,1,1,jfps,jfpe,0,0,0,lfn,'lfn')
 795 call print_chsum(itimestep,ifms,ifme,1,1,jfms,jfme,ifds,ifde,1,1,jfds,jfde,ifps,ifpe,1,1,jfps,jfpe,0,0,0,tign,'tign')
 796
 797
 798 ! fake atm tile bounds
 799 kts=kps
 800 kte=kpe
 801
 802 ! staggered atm patch bounds
 803 ipe1=ifval(ipe.eq.ide,ipe+1,ipe)
 804 jpe1=ifval(jpe.eq.jde,jpe+1,jpe)
 805
 806 ! set up fire tiles & interpolate to fire grid
 807 !$OMP PARALLEL DO PRIVATE(ij,its,ite,jts,jte,ite1,jte1,ifts,ifte,jfts,jfte,msg,id) &
 808 !$OMP SCHEDULE(STATIC)
 809 do ij=1,num_tiles
 810
 811     id = ifval(pid.ne.0,pid+(ij-1)*10000,0) ! for print
 812
 813     ! set up tile bounds
 814     its = i_start(ij)  ! start atmospheric tile in i
 815     ite = i_end(ij)    ! end atmospheric tile in i
 816     jts = j_start(ij)  ! start atmospheric tile in j
 817     jte = j_end(ij)    ! end atmospheric tile in j
 818     ifts= (its-ids)*ir+ifds       ! start fire tile in i
 819     ifte= (ite-ids+1)*ir+ifds-1   ! end fire tile in i
 820     jfts= (jts-jds)*jr+jfds       ! start fire tile in j
 821     jfte= (jte-jds+1)*jr+jfds-1   ! end fire tile in j
 822
 823 ! staggered atm tile bounds
 824     ite1=ifval(ite.eq.ide,ite+1,ite)
 825     jte1=ifval(jte.eq.jde,jte+1,jte)
 826
 827 !$OMP CRITICAL(SFIRE_DRIVER_CRIT)
 828     write(msg,'(a,i3,1x,a,i7,1x,a,i3)')'tile=',ij,' id=',id,' ifun=',ifun
 829     call message(msg)
 830     write(msg,7001)'atm tile   ','its',its,ite,jts,jte
 831     call message(msg)
 832     write(msg,7001)'fire tile  ','ifts',ifts,ifte,jfts,jfte
 833     call message(msg)
 834 !$OMP END CRITICAL(SFIRE_DRIVER_CRIT)
 835
 836     ! check the tiles
 837     call check_mesh_2dim(its,ite,jts,jte,ips,ipe,jps,jpe)                 ! check if atm tile fits in atm patch
 838     call check_mesh_2dim(ifts,ifte,jfts,jfte,ifps,ifpe,jfps,jfpe)         ! check if fire tile fits in fire patch
 839     call check_mesh_2dim(ifts-2,ifte+2,jfts-2,jfte+2,ifms,ifme,jfms,jfme)! check if fire node tile fits in memory
 840
 841
 842 !$OMP CRITICAL(SFIRE_DRIVER_CRIT)
 843     write(msg,'(a,i6,a,2(f15.6,a))')'time step',itimestep,' at',time_start,' duration',dt,'s'
 844     call message(msg)
 845     7001 format(a,' dimensions ',a4,':',i6,' to ',i6,' by ',i6,' to ',i6)
 846     write(msg,'(a,2i9)')'refinement ratio:',ir,jr
 847 !$OMP END CRITICAL(SFIRE_DRIVER_CRIT)
 848
 849     if(run_advance_moisture)then
 850       if(ifun.eq.3)then
 851
 852       ! one timestep of the moisture model
 853           call message('advance_moisture start')
 854           call advance_moisture(    &
 855               itimestep.eq.1,             & ! initialize?
 856               ims,ime,  jms,jme,          & ! memory dimensions
 857               its,ite,  jts,jte,          & ! tile dimensions
 858               nfmc,                       & ! number of moisture fields
 859               dt_moisture,                & ! moisture model time step
 860               fmep_decay_tlag,            & ! moisture extended model assim. diffs decay tlag
 861               rainc, rainnc,              & ! accumulated rain
 862               t2, q2, psfc,               & ! temperature (K), vapor contents (kg/kg), pressure (Pa) at the surface
 863               rain_old,                   & ! previous value of accumulated rain
 864               t2_old, q2_old, psfc_old,   & ! previous values of the atmospheric state at surface
 865               rh_fire,                    & ! relative humidity, diagnostics
 866               fmc_gc,                     & ! fuel moisture fields updated, by class, assumed set to something reasonable
 867               fmep,                       & ! fuel moisture extended model parameters
 868               fmc_equi,                   & ! fuel moisture fields updated, by class equilibrium diagnostic
 869               fmc_lag                     & ! fuel moisture fields updated, by class tendency diagnostic
 870           )
 871           call message('advance_moisture end')
 872       endif
 873     endif
 874
 875     if(fire_run)then
 876
 877      if(ifun.eq.2)then
 878
 879       if(restart)then
 880
 881           call message('restart - interpolation skipped')
 882
 883       else ! restart
 884           call message ('Interpolating static data from atmosphere to fire mesh')
 885         if(kfmc_ndwi > 0 .and. fndwi_from_ndwi .eq.1)then
 886             call print_2d_stats(ips,ipe,jps,jpe,ims,ime,jms,jme,fndwi,'driver:ndwi')
 887             call interpolate_z2fire(id,0,                 & ! for debug output, <= 0 no output, extend strip
 888                 ids,ide,  jds,jde,                    & ! atm grid dimensions
 889                 ims,ime,  jms,jme,                    &
 890                 ips,ipe,jps,jpe,                              &
 891                 its,ite,jts,jte,                              &
 892                 ifds, ifde, jfds, jfde,                       & ! fire grid dimensions
 893                 ifms, ifme, jfms, jfme,                       &
 894                 ifts,ifte,jfts,jfte,                          &
 895                 ir,jr,                                        & ! atm/fire grid ratio
 896                 ndwi,                                       & ! atm grid arrays in
 897                 fndwi)                                      ! fire grid arrays out
 898             call print_2d_stats(ifps,ifpe,jfps,jfpe,ifms,ifme,jfms,jfme,fndwi,'driver:fndwi')
 899          endif
 900 !
 901 !        call print_2d_stats(ips,ipe,jps,jpe,ims,ime,jms,jme,zs,'driver:zs')
 902 !
 903 !        ! interpolate terrain height
 904 !        if(fire_topo_from_atm.eq.1)then
 905 !            call interpolate_z2fire(id,1,                 & ! for debug output, <= 0 no output
 906 !                ids,ide,  jds,jde,                    & ! atm grid dimensions
 907 !                ims,ime,  jms,jme,                    &
 908 !                ips,ipe,jps,jpe,                              &
 909 !                its,ite,jts,jte,                              &
 910 !                ifds, ifde, jfds, jfde,                       & ! fire grid dimensions
 911 !                ifms, ifme, jfms, jfme,                       &
 912 !                ifts,ifte,jfts,jfte,                          &
 913 !                ir,jr,                                        & ! atm/fire grid ratio
 914 !                zs,                                       & ! atm grid arrays in
 915 !                fp%zsf)                                      ! fire grid arrays out
 916 !        else
 917 !!$OMP CRITICAL(SFIRE_DRIVER_CRIT)
 918 !           write(msg,'(a,i3,a)')'fire_topo_from_atm=',fire_topo_from_atm,' assuming ZSF set, interpolation skipped'
 919 !!$OMP END CRITICAL(SFIRE_DRIVER_CRIT)
 920 !        endif
 921
 922 #ifdef DEBUG_OUT
 923          call write_array_m(its,ite,jts,jte,ims,ime,jms,jme,xlat,'xlat',id)
 924          call write_array_m(its,ite,jts,jte,ims,ime,jms,jme,xlong,'xlong',id)
 925 #endif
 926         if (interpolate_long_lat)then
 927          call message('Interpolating node longitude and latitude to fire mesh')
 928          call interpolate_z2fire(id,1,                 & ! for debug output, <= 0 no output
 929             ids,ide,  jds,jde,                    & ! atm grid dimensions
 930             ims,ime,  jms,jme,                    &
 931             ips,ipe,jps,jpe,                              &
 932             its,ite,jts,jte,                              &
 933             ifds, ifde, jfds, jfde,                       & ! fire grid dimensions
 934             ifms, ifme, jfms, jfme,                       &
 935             ifts,ifte,jfts,jfte,                          &
 936             ir,jr,                                        & ! atm/fire grid ratio
 937             xlat,                                       & ! atm grid arrays in
 938             fxlat)                                      ! fire grid arrays out
 939
 940          call interpolate_z2fire(id,1,                 & ! for debug output, <= 0 no output
 941             ids,ide,  jds,jde,                    & ! atm grid dimensions
 942             ims,ime,  jms,jme,                    &
 943             ips,ipe,jps,jpe,                              &
 944             its,ite,jts,jte,                              &
 945             ifds, ifde, jfds, jfde,                       & ! fire grid dimensions
 946             ifms, ifme, jfms, jfme,                       &
 947             ifts,ifte,jfts,jfte,                          &
 948             ir,jr,                                        & ! atm/fire grid ratio
 949             xlong,                                       & ! atm grid arrays in
 950             fxlong)                                      ! fire grid arrays out
 951         endif !interpolate_long_lat
 952
 953         ! cannot initialize moisture model because T2 Q2 PSFC are not set yet
 954
 955         ! not using zsf any more but anyway
 956         call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,fp%zsf,'driver_phys:zsf')
 957
 958      endif !restart
 959
 960     elseif(ifun.eq.3)then  ! interpolate atmosphere variables to the fire grid
 961
 962       for_balbi: if(ibeh.eq.2)then
 963
 964          call message('ibeh=2: interpolating atmosphere input for Balbi rate of spread model')
 965          if(.not.use_atm_vars)call crash('ibeh=2 Balbi model requires atmosphere but use_atm_vars=false')
 966
 967          ! balbi_msglevel = 4
 968
 969          !! halo on t2, psfc, q2 before
 970
 971          call rho2fire( id,       &
 972              ids,ide, jds,jde,                    & ! atm domain dimensions
 973              ims,ime, jms,jme,                    & ! atm memory dimensions
 974              ips,ipe, jps,jpe,                    & ! atm patch dimensions
 975              its,ite, jts,jte,                    & ! atm tile dimensions
 976              ifds, ifde, jfds, jfde,              & ! fire grid dimensions
 977              ifms, ifme, jfms, jfme,              & ! fire memory dimensions
 978              ifts,ifte,jfts,jfte,                 & ! fire tildemensions
 979              ir,jr,                               & ! atm/fire grid ratio
 980              t2, psfc, q2,                        & ! atm surface variables
 981              rhof)                                  ! surface air density
 982
 983          call interpolate_z2fire(id,1,            & ! for debug output, <= 0 no output
 984             ids,ide,  jds,jde,                    & ! atm grid dimensions
 985             ims,ime,  jms,jme,                    &
 986             ips,ipe,jps,jpe,                              &
 987             its,ite,jts,jte,                              &
 988             ifds, ifde, jfds, jfde,                       & ! fire grid dimensions
 989             ifms, ifme, jfms, jfme,                       &
 990             ifts,ifte,jfts,jfte,                          &
 991             ir,jr,                                        & ! atm/fire grid ratio
 992             t2,                                          & ! atm grid arrays in
 993             tempf)                                      ! fire grid arrays out
 994
 995          !! halo on rhof,tempf after
 996
 997          ffwidth(ifts:ifte,jfts:jfte)     = 100.     ! default for now
 998
 999       end if for_balbi
1000
1001       use_atm: if(use_atm_vars)then
1002
1003         call write_array_m(its,ite,jts,jte,ims,ime,jms,jme,z0,'z0',id)
1004         call write_array_m3(its,ite1,kts,kde-1,jts,jte,ims,ime,kms,kme,jms,jme,u,'u_2',id)
1005         call write_array_m3(its,ite,kts,kde-1,jts,jte1,ims,ime,kms,kme,jms,jme,v,'v_2',id)
1006         call write_array_m3(its,ite,kts,kde,jts,jte,ims,ime,kms,kme,jms,jme,ph,'ph_2',id)
1007         call write_array_m3(its,ite,kts,kde,jts,jte,ims,ime,kms,kme,jms,jme,phb,'phb',id)
1008
1009         wind_log: if(fire_wind_log_interp.eq.4)then
1010           call print_2d_stats(its,ite,jts,jte,ims,ime,jms,jme,z0,'driver_phys:z0')
1011           call interpolate_atm2fire(id,                     & ! flag for debug output
1012             ids,ide, kds,kde, jds,jde,                    & ! atm grid dimensions
1013             ims,ime, kms,kme, jms,jme,                    &
1014             ips,ipe, jps,jpe,                             &
1015             its,ite,jts,jte,                              &
1016             ifds, ifde, jfds, jfde,                       & ! fire grid dimensions
1017             ifms, ifme, jfms, jfme,                       &
1018             ifps, ifpe, jfps, jfpe,                       & ! fire patch bounds
1019             ifts, ifte, jfts, jfte,                       &
1020             ir,jr,                                        & ! atm/fire grid ratio
1021             u_frame, v_frame,                             & ! velocity frame correction
1022             u,v,                                          & ! 3D atm grid arrays in
1023             ph,phb,                                       &
1024             z0,zs,                                        & ! 2D atm grid arrays in
1025             uah,vah,                                      & ! 2D atm grid out
1026             fp%vx,fp%vy)                                    ! fire grid arrays out
1027
1028           call apply_windrf(                        &
1029             ifms,ifme,jfms,jfme,                    &
1030             ifts,ifte,jfts,jfte,                    &
1031             nfuel_cat,fp%vx,fp%vy)
1032
1033             wz0 = 0             !Set WRF roughness length to 0 since we dont use this for option 4.
1034
1035         else wind_log
1036
1037           call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,fz0,'driver_phys:fz0')
1038
1039           call interpolate_wind2fire_height(id,       & ! to identify debugging prints and files if needed
1040             ids,ide, kds,kde, jds,jde,                    & ! atm grid dimensions
1041             ims,ime, kms,kme, jms,jme,                    &
1042             ips,ipe,jps,jpe,                              &
1043             its,ite,jts,jte,                              &
1044             ifds, ifde, jfds, jfde,                       & ! fire grid dimensions
1045             ifms, ifme, jfms, jfme,                       &
1046             ifps, ifpe, jfps, jfpe,                       & ! fire patch bounds
1047             ifts,ifte,jfts,jfte,                          &
1048             ir,jr,                                        & ! atm/fire grid ratio
1049             u_frame, v_frame,                             & ! velocity frame correction
1050             u,v,ph,phb,                                   & ! input atmospheric arrays
1051             fz0,fwh,                                      & ! input fire terrain height and wind height
1052             fp%vx,fp%vy)                                          ! output fire arrays
1053
1054           wind_rf: if(fire_use_windrf.eq.1)then             ! wind reduction factors below fwh
1055             call apply_windrf(                      &
1056             ifms,ifme,jfms,jfme,                    &
1057             ifts,ifte,jfts,jfte,                    &
1058             nfuel_cat,fp%vx,fp%vy)
1059
1060
1061           !If we are deploying the Massman parameterization, we need to interpolate our winds
1062           !down to the canopy height instead of the FWH. Here we will need to call the interpolate_wind2can
1063           !subroutine, for our specified canopy height assigned in the namelist (not yet implemented,
1064           !this is just going to be hard coded for now). Once this is done, we can use the massman
1065           !parameterization to scale the winds down to the FWH.
1066
1067           elseif(fire_use_windrf.eq.4)then wind_rf          ! massman interpolation below fwh
1068
1069 !$OMP CRITICAL(SFIRE_DRIVER_CRIT)
1070             call message('fire_use_windrf option 4 has been selected....!')
1071             write(msg, '(A20,I1)') 'fire_can_top_read = ', fire_can_top_read
1072             call message(msg)
1073 !$OMP END CRITICAL(SFIRE_DRIVER_CRIT)
1074
1075             !Lets set the canopy height, using either a look up table, or just explictly setting it
1076             !using WRF input data from LANDFIRE.
1077
1078             if(fire_can_top_read == 0)then
1079
1080                 can_top = 0                             ! Initialize canopy top array on fire mesh
1081
1082                 call find_trees_fmesh(id,             & ! to identify debugging prints and files if needed
1083                 ids,ide, kds,kde, jds,jde,                    & ! atm grid dimensions
1084                 ims,ime, kms,kme, jms,jme,                    &
1085                 ips,ipe,jps,jpe,                              &
1086                 its,ite,jts,jte,                              &
1087                 ifds, ifde, jfds, jfde,                       & ! fire grid dimensions
1088                 ifms, ifme, jfms, jfme,                       &
1089                 ifps, ifpe, jfps, jfpe,                       & ! fire patch bounds
1090                 ifts,ifte,jfts,jfte,                          &
1091                 ir,jr,                                        & ! atm/fire grid ratio
1092                 can_top, nfuel_cat)                             ! tree canopy height info
1093
1094             end if
1095
1096
1097             !Compute canopy top winds. Save this variable to cuf and cvf, which is like
1098             !vf and uf in terms of dimensionality, except that its canopy top winds.
1099
1100             call interpolate_wind2fire_height(id,         & ! to identify debugging prints and files if needed
1101             ids,ide, kds,kde, jds,jde,                    & ! atm grid dimensions
1102             ims,ime, kms,kme, jms,jme,                    &
1103             ips,ipe,jps,jpe,                              &
1104             its,ite,jts,jte,                              &
1105             ifds, ifde, jfds, jfde,                       & ! fire grid dimensions
1106             ifms, ifme, jfms, jfme,                       &
1107             ifps, ifpe, jfps, jfpe,                       & ! fire patch bounds
1108             ifts,ifte,jfts,jfte,                          &
1109             ir,jr,                                        & ! atm/fire grid ratio
1110             u_frame, v_frame,                             & ! velocity frame correction
1111             u,v,ph,phb,                                   & ! input atmospheric arrays
1112             wz0,can_top,                                  & ! input fire arrays
1113             cuf,cvf)                                        ! output fire arrays
1114
1115
1116             !Lets call our massman parameterization and scale winds winds within 'forests'
1117             !within 'forests' using this scheme. uf and vf values that fall within
1118             !forests are modifed accordingly, which will change the fire ROS.
1119
1120             call massman_fwh(id,                          & ! to identify debugging prints and files if needed
1121             ids,ide, kds,kde, jds,jde,                    & ! atm grid dimensions
1122             ims,ime, kms,kme, jms,jme,                    &
1123             ips,ipe,jps,jpe,                              &
1124             its,ite,jts,jte,                              &
1125             ifds, ifde, jfds, jfde,                       & ! fire grid dimensions
1126             ifms, ifme, jfms, jfme,                       &
1127             ifps, ifpe, jfps, jfpe,                       & ! fire patch bounds
1128             ifts,ifte,jfts,jfte,                          &
1129             ir,jr,                                        & ! atm/fire grid ratio
1130             cuf,cvf,                                      & ! canopy top winds (input)
1131             nfuel_cat,                                    & ! fire fuel category (input)
1132             can_top,                                      & ! canopy top (input)
1133             fp%vx,fp%vy)
1134
1135           end if wind_rf
1136         end if wind_log
1137       end if use_atm ! use_atm_vars
1138
1139     elseif(ifun.eq.4)then
1140
1141       ! interpolate and compute weighted average to get the fuel moisture
1142       !! print *,'ifun=4, run_fuel_moisture=',run_fuel_moisture
1143       if(run_fuel_moisture)then
1144         call message('fuel_moisture start')
1145         call fuel_moisture(                &
1146         id,                                     & ! for prints and maybe file names
1147         nfmc,                                &
1148         ids,ide, jds,jde,               & ! atm grid dimensions
1149         ims,ime, jms,jme,           &
1150         ips,ipe,jps,jpe,                &
1151         its,ite,jts,jte,                     &
1152         ifds, ifde, jfds, jfde,         & ! fire grid dimensions
1153         ifms, ifme, jfms, jfme,     &
1154         ifts,ifte,jfts,jfte,                 &
1155         ir,jr,                                   & ! atm/fire grid ratio
1156         nfuel_cat,                         & ! fuel data
1157         fndwi,                             & ! satellite sensing interpolated on fire grid
1158         fmc_gc,                             & ! moisture contents by class on atmospheric grid
1159         fmc_g,                                & ! weighted fuel moisture contents on fire grid
1160         fmc_gc_f                             & ! moisture contents by class on fire grid
1161         )
1162         call message('fuel_moisture end')
1163       endif
1164
1165       call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,fmc_g,'driver_phys:fmc_g')
1166     endif
1167
1168
1169
1170     call sfire_model (id,ifun,restart,replay,  &
1171         run_fuel_moisture,                      & ! if fuel moisture needs to be updated
1172         ifuelread,nfuel_cat0,                   & ! initialize fuel categories
1173         ifds,ifde,jfds,jfde,                    & ! fire domain dims
1174         ifms,ifme,jfms,jfme,                    & ! fire memory dims
1175         ifps,ifpe,jfps,jfpe,                    &
1176         ifts,ifte,jfts,jfte,                    & ! fire patch dims
1177         time_start,dt,                          & ! time and increment
1178         dxf,dyf,                                & ! fire mesh spacing
1179         ignition,hfx,                           & ! description of ignition lines
1180         fxlong,fxlat,                           & ! fire mesh coordinates
1181         fire_hfx,                               & ! given heat flux
1182         tign_in,                                & ! given igntion time
1183         lfn,lfn_out,tign,fuel_frac,             & ! state: level function, ign time, fuel left
1184         fire_area,                              & ! output: fraction of cell burning
1185         fuel_frac_burnt,                        & ! output: fuel fraction burned in this step
1186         fgrnhfx,fgrnqfx,                        & ! output: heat fluxes
1187         ros,flineint,flineint2,                 & ! diagnostic variables
1188         f_ros0,f_rosx,f_rosy,f_ros,             & ! fire risk spread
1189         f_int,f_lineint,f_lineint2,             & ! fire risk intensities
1190         nfuel_cat,                              & ! fuel data per point
1191         fuel_time,fwh,fz0,                      & ! save derived internal data
1192         fp                                      & ! fire coefficients
1193     )
1194
1195      if(ifun.eq.2)then
1196         call setup_wind_log_interpolation(           &
1197             ids,ide,  jds,jde,                    & ! atm grid dimensions
1198             ims,ime,  jms,jme,                    &
1199             ips,ipe,jps,jpe,                              &
1200             its,ite,jts,jte,                              &
1201             ifds, ifde, jfds, jfde,                       & ! fire grid dimensions
1202             ifms, ifme, jfms, jfme,                       &
1203             ifts,ifte,jfts,jfte,                          &
1204             ir,jr,                                        & ! atm/fire grid ratio
1205             z0,                                           & ! atm grid arrays in
1206             nfuel_cat,                                    & ! fuel data per point
1207             fz0,fwh,wz0)                                    ! fire arrays out
1208
1209
1210     elseif(ifun.eq.6)then
1211
1212
1213     ! populate the rate of spread in the 8 directions
1214     do j=jfts,jfte
1215        do i=ifts,ifte
1216           f_ros11(i,j)=fire_rate_of_spread( dxf*(1-2), dyf*(1-2), i,j,fp)
1217           f_ros12(i,j)=fire_rate_of_spread( dxf*(1-2), dyf*(2-2), i,j,fp)
1218           f_ros13(i,j)=fire_rate_of_spread( dxf*(1-2), dyf*(3-2), i,j,fp)
1219           f_ros21(i,j)=fire_rate_of_spread( dxf*(2-2), dyf*(1-2), i,j,fp)
1220           f_ros23(i,j)=fire_rate_of_spread( dxf*(2-2), dyf*(3-2), i,j,fp)
1221           f_ros31(i,j)=fire_rate_of_spread( dxf*(3-2), dyf*(1-2), i,j,fp)
1222           f_ros32(i,j)=fire_rate_of_spread( dxf*(3-2), dyf*(2-2), i,j,fp)
1223           f_ros33(i,j)=fire_rate_of_spread( dxf*(3-2), dyf*(3-2), i,j,fp)
1224        enddo
1225     enddo
1226     call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,f_ros11,'driver_phys:f_ros11')
1227     call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,f_ros12,'driver_phys:f_ros12')
1228     call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,f_ros13,'driver_phys:f_ros13')
1229     call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,f_ros21,'driver_phys:f_ros21')
1230     call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,f_ros23,'driver_phys:f_ros23')
1231     call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,f_ros31,'driver_phys:f_ros31')
1232     call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,f_ros32,'driver_phys:f_ros32')
1233     call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,f_ros33,'driver_phys:f_ros33')
1234
1235     ! heat fluxes into the atmosphere
1236
1237         call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,ros,'sfire_driver:ros')
1238         call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,fgrnhfx,'sfire_driver:fgrnhfx')
1239         call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,fgrnqfx,'sfire_driver:fgrnqfx')
1240
1241         ! sum the fluxes over atm cells
1242         if(use_atm_vars)then
1243           call sum_2d_cells(        &
1244             ifms,ifme,jfms,jfme,  &
1245             ifts,ifte,jfts,jfte,  &
1246             fuel_frac,              &
1247             ims, ime, jms, jme,   &
1248             its,ite,jts,jte,      &
1249             avg_fuel_frac)
1250           call sum_2d_cells(        &
1251             ifms,ifme,jfms,jfme,  &
1252             ifts,ifte,jfts,jfte,  &
1253             fgrnhfx,              &
1254             ims, ime, jms, jme,   &
1255             its,ite,jts,jte,      &
1256             grnhfx)
1257 !comment out the next call to get results as before commit 55fd92051196b796891b60cb7ec1c4bdb8800078
1258           call sum_2d_cells(        &
1259             ifms,ifme,jfms,jfme,  &
1260             ifts,ifte,jfts,jfte,  &
1261             fgrnqfx,              &
1262             ims, ime, jms, jme,   &
1263             its,ite,jts,jte,      &
1264             grnqfx)
1265
1266 !$OMP CRITICAL(SFIRE_DRIVER_CRIT)
1267           write(msg,'(a,f6.3)')'fire-atmosphere feedback scaling ',fire_atm_feedback
1268 !$OMP end CRITICAL(SFIRE_DRIVER_CRIT)
1269           call message(msg)
1270           s = 1./(ir*jr)
1271           do j=jts,jte
1272             do i=its,ite
1273                 ! scale ground fluxes to get the averages
1274                 avg_fuel_frac(i,j)=avg_fuel_frac(i,j)*s
1275                 grnhfx(i,j)=fire_atm_feedback*grnhfx(i,j)*s
1276                 grnqfx(i,j)=fire_atm_feedback*grnqfx(i,j)*s
1277                 ! we do not have canopy fluxes yet...
1278                 canhfx(i,j)=0
1279                 canqfx(i,j)=0
1280             enddo
1281           enddo
1282
1283           call print_2d_stats(its,ite,jts,jte,ims,ime,jms,jme,grnhfx,'fire_driver:grnhfx')
1284           call print_2d_stats(its,ite,jts,jte,ims,ime,jms,jme,grnqfx,'fire_driver:grnqfx')
1285        endif
1286
1287     endif ! ifun=6
1288   endif
1289
1290 enddo ! tiles
1291 !$OMP END PARALLEL DO
1292
1293 #ifdef DEBUG_OUT
1294 if(ifun.eq.1)then
1295     if(pid.ne.0)then
1296         call write_array_m(ips,ipe,jps,jpe,ims,ime,jms,jme,zs,'zs',pid)
1297         call write_array_m(ifps,ifpe,jfps,jfpe,ifms,ifme,jfms,jfme,fp%zsf,'zsf',pid)
1298     endif
1299 endif
1300 #endif
1301
1302 if(ifun.eq.2)then
1303 ! ****      THIS REALLY SHOULD BE DONE ONCE NOT EVERY TIMESTEP
1304 !           we cannot store the ignition structure in grid
1305 !           so we set global ignition instead
1306 !           and we need to have fxlong fxlat for the corners (for checking)
1307 !           which requires MPI calls
1308 !           must be run after the parallel do because of all that
1309             call message('Processing ignition lines')
1310             call fire_ignition_convert (config_flags,ignition,          &
1311                 fxlong, fxlat,                                          &
1312                 ifds,ifde, jfds,jfde,                                   &
1313                 ifms,ifme, jfms,jfme,                                   &
1314                 ifps,ifpe, jfps,jfpe )
1315             call fire_hfx_convert (config_flags,hfx) ! do we even use heat flux lines any more?
1316 endif
1317
1318 if (ifun.eq.3)then
1319  call print_3d_stats_by_slice(ips,ipe,1,moisture_classes,jps,jpe,ims,ime,1,nfmc,jms,jme,fmc_gc,'fmc_gc')
1320  call print_chsum(itimestep,ims,ime,1,nfmc,jms,jme,ids,ide,1,moisture_classes,jds,jde,ips,ipe,1,moisture_classes,jps,jpe,0,0,0,fmc_gc,'fmc_gc')
1321 endif
1322
1323 if (ifun.eq.4)then
1324
1325  call print_chsum(itimestep,ifms,ifme,1,1,jfms,jfme,ifds,ifde,1,1,jfds,jfde,ifps,ifpe,1,1,jfps,jfpe,0,0,0,fmc_g,'fmc_g')
1326  !call print_chsum(itimestep,ims,ime,1,1,jms,jme,ids,ide,1,1,jds,jde,ips,ipe,1,1,jps,jpe,1,0,0,uah,'uah')
1327  !call print_chsum(itimestep,ims,ime,1,1,jms,jme,ids,ide,1,1,jds,jde,ips,ipe,1,1,jps,jpe,0,0,1,vah,'vah')
1328  call print_chsum(itimestep,ifms,ifme,1,1,jfms,jfme,ifds,ifde,1,1,jfds,jfde,ifps,ifpe,1,1,jfps,jfpe,0,0,0,fp%vx,'uf')
1329  call print_chsum(itimestep,ifms,ifme,1,1,jfms,jfme,ifds,ifde,1,1,jfds,jfde,ifps,ifpe,1,1,jfps,jfpe,0,0,0,fp%vy,'vf')
1330  call print_chsum(itimestep,ifms,ifme,1,1,jfms,jfme,ifds,ifde,1,1,jfds,jfde,ifps,ifpe,1,1,jfps,jfpe,0,0,0,lfn,'lfn')
1331  call print_chsum(itimestep,ifms,ifme,1,1,jfms,jfme,ifds,ifde,1,1,jfds,jfde,ifps,ifpe,1,1,jfps,jfpe,0,0,0,tign,'tign')
1332
1333 #ifdef DEBUG_OUT
1334     if(pid.gt.0)then
1335  !       call write_array_m(ips,ipe1,jps,jpe,ims,ime,jms,jme,uah,'uah',pid)
1336  !       call write_array_m(ips,ipe,jps,jpe1,ims,ime,jms,jme,vah,'vah',pid)
1337         call write_array_m(ips,ipe,jps,jpe,ims,ime,jms,jme,grnhfx,'grnhfx',pid)
1338         call write_array_m(ips,ipe,jps,jpe,ims,ime,jms,jme,grnqfx,'grnqfx',pid)
1339         call write_array_m3(ips,ipe1,kds,kde+1,jps,jpe,ims,ime,kms,kme,jms,jme,u,'u',pid)
1340         call write_array_m3(ips,ipe,kds,kde+1,jps,jpe1,ims,ime,kms,kme,jms,jme,v,'v',pid)
1341         call write_array_m(ifps,ifpe,jfps,jfpe,ifms,ifme,jfms,jfme,fp%vx,'uf',pid)
1342         call write_array_m(ifps,ifpe,jfps,jfpe,ifms,ifme,jfms,jfme,fp%vy,'vf',pid)
1343     endif
1344 #endif
1345 endif
1346
1347 if(ifun.eq.5)then
1348     call print_chsum(itimestep,ifms,ifme,1,1,jfms,jfme,ifds,ifde,1,1,jfds,jfde,ifps,ifpe,1,1,jfps,jfpe,0,0,0,lfn,'lfn')
1349     call print_chsum(itimestep,ifms,ifme,1,1,jfms,jfme,ifds,ifde,1,1,jfds,jfde,ifps,ifpe,1,1,jfps,jfpe,0,0,0,tign,'tign')
1350 #ifdef DEBUG_OUT
1351     if(pid.gt.0)then
1352         call write_array_m(ifps,ifpe,jfps,jfpe,ifms,ifme,jfms,jfme,lfn,'lfn',pid)
1353         call write_array_m(ifps,ifpe,jfps,jfpe,ifms,ifme,jfms,jfme,tign,'tign',pid)
1354     endif
1355 #endif
1356 endif
1357
1358 if(ifun.eq.6)then
1359     call print_chsum(itimestep,ifms,ifme,1,1,jfms,jfme,ifds,ifde,1,1,jfds,jfde,ifps,ifpe,1,1,jfps,jfpe,0,0,0,fgrnhfx,'fgrnhfx')
1360     call print_chsum(itimestep,ifms,ifme,1,1,jfms,jfme,ifds,ifde,1,1,jfds,jfde,ifps,ifpe,1,1,jfps,jfpe,0,0,0,fgrnqfx,'fgrnqfx')
1361     call print_chsum(itimestep,ims,ime,1,1,jms,jme,ids,ide,1,1,jds,jde,ips,ipe,1,1,jps,jpe,0,0,0,grnhfx,'grnhfx')
1362     call print_chsum(itimestep,ims,ime,1,1,jms,jme,ids,ide,1,1,jds,jde,ips,ipe,1,1,jps,jpe,0,0,0,grnqfx,'grnqfx')
1363 #ifdef DEBUG_OUT
1364     if(pid.gt.0)then
1365         call write_array_m(ips,ipe,jps,jpe,ims,ime,jms,jme,grnhfx,'grnhfx',pid)
1366         call write_array_m(ips,ipe,jps,jpe,ims,ime,jms,jme,grnqfx,'grnqfx',pid)
1367         call write_array_m(ifps,ifpe,jfps,jfpe,ifms,ifme,jfms,jfme,fuel_frac,'fuel_frac',pid)
1368         call write_array_m(ifps,ifpe,jfps,jfpe,ifms,ifme,jfms,jfme,fgrnhfx,'fgrnhfx',pid)
1369         call write_array_m(ifps,ifpe,jfps,jfpe,ifms,ifme,jfms,jfme,fgrnqfx,'fgrnqfx',pid)
1370     endif
1371 #endif
1372 endif
1373
1374 end subroutine sfire_driver_phys
1375
1376 !
1377 !***
1378 !
1379
1380 subroutine check_fmesh(ids,ide,ifds,ifde,ir,s)
1381 !*** purpose: check if fire and atm meshes line up
1382 implicit none
1383 !*** arguments
1384 integer, intent(in)::ids,ide,ifds,ifde,ir
1385 character(len=*),intent(in)::s
1386 !*** local
1387 character(len=128)msg
1388 !*** executable
1389 if ((ide-ids+1)*ir.ne.(ifde-ifds+1))then
1390 !$OMP CRITICAL(SFIRE_DRIVER_CRIT)
1391     write(msg,1)s,ids,ide,ifds,ifde,ir
1392 1   format('module_fr_sfire_driver: incompatible bounds ',a,' atm ',i5,':',i5,' fire ',i5,':',i5,' ratio ',i3)
1393 !$OMP END CRITICAL(SFIRE_DRIVER_CRIT)
1394     call crash(msg)
1395 endif
1396 end subroutine check_fmesh
1397
1398 subroutine fire_ignition_convert (config_flags,ignition,             &
1399              fxlong, fxlat,                                          &
1400              ifds,ifde, jfds,jfde,                                   &
1401              ifms,ifme, jfms,jfme,                                   &
1402              ifps,ifpe, jfps,jfpe )
1403     implicit none
1404 ! create ignition arrays from scalar flags
1405 !*** arguments
1406     TYPE (grid_config_rec_type) , INTENT(IN)          :: config_flags
1407     TYPE (lines_type), INTENT(OUT):: ignition ! any values from input discarded
1408     integer::ifds,ifde, jfds,jfde,                                   &
1409              ifms,ifme, jfms,jfme,                                   &
1410              ifps,ifpe, jfps,jfpe
1411     real, dimension(ifms:ifme,jfms:jfme):: fxlong,fxlat
1412 !*** local
1413     integer::i,j,ii,jj
1414     logical:: real,ideal
1415     character(len=128)msg
1416     real:: corner_longlat(2,2,2), corner_longlat_1(8), corner_longlat_2(8),lon(2),lat(2)
1417     real, dimension(2,2):: corner_long,corner_lat ! coordinates of fire mesh corner cells
1418
1419 !*** executable
1420
1421
1422     ignition%max_lines=5 ! number of lines that have entries in the namelist
1423     ignition%num_lines=config_flags%fire_num_ignitions
1424
1425     ! this is only until I figure out how to input arrays through the namelist...
1426     if(fire_max_lines.lt.ignition%max_lines)call crash('fire_max_lines too small')
1427
1428     ! figure out which kind of coordinates from the first given
1429     ideal=config_flags%fire_ignition_start_x1 .ne.0. .or. config_flags%fire_ignition_start_y1 .ne. 0.
1430     real=config_flags%fire_ignition_start_lon1 .ne. 0. .or. config_flags%fire_ignition_start_lat1 .ne. 0.
1431     if(ideal)call message('Using ideal ignition coordinates, m from the lower left domain corner')
1432     if(real)call message('Using real ignition coordinates, longitude and latitude')
1433     if(ideal.and.real)call crash('Only one of the ideal or real coordinates may be given')
1434
1435     ignition%longlat=0  ! default, if no ignition
1436     if(ideal)then
1437        ! use values from _x and _y variables
1438        ignition%longlat=0
1439        ignition%line(1)%start_x=config_flags%fire_ignition_start_x1
1440        ignition%line(1)%start_y=config_flags%fire_ignition_start_y1
1441        ignition%line(1)%end_x=config_flags%fire_ignition_end_x1
1442        ignition%line(1)%end_y=config_flags%fire_ignition_end_y1
1443        ignition%line(2)%start_x=config_flags%fire_ignition_start_x2
1444        ignition%line(2)%start_y=config_flags%fire_ignition_start_y2
1445        ignition%line(2)%end_x=config_flags%fire_ignition_end_x2
1446        ignition%line(2)%end_y=config_flags%fire_ignition_end_y2
1447        ignition%line(3)%start_x=config_flags%fire_ignition_start_x3
1448        ignition%line(3)%start_y=config_flags%fire_ignition_start_y3
1449        ignition%line(3)%end_x=config_flags%fire_ignition_end_x3
1450        ignition%line(3)%end_y=config_flags%fire_ignition_end_y3
1451        ignition%line(4)%start_x=config_flags%fire_ignition_start_x4
1452        ignition%line(4)%start_y=config_flags%fire_ignition_start_y4
1453        ignition%line(4)%end_x=config_flags%fire_ignition_end_x4
1454        ignition%line(4)%end_y=config_flags%fire_ignition_end_y4
1455        ignition%line(5)%start_x=config_flags%fire_ignition_start_x5
1456        ignition%line(5)%start_y=config_flags%fire_ignition_start_y5
1457        ignition%line(5)%end_x=config_flags%fire_ignition_end_x5
1458        ignition%line(5)%end_y=config_flags%fire_ignition_end_y5
1459     endif
1460     if(real)then
1461         ! use values from _long and _lat
1462        ignition%longlat=1
1463        ignition%line(1)%start_x=config_flags%fire_ignition_start_lon1
1464        ignition%line(1)%start_y=config_flags%fire_ignition_start_lat1
1465        ignition%line(1)%end_x=config_flags%fire_ignition_end_lon1
1466        ignition%line(1)%end_y=config_flags%fire_ignition_end_lat1
1467        ignition%line(2)%start_x=config_flags%fire_ignition_start_lon2
1468        ignition%line(2)%start_y=config_flags%fire_ignition_start_lat2
1469        ignition%line(2)%end_x=config_flags%fire_ignition_end_lon2
1470        ignition%line(2)%end_y=config_flags%fire_ignition_end_lat2
1471        ignition%line(3)%start_x=config_flags%fire_ignition_start_lon3
1472        ignition%line(3)%start_y=config_flags%fire_ignition_start_lat3
1473        ignition%line(3)%end_x=config_flags%fire_ignition_end_lon3
1474        ignition%line(3)%end_y=config_flags%fire_ignition_end_lat3
1475        ignition%line(4)%start_x=config_flags%fire_ignition_start_lon4
1476        ignition%line(4)%start_y=config_flags%fire_ignition_start_lat4
1477        ignition%line(4)%end_x=config_flags%fire_ignition_end_lon4
1478        ignition%line(4)%end_y=config_flags%fire_ignition_end_lat4
1479        ignition%line(5)%start_x=config_flags%fire_ignition_start_lon5
1480        ignition%line(5)%start_y=config_flags%fire_ignition_start_lat5
1481        ignition%line(5)%end_x=config_flags%fire_ignition_end_lon5
1482        ignition%line(5)%end_y=config_flags%fire_ignition_end_lat5
1483     endif
1484     ! common to both cases
1485        ignition%line(1)%ros=config_flags%fire_ignition_ros1
1486        ignition%line(1)%radius=config_flags%fire_ignition_radius1
1487        ignition%line(1)%start_time=config_flags%fire_ignition_start_time1
1488        ignition%line(1)%end_time=config_flags%fire_ignition_end_time1
1489        ignition%line(2)%ros=config_flags%fire_ignition_ros2
1490        ignition%line(2)%radius=config_flags%fire_ignition_radius2
1491        ignition%line(2)%start_time=config_flags%fire_ignition_start_time2
1492        ignition%line(2)%end_time=config_flags%fire_ignition_end_time2
1493        ignition%line(3)%ros=config_flags%fire_ignition_ros3
1494        ignition%line(3)%radius=config_flags%fire_ignition_radius3
1495        ignition%line(3)%start_time=config_flags%fire_ignition_start_time3
1496        ignition%line(3)%end_time=config_flags%fire_ignition_end_time3
1497        ignition%line(4)%ros=config_flags%fire_ignition_ros4
1498        ignition%line(4)%radius=config_flags%fire_ignition_radius4
1499        ignition%line(4)%start_time=config_flags%fire_ignition_start_time4
1500        ignition%line(4)%end_time=config_flags%fire_ignition_end_time4
1501        ignition%line(5)%ros=config_flags%fire_ignition_ros5
1502        ignition%line(5)%radius=config_flags%fire_ignition_radius5
1503        ignition%line(5)%start_time=config_flags%fire_ignition_start_time5
1504        ignition%line(5)%end_time=config_flags%fire_ignition_end_time5
1505
1506        call postprocess_lines(ignition,'ros',config_flags)
1507
1508 !       get the coordinates of the corner cells
1509         corner_longlat=0.
1510         if(ifds.eq.ifps.and.jfds.eq.jfps)then
1511             corner_longlat(1,1,1)=fxlong(ifps,jfps)
1512             corner_longlat(1,1,2)=fxlat(ifps,jfps)
1513         endif
1514         if(ifds.eq.ifps.and.jfde.eq.jfpe)then
1515             corner_longlat(1,2,1)=fxlong(ifps,jfpe)
1516             corner_longlat(1,2,2)=fxlat(ifps,jfpe)
1517         endif
1518         if(ifde.eq.ifpe.and.jfds.eq.jfps)then
1519             corner_longlat(2,1,1)=fxlong(ifpe,jfps)
1520             corner_longlat(2,1,2)=fxlat(ifpe,jfps)
1521         endif
1522         if(ifde.eq.ifpe.and.jfde.eq.jfpe)then
1523             corner_longlat(2,2,1)=fxlong(ifpe,jfpe)
1524             corner_longlat(2,2,2)=fxlat(ifpe,jfpe)
1525         endif
1526         corner_longlat_1=reshape(corner_longlat,(/8/))
1527 #ifdef DM_PARALLEL
1528         call wrf_dm_sum_reals(corner_longlat_1,corner_longlat_2)
1529 #else
1530         corner_longlat_2=corner_longlat_1
1531 #endif
1532         corner_longlat=reshape(corner_longlat_2,(/2,2,2/))
1533         corner_long=corner_longlat(1:2,1:2,1)
1534         corner_lat=corner_longlat(1:2,1:2,2)
1535         if(fire_print_msg.ge.2)then
1536             do i=1,2
1537                 do j=1,2
1538                      write(msg,'(a,2i2,a,2f14.8)')'corner',i,j,' coordinates ',corner_long(i,j),corner_lat(i,j)
1539                      call message(msg)
1540                 enddo
1541             enddo
1542         endif
1543         lon(1)=(corner_long(1,1)+corner_long(1,2))/2.
1544         lon(2)=(corner_long(2,1)+corner_long(2,2))/2.
1545         lat(1)=(corner_lat(1,1)+corner_lat(2,1))/2.
1546         lat(2)=(corner_lat(1,2)+corner_lat(2,2))/2.
1547         if(fire_print_msg.ge.2)then
1548             write(msg,'(4(a,f14.8))')'coordinates ',lon(1),':',lon(2),',',lat(1),':',lat(2)
1549             call message(msg)
1550         endif
1551
1552        do i=1,ignition%num_lines
1553            call check_ignition_coordinate(ignition%line(i)%start_x,lon(1),lon(2))
1554            call check_ignition_coordinate(ignition%line(i)%start_y,lat(1),lat(2))
1555            call check_ignition_coordinate(ignition%line(i)%end_x,lon(1),lon(2))
1556            call check_ignition_coordinate(ignition%line(i)%end_y,lat(1),lat(2))
1557        enddo
1558
1559        if (fire_ignition_clamp>0) then
1560        do i=1,ignition%num_lines
1561            call clamp_to_grid(ignition%line(i)%start_x,lon(1),lon(2),ifds,ifde,ignition%line(i)%start_x,ii)
1562            call clamp_to_grid(ignition%line(i)%start_y,lat(1),lat(2),jfds,jfde,ignition%line(i)%start_y,jj)
1563            call display_clamp
1564            call clamp_to_grid(ignition%line(i)%end_x,lon(1),lon(2),ifds,ifde,ignition%line(i)%end_x,ii)
1565            call clamp_to_grid(ignition%line(i)%end_y,lat(1),lat(2),jfds,jfde,ignition%line(i)%end_y,jj)
1566            call display_clamp
1567            ! for now, ii jj ignored. In future replace by fxlong(ii,jj), fxlat(ii,jj) to guard against rounding
1568        enddo
1569        endif
1570        contains
1571        subroutine display_clamp
1572         character(len=128)::msg
1573         real::d1,d2
1574            if(ii>=ifps.and.ii<=ifpe.and.jj>=jfps.and.jj<=jfpe)then
1575                write(msg,'(a,2f14.8,a,2i6)')'grid node ',fxlong(ii,jj),fxlat(ii,jj),' index',ii,jj
1576                call message(msg)
1577            endif
1578        end subroutine display_clamp
1579     end subroutine fire_ignition_convert
1580
1581      subroutine check_ignition_coordinate(x,x1,x2)
1582 !***    arguments
1583         real, intent(in)::x,x1,x2
1584         character(len=128)::msg
1585         if (.not.(x>x1 .and. x<x2))then
1586             write(msg,'(a,f14.8,a,2f14.8)')'ignition point coordinate',x,' must be inside the bounds',x1,x2
1587             call crash(msg)
1588         endif
1589      end subroutine check_ignition_coordinate
1590
1591      subroutine clamp_to_grid(x,x1,x2,i1,i2,xout,iout)
1592 !***    round point on uniform mesh to a mesh node
1593 !***    arguments
1594         real, intent(in)::x,x1,x2
1595         integer, intent(in)::i1,i2
1596         real, intent(out)::xout
1597         integer, intent(out)::iout
1598 !***    local
1599         character(len=128)::msg
1600         integer:: i
1601         real::r,dx,xr
1602 !***    executable
1603         dx=(x2-x1)/(i2-i1)
1604         r=i1+(x-x1)/dx
1605         iout=nint(r)
1606         xr=x1+(iout-i1)*dx
1607         if(fire_print_msg.ge.2)then
1608             write(msg,'(a,f14.8,a,f14.8,a,i6)')'coordinate ',x,' clamped to ',xr,' index',iout
1609             call message(msg)
1610         endif
1611         xout=xr
1612     end subroutine clamp_to_grid
1613 !
1614 !***
1615 !
1616      subroutine postprocess_lines(lines,value_name,config_flags)
1617 !      count lines, fill endpoints, print stats
1618 !***   arguments
1619         type(lines_type), intent(inout)::lines
1620         character(len=3), intent(in)::value_name
1621         TYPE (grid_config_rec_type) , INTENT(IN)  :: config_flags ! namelist
1622 !***   local
1623      integer::i,n
1624      real::value
1625      real::lat_ctr,lon_ctr
1626      character(len=128)msg,f2,f3
1627 !***   executable
1628
1629         n=lines%num_lines
1630         do i=1,n
1631             ! find the last line that has positive radius and reset num_lines if needed
1632             if(lines%line(i)%radius.gt.0.)lines%num_lines=i
1633             ! expand ignition data given as zero
1634             if(lines%line(i)%end_x.eq.0.)lines%line(i)%end_x=lines%line(i)%start_x
1635             if(lines%line(i)%end_y.eq.0.)lines%line(i)%end_y=lines%line(i)%start_y
1636             if(lines%line(i)%end_time.eq.0.)lines%line(i)%end_time=lines%line(i)%start_time
1637         enddo
1638
1639     if(lines%longlat .eq. 0)then
1640        ! ideal
1641        !  ignition is in m
1642        lines%unit_fxlong=1.
1643        lines%unit_fxlat=1.
1644        ! will set fire mesh coordinates to uniform mesh below
1645     else
1646        ! real
1647        lat_ctr=config_flags%cen_lat
1648        lon_ctr=config_flags%cen_lon
1649        ! 1 degree in m (approximate OK)
1650        lines%unit_fxlat=pi2/(360.*reradius)  ! earth circumference in m / 360 degrees
1651        lines%unit_fxlong=cos(lat_ctr*pi2/360.)*lines%unit_fxlat  ! latitude
1652     endif
1653
1654     if(fire_print_msg.ge.2)then
1655 !$OMP CRITICAL(SFIRE_DRIVER_CRIT)
1656        if(lines%longlat .eq. 0)then
1657          write(msg,1)'start x','start y','end x','end y','start t','end t',value_name,'radius'
1658        else
1659          write(msg,1)'start lon','start lat','end lon','end lat','start time','end time',value_name,'radius'
1660        endif
1661 1      format(4a10,4a9)
1662        call message(msg)
1663        do i=1,lines%num_lines
1664          select case (value_name)
1665               case ('ros')
1666                   value = lines%line(i)%ros
1667                   f2='(4f10.3,4f9.2)'
1668                   f3='(4f10.5,4f9.2)'
1669               case ('hfx')
1670                   value = lines%line(i)%hfx_value
1671                   f2='(4f10.3,2f9.2,e9.2,f9.2)'
1672                   f3='(4f10.5,2f9.2,e9.2,f9.2)'
1673               case default
1674                   call crash('postprocess_lines: bad value_name '//value_name)
1675          end select
1676          if(lines%longlat .eq. 0)then
1677            write(msg,f2)lines%line(i)%start_x, lines%line(i)%start_y, &
1678                       lines%line(i)%end_x,  lines%line(i)%end_y,        &
1679                       lines%line(i)%start_time, lines%line(i)%end_time, &
1680                       value, lines%line(i)%radius
1681          else
1682            write(msg,f3)lines%line(i)%start_x, lines%line(i)%start_y, &
1683                       lines%line(i)%end_x,  lines%line(i)%end_y,        &
1684                       lines%line(i)%start_time, lines%line(i)%end_time, &
1685                       value, lines%line(i)%radius
1686          endif
1687          call message(msg)
1688          if(lines%line(i)%start_time > lines%line(i)%end_time)then
1689               call crash('start time may not be after end time')
1690          endif
1691        enddo
1692 !$OMP END CRITICAL(SFIRE_DRIVER_CRIT)
1693      endif
1694 end subroutine postprocess_lines
1695
1696
1697 subroutine fire_hfx_convert (config_flags,hfx)
1698     implicit none
1699 ! create heat flux line(s) from scalar flags
1700 !*** arguments
1701     TYPE (grid_config_rec_type) , INTENT(IN)          :: config_flags
1702     TYPE (lines_type), INTENT(OUT):: hfx ! any values from input discarded
1703 !*** local
1704     integer::i
1705     logical:: real,ideal
1706     real::lat_ctr,lon_ctr
1707     character(len=128)msg
1708 !*** executable
1709     ! this is only until I figure out how to input arrays through the namelist...
1710     hfx%num_lines=config_flags%fire_hfx_num_lines
1711     if(fire_max_lines.lt.hfx%num_lines)call crash('fire_max_lines too small')
1712
1713     ! figure out which kind of coordinates from the first given
1714     ideal=config_flags%fire_hfx_start_x1 .ne.0. .or. config_flags%fire_hfx_start_y1 .ne. 0.
1715     real=config_flags%fire_hfx_start_lon1 .ne. 0. .or. config_flags%fire_hfx_start_lat1 .ne. 0.
1716     if(ideal)call message('Using ideal heat flux line coordinates, m from the lower left domain corner')
1717     if(real)call message('Using real heat flux line coordinates, longitude and latitude')
1718     if(ideal.and.real)call crash('Only one of the ideal or real coordinates may be given')
1719
1720     hfx%longlat=0  ! default, if no ignition
1721     if(ideal)then
1722        ! use values from _x and _y variables
1723        hfx%longlat=0
1724        hfx%line(1)%start_x=config_flags%fire_hfx_start_x1
1725        hfx%line(1)%start_y=config_flags%fire_hfx_start_y1
1726        hfx%line(1)%end_x=config_flags%fire_hfx_end_x1
1727        hfx%line(1)%end_y=config_flags%fire_hfx_end_y1
1728     endif
1729     if(real)then
1730         ! use values from _long and _lat
1731        hfx%longlat=1
1732        hfx%line(1)%start_x=config_flags%fire_hfx_start_lon1
1733        hfx%line(1)%start_y=config_flags%fire_hfx_start_lat1
1734        hfx%line(1)%end_x=config_flags%fire_hfx_end_lon1
1735        hfx%line(1)%end_y=config_flags%fire_hfx_end_lat1
1736     endif
1737     ! common to both cases
1738        hfx%line(1)%radius=config_flags%fire_hfx_radius1
1739        hfx%line(1)%start_time=config_flags%fire_hfx_start_time1
1740        hfx%line(1)%end_time=config_flags%fire_hfx_end_time1
1741        hfx%line(1)%trans_time=config_flags%fire_hfx_trans_time1
1742        hfx%line(1)%hfx_value=config_flags%fire_hfx_value1
1743
1744        call postprocess_lines(hfx,'hfx',config_flags)
1745
1746 end subroutine fire_hfx_convert
1747
1748 subroutine set_flags(config_flags)
1749 USE module_configure
1750 implicit none
1751 TYPE (grid_config_rec_type) , INTENT(IN)          :: config_flags
1752 ! copy flags from wrf to module_fr_sfire_util
1753 ! for instructions how to add a flag see the top of module_fr_sfire_util.F
1754
1755 fire_perimeter_time     = config_flags%fire_perimeter_time
1756 fire_tign_in_time       = config_flags%fire_tign_in_time
1757 fire_print_msg          = config_flags%fire_print_msg
1758 fire_print_file         = config_flags%fire_print_file
1759 fuel_left_method        = config_flags%fire_fuel_left_method
1760 fuel_left_irl           = config_flags%fire_fuel_left_irl
1761 fuel_left_jrl           = config_flags%fire_fuel_left_jrl
1762 fire_atm_feedback       = config_flags%fire_atm_feedback
1763 fire_hfx_given          = config_flags%fire_hfx_given
1764 fire_hfx_num_lines      = config_flags%fire_hfx_num_lines
1765 fire_hfx_latent_part    = config_flags%fire_hfx_latent_part
1766 fire_update_fuel_frac   = config_flags%fire_update_fuel_frac
1767 boundary_guard          = config_flags%fire_boundary_guard
1768 fire_back_weight        = config_flags%fire_back_weight
1769 fire_grows_only         = config_flags%fire_grows_only
1770 fire_upwinding          = config_flags%fire_upwinding
1771 if(fire_upwinding .eq.9)fire_upwinding =3  ! replace by old default
1772 fire_viscosity          = config_flags%fire_viscosity
1773 fire_lfn_ext_up         = config_flags%fire_lfn_ext_up
1774 fire_test_steps         = config_flags%fire_test_steps
1775 !fire_topo_from_atm     = config_flags%fire_topo_from_atm
1776 fire_advection          = config_flags%fire_advection
1777 fire_wind_log_interp    = config_flags%fire_wind_log_interp
1778 fire_use_windrf         = config_flags%fire_use_windrf
1779 fire_fmc_read           = config_flags%fire_fmc_read
1780 fire_ignition_clamp     = config_flags%fire_ignition_clamp
1781 kfmc_ndwi               = config_flags%kfmc_ndwi
1782 fndwi_from_ndwi         = config_flags%fndwi_from_ndwi
1783 fire_can_top_read       = config_flags%fire_can_top_read
1784 end subroutine set_flags
1785
1786 !
1787 !*****************************
1788 !
1789
1790
1791 subroutine set_fp_from_grid(grid,fp)
1792     implicit none
1793     type(domain),intent(in)::grid
1794     type(fire_params),intent(out)::fp
1795
1796     ! pointers to be passed to fire spread formulas
1797     fp%vx => grid%uf         ! fire winds
1798     fp%vy => grid%vf         ! fire winds
1799     fp%zsf => grid%zsf       ! terrain height
1800     fp%dzdxf => grid%dzdxf   ! terrain grad
1801     fp%dzdyf => grid%dzdyf   ! terrain grad
1802     fp%bbb => grid%bbb       ! spread formula coeff
1803     fp%phisc => grid%phisc   ! spread formula coeff
1804     fp%phiwc => grid%phiwc   ! spread formula coeff
1805     fp%r_0 => grid%r_0       ! spread formula coeff
1806     fp%fgip => grid%fgip     ! spread formula coeff
1807     fp%ischap => grid%ischap ! spread formula coeff
1808     fp%fuel_time => grid%fuel_time ! time for fuel to burn to 1/e
1809     fp%fmc_g => grid%fmc_g   ! fuel moisture, ground
1810     fp%fmc_gc_f => grid%fmc_gc_f   ! fuel moisture, ground by class
1811     fp%nfuel_cat => grid%nfuel_cat ! fuel category
1812     fp%tempf => grid%tempf     ! temperature
1813     fp%rhof => grid%rhof     ! surface air density
1814     fp%ffwidth => grid%ffwidth     ! fire front width
1815
1816 end subroutine set_fp_from_grid
1817
1818 subroutine check_grid_alloc(grid,config_flags)
1819
1820 ! check if needed arrays in grid are fully allocated
1821 ! they may not be if they are missing in the registry packages
1822 ! leading to memory corruption which is hard to debug
1823
1824     TYPE(domain) , TARGET, INTENT(IN) :: grid                             ! state
1825     TYPE (grid_config_rec_type) , INTENT(IN)  :: config_flags ! namelist
1826
1827     call message('check_grid_alloc: checking if grid arrays used sfire_driver_phys are allocated',level=2)
1828     call check_size(size(grid%u_2),'u_2')
1829     call check_size(size(grid%v_2),'v_2')
1830     call check_size(size(grid%ph_2),'ph_2')
1831     call check_size(size(grid%phb),'phb')
1832     call check_size(size(grid%z0),'z0')
1833     call check_size(size(grid%xlong),'xlong')
1834     call check_size(size(grid%xlat),'xlat')
1835     call check_size(size(grid%tign_in),'tign_in')
1836     call check_size(size(grid%lfn),'lfn')
1837     call check_size(size(grid%tign_g),'tign_g')
1838     call check_size(size(grid%fire_area),'fire_area')
1839     call check_size(size(grid%fuel_frac),'fuel_frac')
1840     call check_size(size(grid%fuel_frac_burnt),'fuel_frac_burnt')
1841     call check_size(size(grid%lfn_out),'lfn_out')
1842     call check_size(size(grid%avg_fuel_frac),'avg_fuel_frac')
1843     call check_size(size(grid%grnhfx),'grnhfx')
1844     call check_size(size(grid%grnqfx),'grnqfx')
1845     call check_size(size(grid%canhfx),'canhfx')
1846     call check_size(size(grid%canqfx),'canqfx')
1847     call check_size(size(grid%fgrnhfx),'fgrnhfx')
1848     call check_size(size(grid%fgrnqfx),'fgrnqfx')
1849     call check_size(size(grid%fcanhfx),'fcanhfx')
1850     call check_size(size(grid%ros),'ros')
1851     call check_size(size(grid%flineint),'flineint')
1852     call check_size(size(grid%flineint2),'flineint2')
1853     call check_size(size(grid%f_ros0),'f_ros0')
1854     call check_size(size(grid%f_rosx),'f_rosx')
1855     call check_size(size(grid%f_rosy),'f_rosy')
1856     call check_size(size(grid%f_ros),'f_ros')
1857     call check_size(size(grid%f_int),'f_int')
1858     call check_size(size(grid%f_lineint),'f_lineint')
1859     call check_size(size(grid%f_lineint2),'f_lineint2')
1860     call check_size(size(grid%f_ros11),'f_ros11')
1861     call check_size(size(grid%f_ros12),'f_ros12')
1862     call check_size(size(grid%f_ros13),'f_ros13')
1863     call check_size(size(grid%f_ros21),'f_ros21')
1864     call check_size(size(grid%f_ros23),'f_ros23')
1865     call check_size(size(grid%f_ros31),'f_ros31')
1866     call check_size(size(grid%f_ros32),'f_ros32')
1867     call check_size(size(grid%f_ros33),'f_ros33')
1868     call check_size(size(grid%fxlong),'fxlong')
1869     call check_size(size(grid%fxlat),'fxlat')
1870     call check_size(size(grid%fire_hfx),'fire_hfx')
1871     call check_size(size(grid%nfuel_cat),'nfuel_cat')
1872     call check_size(size(grid%fuel_time),'fuel_time')
1873     call check_size(size(grid%fz0),'fz0')
1874     call check_size(size(grid%fwh),'fwh')
1875     call check_size(size(grid%uah),'uah')
1876     call check_size(size(grid%vah),'vah')
1877     call check_size(size(grid%rainc),'rainc')
1878     call check_size(size(grid%rainnc),'rainnc')
1879     call check_size(size(grid%t2),'t2')
1880     call check_size(size(grid%q2),'q2')
1881     call check_size(size(grid%psfc),'psfc')
1882     call check_size(size(grid%ndwi),'ndwi')
1883     call check_size(size(grid%fndwi),'fndwi')
1884     call message('check_grid_alloc: checking if grid arrays used in fire params are allocated',level=2)
1885     call check_size(size(grid%uf),'uf')
1886     call check_size(size(grid%vf),'vf')
1887     call check_size(size(grid%zsf),'zsf')
1888     call check_size(size(grid%dzdxf),'dzdxf')
1889     call check_size(size(grid%dzdyf),'dzdyf')
1890     call check_size(size(grid%bbb),'bbb')
1891     call check_size(size(grid%phiwc),'phiwc')
1892     call check_size(size(grid%r_0),'r_0')
1893     call check_size(size(grid%fgip),'fgip')
1894     call check_size(size(grid%ischap),'ischap')
1895     call check_size(size(grid%fmc_g),'fmc_g')
1896     call message('check_grid_alloc: checking if grid arrays used in fuel moisture are allocated',level=2)
1897     if (config_flags%fmoisti_run==1)then
1898         call check_size(size(grid%fmc_gc),'fmc_gc')
1899         call check_size(size(grid%fmep),'fmep')
1900         call check_size(size(grid%fmc_equi),'fmc_equi')
1901         call check_size(size(grid%fmc_lag),'fmc_lag')
1902         call check_size(size(grid%rain_old),'rain_old')
1903         call check_size(size(grid%t2_old),'t2_old')
1904         call check_size(size(grid%q2_old),'q2_old')
1905         call check_size(size(grid%rain_old),'rain_old')
1906         call check_size(size(grid%psfc_old),'psfc_old')
1907         call check_size(size(grid%rh_fire),'rh_fire')
1908         call check_size(size(grid%fmc_g),'fmc_g')
1909         call check_size(size(grid%fmc_gc_f),'fmc_gc_f')
1910     endif
1911     if (config_flags%fmoisti_interp .eq. 1)then
1912         call check_size(size(grid%fmc_gc),'fmc_gc')
1913     endif
1914     call message('check_grid_alloc: checking if grid arrays used in emissions are allocated',level=2)
1915     call check_size(size(grid%chem),'chem',config_flags%chem_opt)
1916     call check_size(size(grid%tracer),'tracer',config_flags%tracer_opt)
1917
1918     if(ibeh.eq.2)then
1919         call check_size(size(grid%t2),'t2')
1920         call check_size(size(grid%psfc),'psfc')
1921         call check_size(size(grid%q2),'q2')
1922         call check_size(size(grid%tempf),'tempf')
1923         call check_size(size(grid%rhof),'rhof')
1924         call check_size(size(grid%ffwidth),'ffwidth')
1925     endif
1926
1927 end subroutine check_grid_alloc
1928
1929 subroutine check_size(sz,name,req)
1930 integer, intent(in)::sz
1931 integer, optional::req
1932 character(len=*), intent(in)::name
1933 character(len=128)::msg
1934 logical::strict = .false.
1935
1936 write(msg,'(a,a10,a,i12)')'array ',name,' has size ',sz
1937 call message(msg,level=3)
1938 if(present(req)) strict=req > 0
1939 if (sz < 2 .and. strict) then
1940     write(msg,'(a,a,a,i2,a)')'array ',name,' is size ',sz,': not allocated, add to registry package'
1941     call crash(trim(msg))
1942 endif
1943 end subroutine check_size
1944
1945
1946 !subroutine print_id
1947 !character(len=128)::id,msg
1948 !include "sfire_id.inc"
1949 !msg=id
1950 !call message(msg,level=1)
1951 !end subroutine print_id
1952
1953 end module module_fr_sfire_driver