Add restart capability to the fire behavior model

CMakeLists.txt

@@ -67,6 +67,7 @@ endif()
 
 # State
 list(APPEND _state_files state/state_mod.F90
+                         state/state_restart_submod.F90
 			 state/ignition_line_mod.F90
 	                 state/tiles_mod.F90)
 

driver/fire_behavior.F90

@@ -8,6 +8,7 @@ program fire_behavior
     use initialize_mod, only : Init_fire_state, Init_atm_state
     use advance_mod, only : Advance_state
     use wrfdata_mod, only : wrfdata_t
+    use stderrout_mod, only : Print_message
     use, intrinsic :: iso_fortran_env, only : ERROR_UNIT, OUTPUT_UNIT
 
     implicit none
@@ -71,16 +72,26 @@ program fire_behavior
 
     end select
 
-    if (DEBUG_LOCAL) write (OUTPUT_UNIT, *) '  Saving fire state...'
-    call grid%Save_state ()
+    if (.not. config_flags%restart) then
+      if (DEBUG_LOCAL) call Print_message ('    Saving fire state')
+      call grid%Save_state ()
+    end if
+
+    call grid%Handle_restart (config_flags, initialize = .true.)
 
     if (DEBUG_LOCAL) write (OUTPUT_UNIT, *) '  Starting temporal loop...'
     do while (grid%datetime_now < grid%datetime_end)
       call Advance_state (grid, config_flags)
+
+      call grid%Check_model_clock ()
+
       ! Save the forcing used for the completed interval, as in the NUOPC cap.
       ! Refresh atmospheric fields afterward for the next fire advance.
+
       call grid%Handle_output (config_flags)
       if (config_flags%ideal_opt == 0) call grid%Handle_wrfdata_update (atm_state, config_flags)
+
+      call grid%Handle_restart (config_flags)
     end do
     if (DEBUG_LOCAL) write (OUTPUT_UNIT, *) '  Completed temporal loop'
 

driver/initialize_mod.F90

@@ -63,46 +63,52 @@ subroutine Init_fire_state (grid, config_flags, wrf)
       if (DEBUG_LOCAL) call Print_message ('  Initialization...')
       if (config_flags%ideal_opt == 0) then
           ! Real world
-        if (DEBUG_LOCAL) call Print_message ('    Reading geogrid file')
+        if (.not. config_flags%restart) then
+          if (DEBUG_LOCAL) call Print_message ('    Reading geogrid file')
 #ifdef DM_PARALLEL
-        if (grid%is_cfbm_comm_set) then
-          call Convert_mpi_comm_to_f08 (grid%cfbm_comm, cfbm_comm_f08)
-          call Mpi_comm_rank (cfbm_comm_f08, rank, ierr)
-        else
-          call Stop_simulation ('The MPI communicator cfbm_comm has not been set')
-        end if
+          if (grid%is_cfbm_comm_set) then
+            call Convert_mpi_comm_to_f08 (grid%cfbm_comm, cfbm_comm_f08)
+            call Mpi_comm_rank (cfbm_comm_f08, rank, ierr)
+          else
+            call Stop_simulation ('The MPI communicator cfbm_comm has not been set')
+          end if
 #else
-        rank = 0
+          rank = 0
 #endif
 
-        if (rank == 0) geogrid = geogrid_t (file_name = 'geo_em.d01.nc')
+          if (rank == 0) geogrid = geogrid_t (file_name = 'geo_em.d01.nc')
 
 #ifdef DM_PARALLEL
-        call MPI_Bcast (geogrid%cen_lon, 1, MPI_REAL, 0, cfbm_comm_f08, ierr)
-        call MPI_Bcast (geogrid%cen_lat, 1, MPI_REAL, 0, cfbm_comm_f08, ierr)
-        call MPI_Bcast (geogrid%dx, 1, MPI_REAL, 0, cfbm_comm_f08, ierr)
-        call MPI_Bcast (geogrid%dy, 1, MPI_REAL, 0, cfbm_comm_f08, ierr)
-        call MPI_Bcast (geogrid%true_lat_1, 1, MPI_REAL, 0, cfbm_comm_f08, ierr)
-        call MPI_Bcast (geogrid%true_lat_2, 1, MPI_REAL, 0, cfbm_comm_f08, ierr)
-        call MPI_Bcast (geogrid%stand_lon, 1, MPI_REAL, 0, cfbm_comm_f08, ierr)
-
-        call MPI_Bcast(geogrid%map_proj, 1, MPI_INTEGER, 0, cfbm_comm_f08, ierr)
-        call MPI_Bcast(geogrid%sr_x, 1, MPI_INTEGER, 0, cfbm_comm_f08, ierr)
-        call MPI_Bcast(geogrid%sr_y, 1, MPI_INTEGER, 0, cfbm_comm_f08, ierr)
-
-        call MPI_Bcast(geogrid%ifds, 1, MPI_INTEGER, 0, cfbm_comm_f08, ierr)
-        call MPI_Bcast(geogrid%ifde, 1, MPI_INTEGER, 0, cfbm_comm_f08, ierr)
-        call MPI_Bcast(geogrid%jfds, 1, MPI_INTEGER, 0, cfbm_comm_f08, ierr)
-        call MPI_Bcast(geogrid%jfde, 1, MPI_INTEGER, 0, cfbm_comm_f08, ierr)
-
-        call MPI_Bcast(geogrid%ids, 1, MPI_INTEGER, 0, cfbm_comm_f08, ierr)
-        call MPI_Bcast(geogrid%ide, 1, MPI_INTEGER, 0, cfbm_comm_f08, ierr)
-        call MPI_Bcast(geogrid%jds, 1, MPI_INTEGER, 0, cfbm_comm_f08, ierr)
-        call MPI_Bcast(geogrid%jde, 1, MPI_INTEGER, 0, cfbm_comm_f08, ierr)
+          call MPI_Bcast (geogrid%cen_lon, 1, MPI_REAL, 0, cfbm_comm_f08, ierr)
+          call MPI_Bcast (geogrid%cen_lat, 1, MPI_REAL, 0, cfbm_comm_f08, ierr)
+          call MPI_Bcast (geogrid%dx, 1, MPI_REAL, 0, cfbm_comm_f08, ierr)
+          call MPI_Bcast (geogrid%dy, 1, MPI_REAL, 0, cfbm_comm_f08, ierr)
+          call MPI_Bcast (geogrid%true_lat_1, 1, MPI_REAL, 0, cfbm_comm_f08, ierr)
+          call MPI_Bcast (geogrid%true_lat_2, 1, MPI_REAL, 0, cfbm_comm_f08, ierr)
+          call MPI_Bcast (geogrid%stand_lon, 1, MPI_REAL, 0, cfbm_comm_f08, ierr)
+
+          call MPI_Bcast(geogrid%map_proj, 1, MPI_INTEGER, 0, cfbm_comm_f08, ierr)
+          call MPI_Bcast(geogrid%sr_x, 1, MPI_INTEGER, 0, cfbm_comm_f08, ierr)
+          call MPI_Bcast(geogrid%sr_y, 1, MPI_INTEGER, 0, cfbm_comm_f08, ierr)
+
+          call MPI_Bcast(geogrid%ifds, 1, MPI_INTEGER, 0, cfbm_comm_f08, ierr)
+          call MPI_Bcast(geogrid%ifde, 1, MPI_INTEGER, 0, cfbm_comm_f08, ierr)
+          call MPI_Bcast(geogrid%jfds, 1, MPI_INTEGER, 0, cfbm_comm_f08, ierr)
+          call MPI_Bcast(geogrid%jfde, 1, MPI_INTEGER, 0, cfbm_comm_f08, ierr)
+
+          call MPI_Bcast(geogrid%ids, 1, MPI_INTEGER, 0, cfbm_comm_f08, ierr)
+          call MPI_Bcast(geogrid%ide, 1, MPI_INTEGER, 0, cfbm_comm_f08, ierr)
+          call MPI_Bcast(geogrid%jds, 1, MPI_INTEGER, 0, cfbm_comm_f08, ierr)
+          call MPI_Bcast(geogrid%jde, 1, MPI_INTEGER, 0, cfbm_comm_f08, ierr)
 #endif
+        end if
 
         if (DEBUG_LOCAL) call Print_message ('    Initializing fire state')
-        call grid%Initialization (config_flags, geogrid)
+        if (config_flags%restart) then
+          call grid%Initialization (config_flags)
+        else
+          call grid%Initialization (config_flags, geogrid)
+        end if
 
       else
           ! Ideal
@@ -112,11 +118,16 @@ subroutine Init_fire_state (grid, config_flags, wrf)
 
       call Init_fire_components (grid, config_flags)
 
-      if (present (wrf)) then
+      if (present (wrf) .and. .not. config_flags%restart) then
         if (DEBUG_LOCAL) call Print_message ('    Initializing atmospheric state')
         call grid%Handle_wrfdata_update (wrf, config_flags)
       end if
 
+      if (config_flags%restart) then
+        if (DEBUG_LOCAL) call Print_message ('    Reading restart state')
+        call grid%Read_restart (config_flags)
+      end if
+
       if (DEBUG_LOCAL) then
           ! print lat/lons
         open (newunit = unit_out, file = 'latlons_c.dat')
@@ -187,6 +198,8 @@ subroutine Init_fire_state_within_wrf (state, config_flags, &
 
       call Init_fire_components (state, config_flags)
 
+      if (config_flags%restart) call state%Read_restart (config_flags)
+
       if (DEBUG_LOCAL) call Print_message ('  Leaving subroutine Init_fire_state_within_wrf...')
 
     end subroutine Init_fire_state_within_wrf

io/namelist_mod.F90

@@ -38,6 +38,10 @@ module namelist_mod
       integer :: interval_output = -1         ! Frequency to save the output [s]
       real :: dt = 2.0                        ! Time step of the fire model [s]
 
+        ! Restart controls
+      logical :: restart = .false.            ! Read initial dynamic fire state from fire_restart_<start datetime>.nc
+      integer :: restart_interval = -1        ! Restart output interval [s]; -1 disables restart writes
+
       integer :: num_tiles = 1                ! Number of tiles in each patch
       integer :: tile_strategy = 0            ! Strategy for the tile decomposition: 0) ...
 
@@ -216,6 +220,8 @@ subroutine Broadcast_nml (this, mpi_comm_cfbm)
 
       call Broadcast_integer (this%ideal_opt)
       call Broadcast_integer (this%devel_opt)
+      call Broadcast_logical (this%restart)
+      call Broadcast_integer (this%restart_interval)
       call Broadcast_integer (this%fuel_opt)
       call Broadcast_integer (this%ros_opt)
       call Broadcast_integer (this%fmc_opt)
@@ -410,6 +416,9 @@ subroutine Check_nml (this, require_atm_interval)
             call Stop_simulation ('perimeter ignition time must align with a fire timestep')
       end if
 
+      if (this%restart_interval /= -1 .and. this%restart_interval <= 0) &
+          call Stop_simulation ('restart_interval must be positive or -1')
+
     end subroutine Check_nml
 
     subroutine Check_time_intervals (this, require_atm_interval)
@@ -829,15 +838,16 @@ subroutine Init_time_block (this, file_name)
 
       integer :: start_year, start_month, start_day, start_hour, start_minute, start_second, &
           end_year, end_month, end_day, end_hour, end_minute, end_second, interval_output, &
-          num_tiles, tile_strategy
+          restart_interval, num_tiles, tile_strategy
       real :: dt
+      logical :: restart
 
       character (len = :), allocatable :: msg
       integer :: unit_nml, io_stat
 
       namelist /time/ start_year, start_month, start_day, start_hour, start_minute, start_second, &
           end_year, end_month, end_day, end_hour, end_minute, end_second, dt, interval_output, &
-          num_tiles, tile_strategy
+          restart, restart_interval, num_tiles, tile_strategy
 
 
         ! Set default values
@@ -855,6 +865,8 @@ subroutine Init_time_block (this, file_name)
       end_second = this%end_second
       dt = this%dt
       interval_output = this%interval_output
+      restart = this%restart
+      restart_interval = this%restart_interval
 
       num_tiles = this%num_tiles
       tile_strategy = this%tile_strategy
@@ -886,6 +898,8 @@ subroutine Init_time_block (this, file_name)
       this%end_second = end_second
       this%dt = dt
       this%interval_output = interval_output
+      this%restart = restart
+      this%restart_interval = restart_interval
 
       this%num_tiles = num_tiles
       this%tile_strategy = tile_strategy