| blob | 136fcec04e94bef2849ebf11258c1b3260351e45 |
2 ! SFIRE - Spread fire model in WRF-Fire
3 !
4 !*** Jan Mandel August 2007 - March 2011
5 !*** email: Jan.Mandel@gmail.com
7 ! For support please subscribe to the wrf-fire mailing list at NCAR at
8 ! http://mailman.ucar.edu/mailman/listinfo/wrf-fire
9 ! or go to http://www.openwfm.org/wiki/WRF-Fire_user_support
11 ! Current drafts of the technical documentation and
12 ! user's guide can be found at
14 ! http://www.openwfm.org/wiki/WRF-Fire_documentation
15 ! http://www.openwfm.org/wiki/WRF-Fire_publications
17 ! This module is the only entry point from WRF-ARW to the wildland
18 ! fire model. The call to sfire_driver advances the fire model by
19 ! one timestep. The fire model inputs the wind and outputs
20 ! temperature and humidity tendencies. The fire model also inputs a
21 ! number of constant arrays (fuel data, topography). Additional
22 ! arguments are model state (for data assimilation) and constant arrays
23 ! the model gives to WRF for safekeeping because it is not allowed
24 ! to save anything.
26 ! This code as of mid-2011 is described in [1]. If you use this code,
27 ! please acknowledge our work by citing [1].
28 ! Thank you.
30 ! Acknowledgements
31 !
32 ! The fire physics code is adapted from the CAWFE code [2].
33 ! The coupling with WRF is adapted from a code by Ned Patton,
34 ! coupling a Fortran 90 port of the CAWFE fire module to WRF [3].
35 ! Support of refined fire grids in WRF was provided by John Michalakes.
36 ! Jonathan D. Beezley has set up and maintained the WRF build and
37 ! execution environment, provided software engineering infrastructure
38 ! including synchronization with the WRF repository, and was responsibe
39 ! for all aspects of WRF modification. UCD students Minjeong Kim and
40 ! Volodymyr Kondratenko have contributed to the implementation of the
41 ! fire propagation by the level set method.
43 ! Refefences
44 !
45 ! [1] Jan Mandel, Jonathan D. Beezley, and Adam K. Kochanski, "Coupled
46 ! atmosphere-wildland fire modeling with WRF 3.3 and SFIRE 2011,
47 ! Geoscientific Model Development (GMD) 4, 591-610, 2011.
48 ! doi:10.5194/gmd-4-591-2011
49 !
50 ! [2] T. L. Clark, J. Coen, and D. Latham, Description of a coupled
51 ! atmosphere-fire model, Intl. J. Wildland Fire, vol. 13, pp. 49-64,
52 ! 2004
53 !
54 ! [3] Edward G. Patton and Janice L. Coen, WRF-Fire: A Coupled
55 ! Atmosphere-Fire Module for WRF, Preprints of Joint MM5/Weather
56 ! Research and Forecasting Model Users' Workshop, Boulder, CO,
57 ! June 22-25, 2004, pp. 221-223, NCAR
58 !
59 ! ---------------------------------------------
60 !
61 ! CURRENT ACTIVITY
62 !
63 ! For current activity and development trends please check out
64 ! http://ccm.ucdenver.edu/wiki/User:Jmandel/blog
65 ! http://www.openwfm.org/wiki/WRF-Fire_development_notes
66 !
68 module module_fr_sfire_driver
69 ! use this module for standalone call, you only need to provide some mock-up wrf modules
71 use module_fr_sfire_model, only: sfire_model
72 use module_fr_sfire_phys, only: fire_params, init_fuel_cats, fuel_moisture, &
73 advance_moisture, moisture_classes, &
74 fire_rate_of_spread
75 use module_fr_sfire_atm, only: apply_windrf,interpolate_wind2fire_height,interpolate_atm2fire, &
76 interpolate_z2fire,setup_wind_log_interpolation
77 use module_fr_sfire_util
78 !use module_fr_sfire_util, only: lines_type,fire_max_lines
79 ! Driver layer modules
80 #ifdef DM_PARALLEL
81 USE module_dm , ONLY : ntasks_x,ntasks_y,local_communicator,mytask,ntasks,wrf_dm_sum_reals
82 USE module_comm_dm , ONLY : halo_fire_fuel_sub, halo_fire_tign_sub, halo_fire_wind_f_sub, &
83 halo_fire_wind_a_sub, halo_fire_ph_sub, halo_fire_zsf_sub, halo_fire_longlat_sub, &
84 halo_fire_phb_sub, halo_fire_z0_sub, halo_fire_lfn_sub, HALO_FIRE_MAG_sub, HALO_FIRE_MFG_sub, &
85 halo_fire_ndwi_sub
86 #endif
87 use module_fr_sfire_atm, only: read_emissions_table, add_fire_emissions
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
96 implicit none
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
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
107 logical:: fmoist_run, fmoist_interp, fire_run ! which kind of model to run overall
109 integer, parameter:: ifun_beg=1, ifun_step=3, ifun_end=6
111 contains
113 ! to write debugging information
114 #define DEBUG_OUT
116 subroutine sfire_driver_em ( grid , config_flags &
117 ,time_step_start,dt &
118 ,fire_ifun_start,fire_ifun_end,tsteps &
119 ,ids,ide, kds,kde, jds,jde &
120 ,ims,ime, kms,kme, jms,jme &
121 ,ips,ipe, kps,kpe, jps,jpe &
122 ,ifds,ifde, jfds,jfde &
123 ,ifms,ifme, jfms,jfme &
124 ,ifps,ifpe, jfps,jfpe &
125 ,rho,z_at_w,dz8w &
126 )
128 !*** purpose: driver from grid structure
132 implicit none
133 !*** arguments
134 TYPE(domain) , TARGET :: grid ! state
135 TYPE (grid_config_rec_type) , INTENT(IN) :: config_flags ! namelist
136 real, intent(in):: time_step_start, dt
137 integer, intent(in):: fire_ifun_start,fire_ifun_end,tsteps ! driver cycle controls
138 integer, intent(in):: &
139 ids,ide, kds,kde, jds,jde, &
140 ims,ime, kms,kme, jms,jme, &
141 ips,ipe, kps,kpe, jps,jpe, &
142 ifds,ifde, jfds,jfde, &
143 ifms,ifme, jfms,jfme, &
144 ifps,ifpe, jfps,jfpe
145 real,dimension(ims:ime, kms:kme, jms:jme),intent(in), optional::rho,z_at_w,dz8w
147 !*** local
148 TYPE(lines_type):: ignition, hfx
149 integer::fire_ifun,ir,jr,istep,itimestep,i,ipe1,kpe1,jpe1,j
150 logical::restart
151 real, dimension(ifms:ifme, jfms:jfme)::lfn_out
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
154 type(fire_params)::fp
155 real:: moisture_time
157 logical:: run_advance_moisture,run_fuel_moisture, moisture_initializing
158 real:: dt_moisture
162 !*** executable
164 call sfire_debug_hook(config_flags%fire_debug_hook_sec)
165 call time_start
166 if(fire_ifun_start.le.1)call print_id ! print id only once, during initialization
169 ! **** THE FOLLOWING REALLY SHOULD BE DONE ONCE NOT EVERY TIMESTEP
171 ! populate our structures from wrf
173 call set_fp_from_grid(grid,fp)
175 ! copy configuration flags to sfire internal structures
176 call set_flags(config_flags)
179 ! get ignition data
180 call fire_ignition_convert (config_flags,ignition, &
181 grid%fxlong, grid%fxlat, &
182 ifds,ifde, jfds,jfde, &
183 ifms,ifme, jfms,jfme, &
184 ifps,ifpe, jfps,jfpe )
185 call fire_hfx_convert (config_flags,hfx)
187 ! store computed mesh units
188 grid%unit_fxlong = ignition%unit_fxlong
189 grid%unit_fxlat = ignition%unit_fxlat
192 #ifndef SFIRE_STANDALONE
194 ! see what we got from wrf
195 !! need to replace ipe by min(ide-1,ipe) and similarly jpe
197 if(fire_print_msg.ge.2 .and. fire_ifun_start .gt. 1)then
199 ipe1=min(ide-1,ipe)
200 jpe1=min(jde-1,jpe)
201 kpe1=kpe-1
203 max_u=fun_real(REAL_AMAX, &
204 ims,ime,kms,kme,jms,jme, & ! memory dims
205 ids,ide,kds,kde,jds,jde, & ! domain dims
206 ips,ipe1,kps,kpe1,jps,jpe1, & ! patch or tile dims
207 1,0,0, & ! staggering
208 grid%u_2,grid%u_2)
210 max_v=fun_real(REAL_AMAX, &
211 ims,ime,kms,kme,jms,jme, & ! memory dims
212 ids,ide,kds,kde,jds,jde, & ! domain dims
213 ips,ipe1,kps,kpe1,jps,jpe1, & ! patch or tile dims
214 0,0,1, & ! staggering
215 grid%v_2,grid%v_2)
217 max_w=fun_real(REAL_AMAX, &
218 ims,ime,kms,kme,jms,jme, & ! memory dims
219 ids,ide,kds,kde,jds,jde, & ! domain dims
220 ips,ipe1,kps,kpe1,jps,jpe1, & ! patch or tile dims
221 0,1,0, & ! staggering
222 grid%w_2,grid%w_2)
224 !write(msg,93)time_step_start,'Maximal u v w wind',max_u,max_v,max_w,'m/s'
225 !call message(msg,0)
226 !write(msg,92)time_step_start,'Min and max rho ',min_rho,max_rho,'kg/m^3'
227 !call message(msg,0)
229 write(msg,91)time_step_start,'Maximal u wind ',max_u,'m/s'
230 call message(msg,0)
231 write(msg,91)time_step_start,'Maximal v wind ',max_v,'m/s'
232 call message(msg,0)
233 write(msg,91)time_step_start,'Maximal w wind ',max_w,'m/s'
234 call message(msg,0)
236 if (present(rho)) then
238 max_rho=fun_real(REAL_MAX, &
239 ims,ime,kms,kme,jms,jme, & ! memory dims
240 ids,ide,kds,kde,jds,jde, & ! domain dims
241 ips,ipe1,kps,kpe1,jps,jpe1, & ! patch or tile dims
242 0,0,0, & ! staggering
243 rho,rho)
245 min_rho=fun_real(REAL_MIN, &
246 ims,ime,kms,kme,jms,jme, & ! memory dims
247 ids,ide,kds,kde,jds,jde, & ! domain dims
248 ips,ipe1,kps,kpe1,jps,jpe1, & ! patch or tile dims
249 0,0,0, & ! staggering
250 rho,rho)
253 write(msg,91)time_step_start,'Minimal rho ',min_rho,'kg/m^3'
254 call message(msg,0)
255 write(msg,91)time_step_start,'Maximal rho ',max_rho,'kg/m^3'
256 call message(msg,0)
258 endif
261 93 format('Time ',f11.3,' s ',a,3e12.3,1x,a)
262 92 format('Time ',f11.3,' s ',a,2e12.3,1x,a)
263 91 format('Time ',f11.3,' s ',a,e12.3,1x,a)
266 endif
267 #endif
270 ! refinement r
271 ir=grid%sr_x ! refinement ratio
272 jr=grid%sr_y
273 write(msg,'(a,2i4)')'fire mesh refinement ratios', ir,jr
274 call message(msg)
275 if(ir.le.0.or.jr.le.0)then
276 call crash('fire mesh refinement ratio must be positive')
277 endif
278 itimestep=grid%itimestep
279 restart=config_flags%restart .or. config_flags%cycling .or. config_flags%fire_restart ! skip state initialization
283 ! **** moisture model
285 ! decide what to run - moisture, interpolation, or fire model itself
286 fmoist_run = config_flags%fmoist_run
287 fmoist_interp = config_flags%fmoist_interp
288 if(fire_fmc_read.ne.0.and.fmoist_run)call crash('fmoist_run=T requires fire_fmc_read=0')
289 fire_run = .not. config_flags%fmoist_only
291 !decide what to run
292 moisture_time = time_step_start
293 run_advance_moisture = .false. ! default
294 run_fuel_moisture = .false. ! default
295 moisture_initializing = fire_ifun_start < 3
299 if(fmoist_run)then
300 if(moisture_initializing)then
301 if(fire_ifun_end>2)call crash('initialization must be run separately')
302 grid%fmoist_lasttime=moisture_time ! initialize the last time the model has run to start of run
303 grid%fmoist_nexttime=moisture_time
304 call message('moisture initialization')
305 run_advance_moisture = .true.
306 else ! regular timestep
307 if(config_flags%fmoist_freq > 0)then ! regular timestep. go by multiples?
308 if(mod(grid%itimestep,config_flags%fmoist_freq) .eq. 0)then
309 write(msg,'(a,i10,a,i10)')'moisture model runs because timestep ',grid%itimestep,' is a multiple of ',config_flags%fmoist_freq
310 call message(msg)
311 run_advance_moisture = .true.
312 endif
313 else
314 if(.not. moisture_time < grid%fmoist_nexttime) then ! no, by time interval
315 write(msg,'(a,f12.2,a)')'moisture model runs because time ',grid%fmoist_nexttime,'s has arrived'
316 call message(msg)
317 run_advance_moisture = .true.
318 endif
319 endif
320 if(run_advance_moisture)then ! decide on timing
321 dt_moisture = moisture_time - grid%fmoist_lasttime ! Time since moisture model run the last time. Should be long.
322 grid%fmoist_lasttime = moisture_time
323 if(config_flags%fmoist_freq > 0)then
324 write(msg,'(a,f12.2,a,i10,a)')'moisture time step is ',dt_moisture,'s running every ',config_flags%fmoist_freq,' steps'
325 call message(msg)
326 else
327 grid%fmoist_nexttime = moisture_time + config_flags%fmoist_dt
328 write(msg,'(a,f12.2,a,f12.2,a)')'moisture time step is ',dt_moisture,'s next run at ',grid%fmoist_nexttime,'s'
329 call message(msg)
330 endif
331 if(fmoist_interp)then
332 call message('moisture interpolation to fuels will run because moisture model does')
333 run_fuel_moisture=.true.
334 endif
335 endif
336 endif
337 elseif(itimestep.eq.1.and.fmoist_interp)then
338 call message('initializing, moisture interpolation to fuels will run from input data')
339 run_fuel_moisture=.true.
340 endif
342 !$OMP CRITICAL(SFIRE_DRIVER_CRIT)
343 write(msg,'(a,i1,a,i1,a,l1)') &
344 'sfire_driver_em: ifun from ',fire_ifun_start,' to ',fire_ifun_end,' restart=',restart
345 !$OMP END CRITICAL(SFIRE_DRIVER_CRIT)
346 call message(msg)
348 do istep=0,tsteps ! istep >0 is for testing only, exit after the first call
349 itimestep = grid%itimestep + istep ! in the first call, do fire_test_steps steps of the fire model
351 do fire_ifun=fire_ifun_start,fire_ifun_end
353 ! 1 = moisture_initialize run pass 1: interpolate height to zsf=terrain
354 ! 2 = initialize run pass 2: set fuel data, terrain gradient
355 ! 3 = initialize timestep: interpolate winds, check for ignition, time step on moisture model
356 ! 4 = do one timestep
357 ! 5 = copy timestep output to input
358 ! 6 = compute output fluxes
360 #ifdef DM_PARALLEL
362 if(fire_run)then
364 if(fire_ifun.eq.1)then
366 ! halo exchange on topography
367 #include "HALO_FIRE_LONGLAT.inc"
368 !! if(fire_topo_from_atm.eq.1)then
369 !!#include "HALO_FIRE_HT.inc"
370 !! endif
371 ! base geopotential and roughness
372 #include "HALO_FIRE_PHB.inc"
373 #include "HALO_FIRE_Z0.inc"
374 if(kfmc_ndwi > 0 .and. fndwi_from_ndwi .eq.1)then
375 #include "HALO_FIRE_NDWI.inc"
376 endif
378 elseif(fire_ifun.eq.2)then
379 ! halo exchange on zsf width 2
380 #include "HALO_FIRE_ZSF.inc"
382 if(config_flags%chem_opt>0 .or. config_flags%tracer_opt > 0)then
383 ! need reading fuel categories first
384 call read_emissions_table(config_flags%chem_opt,config_flags%tracer_opt)
385 endif
388 elseif(fire_ifun.eq.3)then
389 ! halo exchange on atm winds and geopotential, width 1 for interpolation
390 #include "HALO_FIRE_WIND_A.inc"
391 #include "HALO_FIRE_PH.inc"
393 elseif(fire_ifun.eq.4)then
394 ! halo exchange on fire winds width 2 for a 2-step RK method
395 #include "HALO_FIRE_WIND_F.inc"
397 if(run_fuel_moisture)then
398 ! have interpolated to the fire grid
399 #include "HALO_FIRE_MFG.inc"
400 endif
402 elseif(fire_ifun.eq.6)then
403 ! computing fuel_left needs ignition time from neighbors
404 #include "HALO_FIRE_TIGN.inc"
405 call message('halo exchange on lfn width 2')
406 #include "HALO_FIRE_LFN.inc"
408 endif
409 endif
410 #endif
411 ! print *,'dt: ',dt,grid%dt,' diff ', dt-grid%dt
412 ! need domain by 1 smaller, in last row.col winds are not set properly
413 call sfire_driver_phys ( &
414 fire_ifun, &
415 ids,ide-1, kds,kde, jds,jde-1, &
416 ims,ime, kms,kme, jms,jme, &
417 ips,min(ipe,ide-1), kps,kpe, jps,min(jpe,jde-1), &
418 ifds,ifde-ir, jfds,jfde-jr, &
419 ifms,ifme, jfms,jfme, &
420 ifps,min(ifpe,ifde-ir), jfps,min(jfpe,jfde-jr), &
421 ir,jr, & ! atm/fire grid ratio
422 grid%num_tiles, & ! atm grid tiling
423 grid%i_start,min(grid%i_end,ide-1), &
424 grid%j_start,min(grid%j_end,jde-1), &
425 itimestep,restart,config_flags%fire_fuel_read,config_flags%fire_fuel_cat, & ! in scalars
426 time_step_start,dt,grid%dx,grid%dy, &
427 grid%u_frame,grid%v_frame, &
428 config_flags%fire_ext_grnd,config_flags%fire_ext_crwn,config_flags%fire_crwn_hgt, &
429 ignition,hfx, & ! lines
430 grid%u_2,grid%v_2, & ! atm arrays in
431 grid%ph_2,grid%phb, & ! geopotential
432 grid%z0, & ! roughness height
433 grid%ht, & ! terrain height
434 grid%xlong,grid%xlat, & ! coordinates of atm grid centers, for ignition location
435 grid%lfn,grid%tign_g,grid%fuel_frac, & ! state arrays, fire grid
436 grid%fire_area, & ! redundant, for display, fire grid
437 grid%fuel_frac_burnt, &
438 lfn_out, & ! work - one timestep
439 grid%avg_fuel_frac, & ! out redundant arrays, atm grid
440 grid%grnhfx,grid%grnqfx,grid%canhfx,grid%canqfx, & ! out redundant arrays, atm grid
441 grid%uah,grid%vah, &
442 grid%fgrnhfx,grid%fgrnqfx,grid%fcanhfx,grid%fcanqfx, & ! out redundant arrays, atm grid
443 grid%ros,grid%flineint,grid%flineint2, & ! diagnostic variables
444 grid%f_ros0,grid%f_rosx,grid%f_rosy,grid%f_ros,& ! fire risk spread
445 grid%f_int,grid%f_lineint,grid%f_lineint2, & ! fire risk intensities
446 grid%f_ros11,grid%f_ros12,grid%f_ros13,grid%f_ros21, & ! fire spread in nodal directions
447 grid%f_ros23,grid%f_ros31,grid%f_ros32,grid%f_ros33, & ! fire spread in nodal directions
448 grid%fxlong,grid%fxlat, &
449 grid%fire_hfx, & !
450 grid%nfuel_cat, & ! input, or internal for safekeeping
451 grid%fuel_time, &
452 grid%fz0, grid%fwh, &
453 fp, & ! structure with pointers passed to spread rate calculation
454 config_flags%nfmc, & ! moisture model variables start
455 run_advance_moisture,run_fuel_moisture,dt_moisture, & ! moisture model control
456 config_flags%fmep_decay_tlag, & ! moisture extended model assim. diffs decay time lag
457 grid%rainc, grid%rainnc, & ! accumulated rain from different sources
458 grid%t2, grid%q2, grid%psfc, & ! temperature (K), vapor contents (kg/kg), pressure (Pa) at the surface
459 grid%rain_old, & ! previous value of accumulated rain
460 grid%t2_old, grid%q2_old, grid%psfc_old, & ! previous values of the atmospheric state at surface
461 grid%rh_fire, & ! relative humidity, diagnostics
462 grid%fmc_gc, & ! fuel moisture fields updated, by class, assumed set to something reasonable
463 grid%fmep, & ! fuel moisture extended model parameters
464 grid%fmc_equi, & ! fuel moisture fields updated, by class, equilibrium diagnostic
465 grid%fmc_lag, & ! fuel moisture fields updated, by class, tendency diagnostic
466 fp%fmc_g, & ! write-only alias. need to exit before using fp again
467 grid%ndwi, &
468 grid%fndwi)
471 #ifdef DM_PARALLEL
472 if(fire_run)then
473 if(fire_ifun.eq.2)then
474 ! halo exchange on all fuel data width 2
475 #include "HALO_FIRE_FUEL.inc"
476 ! fire state was initialized
477 call message('halo exchange on lfn width 2')
478 #include "HALO_FIRE_LFN.inc"
479 endif
480 if(run_fuel_moisture)then
481 if(fire_ifun.eq.3)then
482 ! prepare for interpolation to the fire grid
483 #include "HALO_FIRE_MAG.inc"
484 endif
485 endif
486 endif
487 #endif
489 if(fire_ifun.eq.6)then
490 if(config_flags%chem_opt>0 .or. config_flags%tracer_opt>0)then
491 if(.not.(present(rho).and.present(dz8w)))then
492 call crash('sfire_driver_em: must have rho and dz8w to call add_fire_emissions')
493 endif
494 call add_fire_emissions( &
495 config_flags%chem_opt,config_flags%tracer_opt,dt,grid%dx,grid%dy, &
496 ifms,ifme,jfms,jfme, &
497 ifps,ifpe,jfps,jfpe, & ! use patch instead of tile
498 ids,ide,kds,kde,jds,jde, &
499 ims,ime,kms,kme,jms,jme, &
500 ips,ipe,kps,kpe,jps,jpe, &
501 rho,dz8w, & ! from atmosphere state
502 grid%fgip, grid%fuel_frac_burnt, grid%nfuel_cat, & ! from fire state
503 grid%chem,grid%tracer) ! update/output
504 endif
505 endif
508 enddo
509 enddo
511 if(tsteps>0)call crash('sfire_driver_em: test run of uncoupled fire model completed')
512 call time_end('sfire')
514 end subroutine sfire_driver_em
516 !
517 !*******************
518 !
520 subroutine sfire_driver_phys (ifun, &
521 ids,ide, kds,kde, jds,jde, & ! atm grid dimensions
522 ims,ime, kms,kme, jms,jme, &
523 ips,ipe, kps,kpe, jps,jpe, &
524 ifds, ifde, jfds, jfde, & ! fire grid dimensions
525 ifms, ifme, jfms, jfme, &
526 ifps, ifpe, jfps, jfpe, & ! fire patch in - will use smaller
527 ir,jr, & ! atm/fire grid ratio
528 num_tiles,i_start,i_end,j_start,j_end, & ! atm grid tiling
529 itimestep,restart,ifuelread,nfuel_cat0, & ! in scalars
530 time_step_start,dt,dx,dy, & ! in scalars
531 u_frame,v_frame, &
532 fire_ext_grnd,fire_ext_crwn,fire_crwn_hgt, &
533 ignition,hfx, & ! lines
534 u,v, & ! in arrays, atm grid
535 ph,phb, &
536 z0,zs, &
537 xlong,xlat, &
538 lfn,tign,fuel_frac, & ! state arrays, fire grid
539 fire_area, & ! redundant state, fire grid
540 fuel_frac_burnt, &
541 lfn_out, & ! out level set function
542 avg_fuel_frac, &
543 grnhfx,grnqfx,canhfx,canqfx, & ! out redundant arrays, atm grid
544 uah,vah, & ! out atm grid
545 fgrnhfx,fgrnqfx,fcanhfx,fcanqfx, & ! out redundant arrays, fire grid
546 ros,flineint,flineint2, & ! diagnostic variables
547 f_ros0,f_rosx,f_rosy,f_ros, & ! fire risk spread
548 f_int,f_lineint,f_lineint2, & ! fire risk intensities
549 f_ros11,f_ros12,f_ros13,f_ros21, & ! fire spread in nodal directions
550 f_ros23,f_ros31,f_ros32,f_ros33, & ! fire spread in nodal directions
551 fxlong,fxlat, & !
552 fire_hfx, & !
553 nfuel_cat, & ! in array, data, fire grid, or constant internal
554 fuel_time, & ! save constant internal data, fire grid
555 fz0,fwh, &
556 fp, & ! fire params
557 nfmc, & ! number of fuel moisture classes
558 run_advance_moisture,run_fuel_moisture,dt_moisture,& ! moisture model control
559 fmep_decay_tlag, & ! moist. extended model assim. diffs time lag
560 rainc,rainnc, & ! accumulated rain from different sources
561 t2, q2, psfc, & ! temperature (K), vapor contents (kg/kg), pressure (Pa) at the surface
562 rain_old, & ! previous value of accumulated rain
563 t2_old, q2_old, psfc_old, & ! previous values of the atmospheric state at surface
564 rh_fire, & ! relative humidity, diagnostics
565 fmc_gc, & ! fuel moisture fields updated, by class, assumed set to something reasonable
566 fmep, & ! fuel moisture extended model parameters
567 fmc_equi, & ! fuel moisture fields updated, by class equilibrium diagnostic
568 fmc_lag, & ! fuel moisture fields updated, by class tendency diagnostic
569 fmc_g, & ! fuel moisture, alias of fp%fmc_g
570 ndwi, &
571 fndwi) ! ndwi on fire grid
574 implicit none
576 !*** arguments
578 integer, intent(in)::ifun, &
579 ids,ide, kds,kde, jds,jde, & ! atm domain bounds
580 ims,ime, kms,kme, jms,jme, & ! atm memory bounds
581 ips,ipe, kps,kpe, jps,jpe, & ! atm patch bounds
582 ifds, ifde, jfds, jfde, & ! fire domain bounds
583 ifms, ifme, jfms, jfme, & ! fire memory bounds
584 ifps, ifpe, jfps, jfpe, & ! fire patch bounds
585 ir,jr, & ! atm/fire grid refinement ratio
586 nfmc, & ! number of fuel moisture classes
587 itimestep, & ! number of this timestep
588 ifuelread, & ! how to initialize nfuel_cat:
589 ! -1=not at all, done outside
590 ! 0=from nfuel_cat0
591 ! 1=from altitude
592 ! 2=from file
593 nfuel_cat0, & ! fuel category to initialize everything to
594 num_tiles ! number of tiles
596 logical, intent(in)::restart
598 integer,dimension(num_tiles),intent(in) :: i_start,i_end,j_start,j_end ! atm grid tiling
600 real, intent(in):: &
601 time_step_start, & ! time step start
602 dt, & ! time step length
603 dx,dy, & ! atm grid step
604 u_frame,v_frame, & ! velocity offset
605 fire_crwn_hgt, & ! lowest height crown fire heat is released (m)
606 fire_ext_grnd, & ! extinction depth of ground fire heat (m)
607 fire_ext_crwn ! and for the canopy (m)
610 TYPE (lines_type), intent(inout):: ignition,hfx
612 real,intent(in),dimension(ims:ime,kms:kme,jms:jme)::u,v, & ! wind velocity (m/s) (staggered atm grid)
613 ph, phb ! geopotential (w-points atm grid)
614 real,intent(in),dimension(ims:ime, jms:jme):: z0, & ! roughness height
615 zs ! terrain height
616 real,intent(out),dimension(ims:ime,jms:jme)::&
617 uah, & ! atm wind at fire_wind_height, diagnostics
618 vah ! atm wind at fire_wind_height, diagnostics
620 real, dimension(ims:ime, jms:jme), intent(inout)::xlong, xlat, ndwi ! inout because of extension at bdry
622 real, intent(inout), dimension(ifms:ifme,jfms:jfme):: & ! fuel data; can be also set inside (cell based, fire grid)
623 fz0,fwh, &
624 nfuel_cat,fndwi
626 real, intent(inout), dimension(ifms:ifme, jfms:jfme):: &
627 lfn,tign,fuel_frac, & ! state: level function, ign time, fuel left
628 lfn_out ! fire wind velocities
630 real, intent(out), dimension(ifms:ifme, jfms:jfme):: &
631 fire_area, & ! fraction of each cell burning
632 fuel_frac_burnt
634 real, intent(out), dimension(ims:ime, jms:jme):: & ! redundant arrays, for display purposes only (atm grid)
635 avg_fuel_frac, & ! average fuel fraction
636 grnhfx, & ! heat flux from ground fire (W/m^2)
637 grnqfx, & ! moisture flux from ground fire (W/m^2)
638 canhfx, & ! heat flux from crown fire (W/m^2)
639 canqfx ! moisture flux from crown fire (W/m^2)
641 real, intent(out), dimension(ifms:ifme, jfms:jfme):: & ! redundant arrays, for display only, fire grid
642 fgrnhfx, & ! heat flux from ground fire (W/m^2)
643 fgrnqfx, & ! moisture flux from ground fire (W/m^2)
644 fcanhfx, & ! heat flux from crown fire (W/m^2)
645 fcanqfx, & ! moisture flux from crown fire (W/m^2)
646 ros,flineint,flineint2, & ! diagnostic variables
647 f_ros0,f_rosx,f_rosy,f_ros, & ! fire risk spread
648 f_int,f_lineint,f_lineint2, & ! fire risk intensities
649 f_ros11,f_ros12,f_ros13,f_ros21, & ! fire spread in nodal directions
650 f_ros23,f_ros31,f_ros32,f_ros33 ! fire spread in nodal directions
653 ! moisture model arguments
654 logical, intent(in)::run_advance_moisture,run_fuel_moisture
655 real, intent(in)::dt_moisture
656 real, intent(in)::fmep_decay_tlag
657 real, intent(in), dimension(ims:ime,jms:jme):: t2, q2, psfc, rainc, rainnc
658 real, intent(inout), dimension(ims:ime,jms:jme):: t2_old, q2_old, psfc_old, rain_old
659 real, intent(out),dimension(ims:ime,jms:jme):: rh_fire
660 real, intent(inout), dimension(ims:ime,nfmc,jms:jme):: fmc_gc
661 real, intent(inout), dimension(ims:ime,2,jms:jme):: fmep
662 real, intent(out), dimension(ims:ime,nfmc,jms:jme):: fmc_equi,fmc_lag
663 real, intent(inout), dimension(ifms:ifme,jfms:jfme):: fmc_g
667 ! ***** data (constant in time) *****
669 real, dimension(ifms:ifme, jfms:jfme), intent(inout)::fxlong,fxlat, & ! fire mesh coordinates
670 fire_hfx
671 real, intent(out), dimension(ifms:ifme, jfms:jfme)::fuel_time ! fire params arrays
673 type(fire_params),intent(inout)::fp
675 !*** local
676 real :: dxf,dyf,time_start,latm, s
677 integer :: its,ite,jts,jte,kts,kte, & ! tile
678 ij,i,j,k,id,pid,ipe1,jpe1,ite1,jte1, &
679 ii,jj, &
680 ifts,ifte,jfts,jfte ! fire tile
681 character(len=128)::msg
682 character(len=3)::kk
683 real, parameter::zero=0.
685 !*** executable
687 ! time - assume dt does not change
688 ! time_start = (itimestep-1) * dt ! timestep 1 starts at 0
689 ! print *,'time_start: ',time_start,time_step_start,' diff ', time_start-time_step_start
690 time_start = time_step_start ! use the time passed from wrf
692 ! fire mesh step
693 dxf=dx/ir
694 dyf=dy/jr
697 !$OMP CRITICAL(SFIRE_DRIVER_CRIT)
698 write(msg,'(a,i5)')'sfire_driver_phys stage ',ifun
699 call message(msg)
700 write(msg,'(a,2f15.6)')'atmosphere mesh step:',dx,dy
701 call message(msg)
702 write(msg,'(a,2f15.6)')'fire mesh step: ',dxf,dyf
703 call message(msg)
704 write(msg,7001)'atm domain ','ids',ids,ide,jds,jde
705 call message(msg)
706 write(msg,7001)'atm memory ','ims',ims,ime,jms,jme
707 call message(msg)
708 write(msg,7001)'atm patch ','ips',ips,ipe,jps,jpe
709 call message(msg)
710 write(msg,7001)'fire domain ','ifds',ifds,ifde,jfds,jfde
711 call message(msg)
712 write(msg,7001)'fire memory ','ifms',ifms,ifme,jfms,jfme
713 call message(msg)
714 write(msg,7001)'fire patch ','ifps',ifps,ifpe,jfps,jfpe
715 call message(msg)
716 !$OMP END CRITICAL(SFIRE_DRIVER_CRIT)
718 ! check mesh dimensions
719 call check_fmesh(ids,ide,ifds,ifde,ir,'id') ! check if atm and fire grids line up
720 call check_fmesh(jds,jde,jfds,jfde,jr,'jd')
721 call check_fmesh(ips,ipe,ifps,ifpe,ir,'ip')
722 call check_fmesh(jps,jpe,jfps,jfpe,jr,'jp')
723 call check_mesh_2dim(ips,ipe,jps,jpe,ims,ime,jms,jme) ! check if atm patch fits in atm array
724 call check_mesh_2dim(ifps,ifpe,jfps,jfpe,ifms,ifme,jfms,jfme) ! check if fire patch fits in fire array
725 call check_mesh_2dim(ips,ipe,jps,jpe,ids,ide,jds,jde) ! check if atm patch fits in atm domain
726 call check_mesh_2dim(ifps,ifpe,jfps,jfpe,ifds,ifde,jfds,jfde) ! check if fire patch fits in fire domain
728 pid=0
729 if(fire_print_file.gt.0)then
730 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
731 endif
734 if(ifun.eq.1)then
735 call init_fuel_cats(fmoist_run .or. fmoist_interp) ! properties of fuel categories and moisture classes from namelist.fire
736 endif
738 if(ifun.eq.3)then
739 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')
740 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')
741 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')
742 endif
744 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')
745 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')
748 ! fake atm tile bounds
749 kts=kps
750 kte=kpe
752 ! staggered atm patch bounds
753 ipe1=ifval(ipe.eq.ide,ipe+1,ipe)
754 jpe1=ifval(jpe.eq.jde,jpe+1,jpe)
756 ! set up fire tiles & interpolate to fire grid
757 !$OMP PARALLEL DO PRIVATE(ij,its,ite,jts,jte,ite1,jte1,ifts,ifte,jfts,jfte,msg,id) &
758 !$OMP SCHEDULE(STATIC)
759 do ij=1,num_tiles
761 id = ifval(pid.ne.0,pid+ij*10000,0) ! for print
763 ! set up tile bounds
764 its = i_start(ij) ! start atmospheric tile in i
765 ite = i_end(ij) ! end atmospheric tile in i
766 jts = j_start(ij) ! start atmospheric tile in j
767 jte = j_end(ij) ! end atmospheric tile in j
768 ifts= (its-ids)*ir+ifds ! start fire tile in i
769 ifte= (ite-ids+1)*ir+ifds-1 ! end fire tile in i
770 jfts= (jts-jds)*jr+jfds ! start fire tile in j
771 jfte= (jte-jds+1)*jr+jfds-1 ! end fire tile in j
773 ! staggered atm tile bounds
774 ite1=ifval(ite.eq.ide,ite+1,ite)
775 jte1=ifval(jte.eq.jde,jte+1,jte)
777 !$OMP CRITICAL(SFIRE_DRIVER_CRIT)
778 write(msg,'(a,i3,1x,a,i7,1x,a,i3)')'tile=',ij,' id=',id,' ifun=',ifun
779 call message(msg)
780 write(msg,7001)'atm tile ','its',its,ite,jts,jte
781 call message(msg)
782 write(msg,7001)'fire tile ','ifts',ifts,ifte,jfts,jfte
783 call message(msg)
784 !$OMP END CRITICAL(SFIRE_DRIVER_CRIT)
786 ! check the tiles
787 call check_mesh_2dim(its,ite,jts,jte,ips,ipe,jps,jpe) ! check if atm tile fits in atm patch
788 call check_mesh_2dim(ifts,ifte,jfts,jfte,ifps,ifpe,jfps,jfpe) ! check if fire tile fits in fire patch
789 call check_mesh_2dim(ifts-2,ifte+2,jfts-2,jfte+2,ifms,ifme,jfms,jfme)! check if fire node tile fits in memory
792 !$OMP CRITICAL(SFIRE_DRIVER_CRIT)
793 write(msg,'(a,i6,a,2(f15.6,a))')'time step',itimestep,' at',time_start,' duration',dt,'s'
794 call message(msg)
795 7001 format(a,' dimensions ',a4,':',i6,' to ',i6,' by ',i6,' to ',i6)
796 write(msg,'(a,2i9)')'refinement ratio:',ir,jr
797 !$OMP END CRITICAL(SFIRE_DRIVER_CRIT)
799 if(run_advance_moisture)then
800 if(ifun.eq.3)then
802 ! one timestep of the moisture model
803 call message('advance_moisture start')
804 call advance_moisture( &
805 itimestep.eq.1, & ! initialize?
806 ims,ime, jms,jme, & ! memory dimensions
807 its,ite, jts,jte, & ! tile dimensions
808 nfmc, & ! number of moisture fields
809 dt_moisture, & ! moisture model time step
810 fmep_decay_tlag, & ! moisture extended model assim. diffs decay tlag
811 rainc, rainnc, & ! accumulated rain
812 t2, q2, psfc, & ! temperature (K), vapor contents (kg/kg), pressure (Pa) at the surface
813 rain_old, & ! previous value of accumulated rain
814 t2_old, q2_old, psfc_old, & ! previous values of the atmospheric state at surface
815 rh_fire, & ! relative humidity, diagnostics
816 fmc_gc, & ! fuel moisture fields updated, by class, assumed set to something reasonable
817 fmep, & ! fuel moisture extended model parameters
818 fmc_equi, & ! fuel moisture fields updated, by class equilibrium diagnostic
819 fmc_lag & ! fuel moisture fields updated, by class tendency diagnostic
820 )
821 call message('advance_moisture end')
822 endif
823 endif
825 if(fire_run)then
827 if(ifun.eq.2)then ! interpolate
829 if(restart)then
831 call message('restart - interpolation skipped')
833 else
834 if(kfmc_ndwi > 0 .and. fndwi_from_ndwi .eq.1)then
835 call print_2d_stats(ips,ipe,jps,jpe,ims,ime,jms,jme,fndwi,'driver:ndwi')
836 call interpolate_z2fire(id,0, & ! for debug output, <= 0 no output, extend strip
837 ids,ide, jds,jde, & ! atm grid dimensions
838 ims,ime, jms,jme, &
839 ips,ipe,jps,jpe, &
840 its,ite,jts,jte, &
841 ifds, ifde, jfds, jfde, & ! fire grid dimensions
842 ifms, ifme, jfms, jfme, &
843 ifts,ifte,jfts,jfte, &
844 ir,jr, & ! atm/fire grid ratio
845 ndwi, & ! atm grid arrays in
846 fndwi) ! fire grid arrays out
847 call print_2d_stats(ifps,ifpe,jfps,jfpe,ifms,ifme,jfms,jfme,fndwi,'driver:fndwi')
848 endif
849 !
850 ! call print_2d_stats(ips,ipe,jps,jpe,ims,ime,jms,jme,zs,'driver:zs')
851 !
852 ! ! interpolate terrain height
853 ! if(fire_topo_from_atm.eq.1)then
854 ! call interpolate_z2fire(id,1, & ! for debug output, <= 0 no output
855 ! ids,ide, jds,jde, & ! atm grid dimensions
856 ! ims,ime, jms,jme, &
857 ! ips,ipe,jps,jpe, &
858 ! its,ite,jts,jte, &
859 ! ifds, ifde, jfds, jfde, & ! fire grid dimensions
860 ! ifms, ifme, jfms, jfme, &
861 ! ifts,ifte,jfts,jfte, &
862 ! ir,jr, & ! atm/fire grid ratio
863 ! zs, & ! atm grid arrays in
864 ! fp%zsf) ! fire grid arrays out
865 ! else
866 !!$OMP CRITICAL(SFIRE_DRIVER_CRIT)
867 ! write(msg,'(a,i3,a)')'fire_topo_from_atm=',fire_topo_from_atm,' assuming ZSF set, interpolation skipped'
868 !!$OMP END CRITICAL(SFIRE_DRIVER_CRIT)
869 ! endif
871 if(ignition%longlat .eq.0)then
872 ! set ideal fire mesh coordinates - used for ignition only
873 ! do not forget to set unit_fxlong, unit_fxlat outside of parallel loop
874 !call set_ideal_coord( dxf,dyf, &
875 ! ifds,ifde,jfds,jfde, &
876 ! ifms,ifme,jfms,jfme, &
877 ! ifts,ifte,jfts,jfte, &
878 ! fxlong,fxlat )
879 !call set_ideal_coord( dx,dy, &
880 ! ids,ide,jds,jde, &
881 ! ims,ime,jms,jme, &
882 ! its,ite,jts,jte, &
883 ! xlong,xlat )
884 elseif(use_atm_vars)then
885 ! assume halo xlong xlat
886 ! interpolate nodal coordinates
888 #ifdef DEBUG_OUT
889 call write_array_m(its,ite,jts,jte,ims,ime,jms,jme,xlat,'xlat',id)
890 call write_array_m(its,ite,jts,jte,ims,ime,jms,jme,xlong,'xlong',id)
891 #endif
893 if (interpolate_long_lat)then
894 call message('Intepolating node longitude and latitude to fire mesh')
895 call interpolate_z2fire(id,1, & ! for debug output, <= 0 no output
896 ids,ide, jds,jde, & ! atm grid dimensions
897 ims,ime, jms,jme, &
898 ips,ipe,jps,jpe, &
899 its,ite,jts,jte, &
900 ifds, ifde, jfds, jfde, & ! fire grid dimensions
901 ifms, ifme, jfms, jfme, &
902 ifts,ifte,jfts,jfte, &
903 ir,jr, & ! atm/fire grid ratio
904 xlat, & ! atm grid arrays in
905 fxlat) ! fire grid arrays out
907 call interpolate_z2fire(id,1, & ! for debug output, <= 0 no output
908 ids,ide, jds,jde, & ! atm grid dimensions
909 ims,ime, jms,jme, &
910 ips,ipe,jps,jpe, &
911 its,ite,jts,jte, &
912 ifds, ifde, jfds, jfde, & ! fire grid dimensions
913 ifms, ifme, jfms, jfme, &
914 ifts,ifte,jfts,jfte, &
915 ir,jr, & ! atm/fire grid ratio
916 xlong, & ! atm grid arrays in
917 fxlong) ! fire grid arrays out
918 endif
920 ! after the loop where zsf created exited and all synced
921 call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,fp%zsf,'driver_phys:zsf')
923 ! cannot initialize moisture model because T2 Q2 PSFC are not set yet
924 endif
925 endif
927 elseif(ifun.eq.3)then ! interpolate winds to the fire grid
929 if(use_atm_vars)then
931 call write_array_m(its,ite,jts,jte,ims,ime,jms,jme,z0,'z0',id)
932 call write_array_m3(its,ite1,kts,kde-1,jts,jte,ims,ime,kms,kme,jms,jme,u,'u_2',id)
933 call write_array_m3(its,ite,kts,kde-1,jts,jte1,ims,ime,kms,kme,jms,jme,v,'v_2',id)
934 call write_array_m3(its,ite,kts,kde,jts,jte,ims,ime,kms,kme,jms,jme,ph,'ph_2',id)
935 call write_array_m3(its,ite,kts,kde,jts,jte,ims,ime,kms,kme,jms,jme,phb,'phb',id)
937 if(fire_wind_log_interp.eq.4)then
938 call print_2d_stats(its,ite,jts,jte,ims,ime,jms,jme,z0,'driver_phys:z0')
939 call interpolate_atm2fire(id, & ! flag for debug output
940 ids,ide, kds,kde, jds,jde, & ! atm grid dimensions
941 ims,ime, kms,kme, jms,jme, &
942 ips,ipe, jps,jpe, &
943 its,ite,jts,jte, &
944 ifds, ifde, jfds, jfde, & ! fire grid dimensions
945 ifms, ifme, jfms, jfme, &
946 ifps, ifpe, jfps, jfpe, & ! fire patch bounds
947 ifts, ifte, jfts, jfte, &
948 ir,jr, & ! atm/fire grid ratio
949 u_frame, v_frame, & ! velocity frame correction
950 u,v, & ! 3D atm grid arrays in
951 ph,phb, &
952 z0,zs, & ! 2D atm grid arrays in
953 uah,vah, & ! 2D atm grid out
954 fp%vx,fp%vy) ! fire grid arrays out
956 call apply_windrf( &
957 ifms,ifme,jfms,jfme, &
958 ifts,ifte,jfts,jfte, &
959 nfuel_cat,fp%vx,fp%vy)
961 else
962 call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,fz0,'driver_phys:fz0')
963 call interpolate_wind2fire_height(id, & ! to identify debugging prints and files if needed
964 ids,ide, kds,kde, jds,jde, & ! atm grid dimensions
965 ims,ime, kms,kme, jms,jme, &
966 ips,ipe,jps,jpe, &
967 its,ite,jts,jte, &
968 ifds, ifde, jfds, jfde, & ! fire grid dimensions
969 ifms, ifme, jfms, jfme, &
970 ifps, ifpe, jfps, jfpe, & ! fire patch bounds
971 ifts,ifte,jfts,jfte, &
972 ir,jr, & ! atm/fire grid ratio
973 u_frame, v_frame, & ! velocity frame correction
974 u,v,ph,phb, & ! input atmospheric arrays
975 fz0,fwh, & ! input fire arrays
976 fp%vx,fp%vy) ! output fire arrays
978 if(fire_use_windrf.eq.1)then
979 call apply_windrf( &
980 ifms,ifme,jfms,jfme, &
981 ifts,ifte,jfts,jfte, &
982 nfuel_cat,fp%vx,fp%vy)
983 endif
985 endif
988 endif
990 elseif(ifun.eq.4)then
992 ! interpolate and compute weighted average to get the fuel moisture
993 !! print *,'ifun=4, run_fuel_moisture=',run_fuel_moisture
994 if(run_fuel_moisture)then
995 call message('fuel_moisture start')
996 call fuel_moisture( &
997 id, & ! for prints and maybe file names
998 nfmc, &
999 ids,ide, jds,jde, & ! atm grid dimensions
1000 ims,ime, jms,jme, &
1001 ips,ipe,jps,jpe, &
1002 its,ite,jts,jte, &
1003 ifds, ifde, jfds, jfde, & ! fire grid dimensions
1004 ifms, ifme, jfms, jfme, &
1005 ifts,ifte,jfts,jfte, &
1006 ir,jr, & ! atm/fire grid ratio
1007 nfuel_cat, & ! fuel data
1008 fndwi, & ! satellite sensing interpolated on fire grid
1009 fmc_gc, & ! moisture contents by class on atmospheric grid
1010 fmc_g & ! weighted fuel moisture contents on fire grid
1011 )
1012 call message('fuel_moisture end')
1013 endif
1015 call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,fmc_g,'driver_phys:fmc_g')
1016 endif
1020 call sfire_model (id,ifun,restart, &
1021 run_fuel_moisture, & ! if fuel moisture needs to be updated
1022 ifuelread,nfuel_cat0, & ! initialize fuel categories
1023 ifds,ifde,jfds,jfde, & ! fire domain dims
1024 ifms,ifme,jfms,jfme, & ! fire memory dims
1025 ifps,ifpe,jfps,jfpe, &
1026 ifts,ifte,jfts,jfte, & ! fire patch dims
1027 time_start,dt, & ! time and increment
1028 dxf,dyf, & ! fire mesh spacing
1029 ignition,hfx, & ! description of ignition lines
1030 fxlong,fxlat, & ! fire mesh coordinates
1031 fire_hfx, & ! given heat flux
1032 lfn,lfn_out,tign,fuel_frac, & ! state: level function, ign time, fuel left
1033 fire_area, & ! output: fraction of cell burning
1034 fuel_frac_burnt, & ! output: fuel fraction burned in this step
1035 fgrnhfx,fgrnqfx, & ! output: heat fluxes
1036 ros,flineint,flineint2, & ! diagnostic variables
1037 f_ros0,f_rosx,f_rosy,f_ros, & ! fire risk spread
1038 f_int,f_lineint,f_lineint2, & ! fire risk intensities
1039 nfuel_cat, & ! fuel data per point
1040 fuel_time,fwh,fz0, & ! save derived internal data
1041 fp & ! fire coefficients
1042 )
1044 if(ifun.eq.2)then
1045 call setup_wind_log_interpolation( &
1046 ids,ide, jds,jde, & ! atm grid dimensions
1047 ims,ime, jms,jme, &
1048 ips,ipe,jps,jpe, &
1049 its,ite,jts,jte, &
1050 ifds, ifde, jfds, jfde, & ! fire grid dimensions
1051 ifms, ifme, jfms, jfme, &
1052 ifts,ifte,jfts,jfte, &
1053 ir,jr, & ! atm/fire grid ratio
1054 z0, & ! atm grid arrays in
1055 nfuel_cat, & ! fuel data per point
1056 fz0,fwh) ! fire arrays out
1059 elseif(ifun.eq.6)then
1061 ! populate the rate of spread in the 8 directions
1062 do j=jfts,jfte
1063 do i=ifts,ifte
1064 f_ros11(i,j)=fire_rate_of_spread( dxf*(1-2), dyf*(1-2), i,j,fp)
1065 f_ros12(i,j)=fire_rate_of_spread( dxf*(1-2), dyf*(2-2), i,j,fp)
1066 f_ros13(i,j)=fire_rate_of_spread( dxf*(1-2), dyf*(3-2), i,j,fp)
1067 f_ros21(i,j)=fire_rate_of_spread( dxf*(2-2), dyf*(1-2), i,j,fp)
1068 f_ros23(i,j)=fire_rate_of_spread( dxf*(2-2), dyf*(3-2), i,j,fp)
1069 f_ros31(i,j)=fire_rate_of_spread( dxf*(3-2), dyf*(1-2), i,j,fp)
1070 f_ros32(i,j)=fire_rate_of_spread( dxf*(3-2), dyf*(2-2), i,j,fp)
1071 f_ros33(i,j)=fire_rate_of_spread( dxf*(3-2), dyf*(3-2), i,j,fp)
1072 enddo
1073 enddo
1074 call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,f_ros11,'driver_phys:f_ros11')
1075 call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,f_ros12,'driver_phys:f_ros12')
1076 call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,f_ros13,'driver_phys:f_ros13')
1077 call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,f_ros21,'driver_phys:f_ros21')
1078 call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,f_ros23,'driver_phys:f_ros23')
1079 call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,f_ros31,'driver_phys:f_ros31')
1080 call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,f_ros32,'driver_phys:f_ros32')
1081 call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,f_ros33,'driver_phys:f_ros33')
1083 ! heat fluxes into the atmosphere
1085 call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,fgrnhfx,'fire_driver:fgrnhfx')
1086 call print_2d_stats(ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,fgrnqfx,'fire_driver:fgrnqfx')
1088 ! sum the fluxes over atm cells
1089 if(use_atm_vars)then
1090 call sum_2d_cells( &
1091 ifms,ifme,jfms,jfme, &
1092 ifts,ifte,jfts,jfte, &
1093 fuel_frac, &
1094 ims, ime, jms, jme, &
1095 its,ite,jts,jte, &
1096 avg_fuel_frac)
1097 call sum_2d_cells( &
1098 ifms,ifme,jfms,jfme, &
1099 ifts,ifte,jfts,jfte, &
1100 fgrnhfx, &
1101 ims, ime, jms, jme, &
1102 its,ite,jts,jte, &
1103 grnhfx)
1104 !comment out the next call to get results as before commit 55fd92051196b796891b60cb7ec1c4bdb8800078
1105 call sum_2d_cells( &
1106 ifms,ifme,jfms,jfme, &
1107 ifts,ifte,jfts,jfte, &
1108 fgrnqfx, &
1109 ims, ime, jms, jme, &
1110 its,ite,jts,jte, &
1111 grnqfx)
1113 !$OMP CRITICAL(SFIRE_DRIVER_CRIT)
1114 write(msg,'(a,f6.3)')'fire-atmosphere feedback scaling ',fire_atm_feedback
1115 !$OMP end CRITICAL(SFIRE_DRIVER_CRIT)
1116 call message(msg)
1117 s = 1./(ir*jr)
1118 do j=jts,jte
1119 do i=its,ite
1120 ! scale ground fluxes to get the averages
1121 avg_fuel_frac(i,j)=avg_fuel_frac(i,j)*s
1122 grnhfx(i,j)=fire_atm_feedback*grnhfx(i,j)*s
1123 grnqfx(i,j)=fire_atm_feedback*grnqfx(i,j)*s
1124 ! we do not have canopy fluxes yet...
1125 canhfx(i,j)=0
1126 canqfx(i,j)=0
1127 enddo
1128 enddo
1130 call print_2d_stats(its,ite,jts,jte,ims,ime,jms,jme,grnhfx,'fire_driver:grnhfx')
1131 call print_2d_stats(its,ite,jts,jte,ims,ime,jms,jme,grnqfx,'fire_driver:grnqfx')
1132 endif
1134 endif ! ifun=6
1135 endif
1137 enddo ! tiles
1138 !$OMP END PARALLEL DO
1140 #ifdef DEBUG_OUT
1141 if(ifun.eq.1)then
1142 if(pid.ne.0)then
1143 call write_array_m(ips,ipe,jps,jpe,ims,ime,jms,jme,zs,'zs',pid)
1144 call write_array_m(ifps,ifpe,jfps,jfpe,ifms,ifme,jfms,jfme,fp%zsf,'zsf',pid)
1145 endif
1146 endif
1147 #endif
1149 if (ifun.eq.3)then
1150 call print_3d_stats_by_slice(ips,ipe,1,moisture_classes,jps,jpe,ims,ime,1,nfmc,jms,jme,fmc_gc,'fmc_gc')
1151 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')
1152 endif
1154 if (ifun.eq.4)then
1156 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')
1157 !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')
1158 !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')
1159 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')
1160 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')
1161 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')
1162 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')
1164 #ifdef DEBUG_OUT
1165 if(pid.gt.0)then
1166 ! call write_array_m(ips,ipe1,jps,jpe,ims,ime,jms,jme,uah,'uah',pid)
1167 ! call write_array_m(ips,ipe,jps,jpe1,ims,ime,jms,jme,vah,'vah',pid)
1168 call write_array_m(ips,ipe,jps,jpe,ims,ime,jms,jme,grnhfx,'grnhfx',pid)
1169 call write_array_m(ips,ipe,jps,jpe,ims,ime,jms,jme,grnqfx,'grnqfx',pid)
1170 call write_array_m3(ips,ipe1,kds,kde+1,jps,jpe,ims,ime,kms,kme,jms,jme,u,'u',pid)
1171 call write_array_m3(ips,ipe,kds,kde+1,jps,jpe1,ims,ime,kms,kme,jms,jme,v,'v',pid)
1172 call write_array_m(ifps,ifpe,jfps,jfpe,ifms,ifme,jfms,jfme,fp%vx,'uf',pid)
1173 call write_array_m(ifps,ifpe,jfps,jfpe,ifms,ifme,jfms,jfme,fp%vy,'vf',pid)
1174 endif
1175 #endif
1176 endif
1178 if(ifun.eq.5)then
1179 #ifdef DEBUG_OUT
1180 if(pid.gt.0)then
1181 call write_array_m(ifps,ifpe,jfps,jfpe,ifms,ifme,jfms,jfme,lfn,'lfn',pid)
1182 call write_array_m(ifps,ifpe,jfps,jfpe,ifms,ifme,jfms,jfme,tign,'tign',pid)
1183 endif
1184 #endif
1185 endif
1187 if(ifun.eq.6)then
1188 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')
1189 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')
1190 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')
1191 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')
1192 #ifdef DEBUG_OUT
1193 if(pid.gt.0)then
1194 call write_array_m(ips,ipe,jps,jpe,ims,ime,jms,jme,grnhfx,'grnhfx',pid)
1195 call write_array_m(ips,ipe,jps,jpe,ims,ime,jms,jme,grnqfx,'grnqfx',pid)
1196 call write_array_m(ifps,ifpe,jfps,jfpe,ifms,ifme,jfms,jfme,fuel_frac,'fuel_frac',pid)
1197 call write_array_m(ifps,ifpe,jfps,jfpe,ifms,ifme,jfms,jfme,fgrnhfx,'fgrnhfx',pid)
1198 call write_array_m(ifps,ifpe,jfps,jfpe,ifms,ifme,jfms,jfme,fgrnqfx,'fgrnqfx',pid)
1199 endif
1200 #endif
1201 endif
1203 end subroutine sfire_driver_phys
1205 !
1206 !***
1207 !
1209 subroutine check_fmesh(ids,ide,ifds,ifde,ir,s)
1210 !*** purpose: check if fire and atm meshes line up
1211 implicit none
1212 !*** arguments
1213 integer, intent(in)::ids,ide,ifds,ifde,ir
1214 character(len=*),intent(in)::s
1215 !*** local
1216 character(len=128)msg
1217 !*** executable
1218 if ((ide-ids+1)*ir.ne.(ifde-ifds+1))then
1219 !$OMP CRITICAL(SFIRE_DRIVER_CRIT)
1220 write(msg,1)s,ids,ide,ifds,ifde,ir
1221 1 format('module_fr_sfire_driver: incompatible bounds ',a,' atm ',i5,':',i5,' fire ',i5,':',i5,' ratio ',i3)
1222 !$OMP END CRITICAL(SFIRE_DRIVER_CRIT)
1223 call crash(msg)
1224 endif
1225 end subroutine check_fmesh
1227 subroutine fire_ignition_convert (config_flags,ignition, &
1228 fxlong, fxlat, &
1229 ifds,ifde, jfds,jfde, &
1230 ifms,ifme, jfms,jfme, &
1231 ifps,ifpe, jfps,jfpe )
1232 implicit none
1233 ! create ignition arrays from scalar flags
1234 !*** arguments
1235 TYPE (grid_config_rec_type) , INTENT(IN) :: config_flags
1236 TYPE (lines_type), INTENT(OUT):: ignition ! any values from input discarded
1237 integer::ifds,ifde, jfds,jfde, &
1238 ifms,ifme, jfms,jfme, &
1239 ifps,ifpe, jfps,jfpe
1240 real, dimension(ifms:ifme,jfms:jfme):: fxlong,fxlat
1241 !*** local
1242 integer::i,j,ii,jj
1243 logical:: real,ideal
1244 character(len=128)msg
1245 real:: corner_longlat(2,2,2), corner_longlat_1(8), corner_longlat_2(8),lon(2),lat(2)
1246 real, dimension(2,2):: corner_long,corner_lat ! coordinates of fire mesh corner cells
1248 !*** executable
1251 ignition%max_lines=5 ! number of lines that have entries in the namelist
1252 ignition%num_lines=config_flags%fire_num_ignitions
1254 ! this is only until I figure out how to input arrays through the namelist...
1255 if(fire_max_lines.lt.ignition%max_lines)call crash('fire_max_lines too small')
1257 ! figure out which kind of coordinates from the first given
1258 ideal=config_flags%fire_ignition_start_x1 .ne.0. .or. config_flags%fire_ignition_start_y1 .ne. 0.
1259 real=config_flags%fire_ignition_start_lon1 .ne. 0. .or. config_flags%fire_ignition_start_lat1 .ne. 0.
1260 if(ideal)call message('Using ideal ignition coordinates, m from the lower left domain corner')
1261 if(real)call message('Using real ignition coordinates, longitude and latitude')
1262 if(ideal.and.real)call crash('Only one of the ideal or real coordinates may be given')
1264 ignition%longlat=0 ! default, if no ignition
1265 if(ideal)then
1266 ! use values from _x and _y variables
1267 ignition%longlat=0
1268 ignition%line(1)%start_x=config_flags%fire_ignition_start_x1
1269 ignition%line(1)%start_y=config_flags%fire_ignition_start_y1
1270 ignition%line(1)%end_x=config_flags%fire_ignition_end_x1
1271 ignition%line(1)%end_y=config_flags%fire_ignition_end_y1
1272 ignition%line(2)%start_x=config_flags%fire_ignition_start_x2
1273 ignition%line(2)%start_y=config_flags%fire_ignition_start_y2
1274 ignition%line(2)%end_x=config_flags%fire_ignition_end_x2
1275 ignition%line(2)%end_y=config_flags%fire_ignition_end_y2
1276 ignition%line(3)%start_x=config_flags%fire_ignition_start_x3
1277 ignition%line(3)%start_y=config_flags%fire_ignition_start_y3
1278 ignition%line(3)%end_x=config_flags%fire_ignition_end_x3
1279 ignition%line(3)%end_y=config_flags%fire_ignition_end_y3
1280 ignition%line(4)%start_x=config_flags%fire_ignition_start_x4
1281 ignition%line(4)%start_y=config_flags%fire_ignition_start_y4
1282 ignition%line(4)%end_x=config_flags%fire_ignition_end_x4
1283 ignition%line(4)%end_y=config_flags%fire_ignition_end_y4
1284 ignition%line(5)%start_x=config_flags%fire_ignition_start_x5
1285 ignition%line(5)%start_y=config_flags%fire_ignition_start_y5
1286 ignition%line(5)%end_x=config_flags%fire_ignition_end_x5
1287 ignition%line(5)%end_y=config_flags%fire_ignition_end_y5
1288 endif
1289 if(real)then
1290 ! use values from _long and _lat
1291 ignition%longlat=1
1292 ignition%line(1)%start_x=config_flags%fire_ignition_start_lon1
1293 ignition%line(1)%start_y=config_flags%fire_ignition_start_lat1
1294 ignition%line(1)%end_x=config_flags%fire_ignition_end_lon1
1295 ignition%line(1)%end_y=config_flags%fire_ignition_end_lat1
1296 ignition%line(2)%start_x=config_flags%fire_ignition_start_lon2
1297 ignition%line(2)%start_y=config_flags%fire_ignition_start_lat2
1298 ignition%line(2)%end_x=config_flags%fire_ignition_end_lon2
1299 ignition%line(2)%end_y=config_flags%fire_ignition_end_lat2
1300 ignition%line(3)%start_x=config_flags%fire_ignition_start_lon3
1301 ignition%line(3)%start_y=config_flags%fire_ignition_start_lat3
1302 ignition%line(3)%end_x=config_flags%fire_ignition_end_lon3
1303 ignition%line(3)%end_y=config_flags%fire_ignition_end_lat3
1304 ignition%line(4)%start_x=config_flags%fire_ignition_start_lon4
1305 ignition%line(4)%start_y=config_flags%fire_ignition_start_lat4
1306 ignition%line(4)%end_x=config_flags%fire_ignition_end_lon4
1307 ignition%line(4)%end_y=config_flags%fire_ignition_end_lat4
1308 ignition%line(5)%start_x=config_flags%fire_ignition_start_lon5
1309 ignition%line(5)%start_y=config_flags%fire_ignition_start_lat5
1310 ignition%line(5)%end_x=config_flags%fire_ignition_end_lon5
1311 ignition%line(5)%end_y=config_flags%fire_ignition_end_lat5
1312 endif
1313 ! common to both cases
1314 ignition%line(1)%ros=config_flags%fire_ignition_ros1
1315 ignition%line(1)%radius=config_flags%fire_ignition_radius1
1316 ignition%line(1)%start_time=config_flags%fire_ignition_start_time1
1317 ignition%line(1)%end_time=config_flags%fire_ignition_end_time1
1318 ignition%line(2)%ros=config_flags%fire_ignition_ros2
1319 ignition%line(2)%radius=config_flags%fire_ignition_radius2
1320 ignition%line(2)%start_time=config_flags%fire_ignition_start_time2
1321 ignition%line(2)%end_time=config_flags%fire_ignition_end_time2
1322 ignition%line(3)%ros=config_flags%fire_ignition_ros3
1323 ignition%line(3)%radius=config_flags%fire_ignition_radius3
1324 ignition%line(3)%start_time=config_flags%fire_ignition_start_time3
1325 ignition%line(3)%end_time=config_flags%fire_ignition_end_time3
1326 ignition%line(4)%ros=config_flags%fire_ignition_ros4
1327 ignition%line(4)%radius=config_flags%fire_ignition_radius4
1328 ignition%line(4)%start_time=config_flags%fire_ignition_start_time4
1329 ignition%line(4)%end_time=config_flags%fire_ignition_end_time4
1330 ignition%line(5)%ros=config_flags%fire_ignition_ros5
1331 ignition%line(5)%radius=config_flags%fire_ignition_radius5
1332 ignition%line(5)%start_time=config_flags%fire_ignition_start_time5
1333 ignition%line(5)%end_time=config_flags%fire_ignition_end_time5
1335 call postprocess_lines(ignition,'ros',config_flags)
1337 ! get the coordinates of the corner cells
1338 corner_longlat=0.
1339 if(ifds.eq.ifps.and.jfds.eq.jfps)then
1340 corner_longlat(1,1,1)=fxlong(ifps,jfps)
1341 corner_longlat(1,1,2)=fxlat(ifps,jfps)
1342 endif
1343 if(ifds.eq.ifps.and.jfde.eq.jfpe)then
1344 corner_longlat(1,2,1)=fxlong(ifps,jfpe)
1345 corner_longlat(1,2,2)=fxlat(ifps,jfpe)
1346 endif
1347 if(ifde.eq.ifpe.and.jfds.eq.jfps)then
1348 corner_longlat(2,1,1)=fxlong(ifpe,jfps)
1349 corner_longlat(2,1,2)=fxlat(ifpe,jfps)
1350 endif
1351 if(ifde.eq.ifpe.and.jfde.eq.jfpe)then
1352 corner_longlat(2,2,1)=fxlong(ifpe,jfpe)
1353 corner_longlat(2,2,2)=fxlat(ifpe,jfpe)
1354 endif
1355 corner_longlat_1=reshape(corner_longlat,(/8/))
1356 #ifdef DM_PARALLEL
1357 call wrf_dm_sum_reals(corner_longlat_1,corner_longlat_2)
1358 #else
1359 corner_longlat_2=corner_longlat_1
1360 #endif
1361 corner_longlat=reshape(corner_longlat_2,(/2,2,2/))
1362 corner_long=corner_longlat(1:2,1:2,1)
1363 corner_lat=corner_longlat(1:2,1:2,2)
1364 if(fire_print_msg.ge.2)then
1365 do i=1,2
1366 do j=1,2
1367 write(msg,'(a,2i2,a,2f14.8)')'corner',i,j,' coordinates ',corner_long(i,j),corner_lat(i,j)
1368 call message(msg)
1369 enddo
1370 enddo
1371 endif
1372 lon(1)=(corner_long(1,1)+corner_long(1,2))/2.
1373 lon(2)=(corner_long(2,1)+corner_long(2,2))/2.
1374 lat(1)=(corner_lat(1,1)+corner_lat(2,1))/2.
1375 lat(2)=(corner_lat(1,2)+corner_lat(2,2))/2.
1376 if(fire_print_msg.ge.2)then
1377 write(msg,'(4(a,f14.8))')'coordinates ',lon(1),':',lon(2),',',lat(1),':',lat(2)
1378 call message(msg)
1379 endif
1381 do i=1,ignition%num_lines
1382 call check_ignition_coordinate(ignition%line(i)%start_x,lon(1),lon(2))
1383 call check_ignition_coordinate(ignition%line(i)%start_y,lat(1),lat(2))
1384 call check_ignition_coordinate(ignition%line(i)%end_x,lon(1),lon(2))
1385 call check_ignition_coordinate(ignition%line(i)%end_y,lat(1),lat(2))
1386 enddo
1388 if (fire_ignition_clamp>0) then
1389 do i=1,ignition%num_lines
1390 call clamp_to_grid(ignition%line(i)%start_x,lon(1),lon(2),ifds,ifde,ignition%line(i)%start_x,ii)
1391 call clamp_to_grid(ignition%line(i)%start_y,lat(1),lat(2),jfds,jfde,ignition%line(i)%start_y,jj)
1392 call display_clamp
1393 call clamp_to_grid(ignition%line(i)%end_x,lon(1),lon(2),ifds,ifde,ignition%line(i)%end_x,ii)
1394 call clamp_to_grid(ignition%line(i)%end_y,lat(1),lat(2),jfds,jfde,ignition%line(i)%end_y,jj)
1395 call display_clamp
1396 ! for now, ii jj ignored. In future replace by fxlong(ii,jj), fxlat(ii,jj) to guard against rounding
1397 enddo
1398 endif
1399 contains
1400 subroutine display_clamp
1401 character(len=128)::msg
1402 real::d1,d2
1403 if(ii>=ifps.and.ii<=ifpe.and.jj>=jfps.and.jj<=jfpe)then
1404 write(msg,'(a,2f14.8,a,2i6)')'grid node ',fxlong(ii,jj),fxlat(ii,jj),' index',ii,jj
1405 call message(msg)
1406 endif
1407 end subroutine display_clamp
1408 end subroutine fire_ignition_convert
1410 subroutine check_ignition_coordinate(x,x1,x2)
1411 !*** arguments
1412 real, intent(in)::x,x1,x2
1413 character(len=128)::msg
1414 if (.not.(x>x1 .and. x<x2))then
1415 write(msg,'(a,f14.8,a,2f14.8)')'ignition point coordinate',x,' must be inside the bounds',x1,x2
1416 call crash(msg)
1417 endif
1418 end subroutine check_ignition_coordinate
1420 subroutine clamp_to_grid(x,x1,x2,i1,i2,xout,iout)
1421 !*** round point on uniform mesh to a mesh node
1422 !*** arguments
1423 real, intent(in)::x,x1,x2
1424 integer, intent(in)::i1,i2
1425 real, intent(out)::xout
1426 integer, intent(out)::iout
1427 !*** local
1428 character(len=128)::msg
1429 integer:: i
1430 real::r,dx,xr
1431 !*** executable
1432 dx=(x2-x1)/(i2-i1)
1433 r=i1+(x-x1)/dx
1434 iout=nint(r)
1435 xr=x1+(iout-i1)*dx
1436 if(fire_print_msg.ge.2)then
1437 write(msg,'(a,f14.8,a,f14.8,a,i6)')'coordinate ',x,' clamped to ',xr,' index',iout
1438 call message(msg)
1439 endif
1440 xout=xr
1441 end subroutine clamp_to_grid
1442 !
1443 !***
1444 !
1445 subroutine postprocess_lines(lines,value_name,config_flags)
1446 ! count lines, fill endpoints, print stats
1447 !*** arguments
1448 type(lines_type), intent(inout)::lines
1449 character(len=3), intent(in)::value_name
1450 TYPE (grid_config_rec_type) , INTENT(IN) :: config_flags ! namelist
1451 !*** local
1452 integer::i,n
1453 real::value
1454 real::lat_ctr,lon_ctr
1455 character(len=128)msg,f2,f3
1456 !*** executable
1458 n=lines%num_lines
1459 do i=1,n
1460 ! find the last line that has positive radius and reset num_lines if needed
1461 if(lines%line(i)%radius.gt.0.)lines%num_lines=i
1462 ! expand ignition data given as zero
1463 if(lines%line(i)%end_x.eq.0.)lines%line(i)%end_x=lines%line(i)%start_x
1464 if(lines%line(i)%end_y.eq.0.)lines%line(i)%end_y=lines%line(i)%start_y
1465 if(lines%line(i)%end_time.eq.0.)lines%line(i)%end_time=lines%line(i)%start_time
1466 enddo
1468 if(lines%longlat .eq. 0)then
1469 ! ideal
1470 ! ignition is in m
1471 lines%unit_fxlong=1.
1472 lines%unit_fxlat=1.
1473 ! will set fire mesh coordinates to uniform mesh below
1474 else
1475 ! real
1476 lat_ctr=config_flags%cen_lat
1477 lon_ctr=config_flags%cen_lon
1478 ! 1 degree in m (approximate OK)
1479 lines%unit_fxlat=pi2/(360.*reradius) ! earth circumference in m / 360 degrees
1480 lines%unit_fxlong=cos(lat_ctr*pi2/360.)*lines%unit_fxlat ! latitude
1481 endif
1483 if(fire_print_msg.ge.2)then
1484 !$OMP CRITICAL(SFIRE_DRIVER_CRIT)
1485 if(lines%longlat .eq. 0)then
1486 write(msg,1)'start x','start y','end x','end y','start t','end t',value_name,'radius'
1487 else
1488 write(msg,1)'start lon','start lat','end lon','end lat','start time','end time',value_name,'radius'
1489 endif
1490 1 format(4a10,4a9)
1491 call message(msg)
1492 do i=1,lines%num_lines
1493 select case (value_name)
1494 case ('ros')
1495 value = lines%line(i)%ros
1496 f2='(4f10.3,4f9.2)'
1497 f3='(4f10.5,4f9.2)'
1498 case ('hfx')
1499 value = lines%line(i)%hfx_value
1500 f2='(4f10.3,2f9.2,e9.2,f9.2)'
1501 f3='(4f10.5,2f9.2,e9.2,f9.2)'
1502 case default
1503 call crash('postprocess_lines: bad value_name '//value_name)
1504 end select
1505 if(lines%longlat .eq. 0)then
1506 write(msg,f2)lines%line(i)%start_x, lines%line(i)%start_y, &
1507 lines%line(i)%end_x, lines%line(i)%end_y, &
1508 lines%line(i)%start_time, lines%line(i)%end_time, &
1509 value, lines%line(i)%radius
1510 else
1511 write(msg,f3)lines%line(i)%start_x, lines%line(i)%start_y, &
1512 lines%line(i)%end_x, lines%line(i)%end_y, &
1513 lines%line(i)%start_time, lines%line(i)%end_time, &
1514 value, lines%line(i)%radius
1515 endif
1516 call message(msg)
1517 if(lines%line(i)%start_time > lines%line(i)%end_time)then
1518 call crash('start time may not be after end time')
1519 endif
1520 enddo
1521 !$OMP END CRITICAL(SFIRE_DRIVER_CRIT)
1522 endif
1523 end subroutine postprocess_lines
1526 subroutine fire_hfx_convert (config_flags,hfx)
1527 implicit none
1528 ! create heat flux line(s) from scalar flags
1529 !*** arguments
1530 TYPE (grid_config_rec_type) , INTENT(IN) :: config_flags
1531 TYPE (lines_type), INTENT(OUT):: hfx ! any values from input discarded
1532 !*** local
1533 integer::i
1534 logical:: real,ideal
1535 real::lat_ctr,lon_ctr
1536 character(len=128)msg
1537 !*** executable
1538 ! this is only until I figure out how to input arrays through the namelist...
1539 hfx%num_lines=config_flags%fire_hfx_num_lines
1540 if(fire_max_lines.lt.hfx%num_lines)call crash('fire_max_lines too small')
1542 ! figure out which kind of coordinates from the first given
1543 ideal=config_flags%fire_hfx_start_x1 .ne.0. .or. config_flags%fire_hfx_start_y1 .ne. 0.
1544 real=config_flags%fire_hfx_start_lon1 .ne. 0. .or. config_flags%fire_hfx_start_lat1 .ne. 0.
1545 if(ideal)call message('Using ideal heat flux line coordinates, m from the lower left domain corner')
1546 if(real)call message('Using real heat flux line coordinates, longitude and latitude')
1547 if(ideal.and.real)call crash('Only one of the ideal or real coordinates may be given')
1549 hfx%longlat=0 ! default, if no ignition
1550 if(ideal)then
1551 ! use values from _x and _y variables
1552 hfx%longlat=0
1553 hfx%line(1)%start_x=config_flags%fire_hfx_start_x1
1554 hfx%line(1)%start_y=config_flags%fire_hfx_start_y1
1555 hfx%line(1)%end_x=config_flags%fire_hfx_end_x1
1556 hfx%line(1)%end_y=config_flags%fire_hfx_end_y1
1557 endif
1558 if(real)then
1559 ! use values from _long and _lat
1560 hfx%longlat=1
1561 hfx%line(1)%start_x=config_flags%fire_hfx_start_lon1
1562 hfx%line(1)%start_y=config_flags%fire_hfx_start_lat1
1563 hfx%line(1)%end_x=config_flags%fire_hfx_end_lon1
1564 hfx%line(1)%end_y=config_flags%fire_hfx_end_lat1
1565 endif
1566 ! common to both cases
1567 hfx%line(1)%radius=config_flags%fire_hfx_radius1
1568 hfx%line(1)%start_time=config_flags%fire_hfx_start_time1
1569 hfx%line(1)%end_time=config_flags%fire_hfx_end_time1
1570 hfx%line(1)%trans_time=config_flags%fire_hfx_trans_time1
1571 hfx%line(1)%hfx_value=config_flags%fire_hfx_value1
1573 call postprocess_lines(hfx,'hfx',config_flags)
1575 end subroutine fire_hfx_convert
1577 subroutine set_flags(config_flags)
1578 USE module_configure
1579 implicit none
1580 TYPE (grid_config_rec_type) , INTENT(IN) :: config_flags
1581 ! copy flags from wrf to module_fr_sfire_util
1582 ! for instructions how to add a flag see the top of module_fr_sfire_util.F
1584 fire_perimeter_time = config_flags%fire_perimeter_time
1585 fire_print_msg = config_flags%fire_print_msg
1586 fire_print_file = config_flags%fire_print_file
1587 fuel_left_method = config_flags%fire_fuel_left_method
1588 fuel_left_irl = config_flags%fire_fuel_left_irl
1589 fuel_left_jrl = config_flags%fire_fuel_left_jrl
1590 fire_atm_feedback = config_flags%fire_atm_feedback
1591 fire_hfx_given = config_flags%fire_hfx_given
1592 fire_hfx_num_lines = config_flags%fire_hfx_num_lines
1593 fire_hfx_latent_part = config_flags%fire_hfx_latent_part
1594 boundary_guard = config_flags%fire_boundary_guard
1595 fire_back_weight = config_flags%fire_back_weight
1596 fire_grows_only = config_flags%fire_grows_only
1597 fire_upwinding = config_flags%fire_upwinding
1598 fire_viscosity = config_flags%fire_viscosity
1599 fire_lfn_ext_up = config_flags%fire_lfn_ext_up
1600 fire_test_steps = config_flags%fire_test_steps
1601 !fire_topo_from_atm = config_flags%fire_topo_from_atm
1602 fire_advection = config_flags%fire_advection
1603 fire_wind_log_interp = config_flags%fire_wind_log_interp
1604 fire_use_windrf = config_flags%fire_use_windrf
1605 fire_fmc_read = config_flags%fire_fmc_read
1606 fire_ignition_clamp = config_flags%fire_ignition_clamp
1607 kfmc_ndwi = config_flags%kfmc_ndwi
1608 fndwi_from_ndwi = config_flags%fndwi_from_ndwi
1610 end subroutine set_flags
1612 !
1613 !*****************************
1614 !
1617 subroutine set_fp_from_grid(grid,fp)
1618 implicit none
1619 type(domain),intent(in)::grid
1620 type(fire_params),intent(out)::fp
1622 ! pointers to be passed to fire spread formulas
1623 fp%vx => grid%uf ! fire winds
1624 fp%vy => grid%vf ! fire winds
1625 fp%zsf => grid%zsf ! terrain height
1626 fp%dzdxf => grid%dzdxf ! terrain grad
1627 fp%dzdyf => grid%dzdyf ! terrain grad
1628 fp%bbb => grid%bbb ! spread formula coeff
1629 fp%phisc => grid%phisc ! spread formula coeff
1630 fp%phiwc => grid%phiwc ! spread formula coeff
1631 fp%r_0 => grid%r_0 ! spread formula coeff
1632 fp%fgip => grid%fgip ! spread formula coeff
1633 fp%ischap => grid%ischap ! spread formula coeff
1634 fp%fuel_time => grid%fuel_time ! time for fuel to burn to 1/e
1635 fp%fmc_g => grid%fmc_g ! fuel moisture, ground
1636 fp%nfuel_cat => grid%nfuel_cat ! fuel category
1639 end subroutine set_fp_from_grid
1642 subroutine print_id
1643 character(len=128)::id,msg
1644 #include "sfire_id.inc"
1645 msg=id
1646 call message(msg,level=1)
1647 end subroutine print_id
1649 end module module_fr_sfire_driver