m_data_input.f90

Go to the documentation of this file.

!>
!! @file
!> @brief Contains module m_data_input
 
!> @brief Reads raw simulation grid and conservative-variable data for a given time-step and fills buffer regions
module m_data_input
 
#ifdef MFC_MPI
    use mpi                     !< message passing interface (mpi) module
#endif
 
    use m_derived_types         !< definitions of the derived types
 
    use m_global_parameters     !< global parameters for the code
 
    use m_mpi_proxy             !< message passing interface (mpi) module proxy
 
    use m_mpi_common
 
    use m_compile_specific
 
    use m_boundary_common
 
    use m_helper
 
    implicit none
 
    private; public :: s_initialize_data_input_module, &
 s_read_data_files, &
 s_read_serial_data_files, &
 s_read_parallel_data_files, &
 s_finalize_data_input_module
 
    abstract interface
 
        !> Subroutine for reading data files
        !!  @param t_step Current time-step to input
        impure subroutine s_read_abstract_data_files(t_step)
 
            implicit none
 
            integer, intent(in) :: t_step
 
        end subroutine s_read_abstract_data_files
 
    end interface
 
    type(scalar_field), allocatable, dimension(:), public :: q_cons_vf !<
    !! Conservative variables
 
    type(scalar_field), allocatable, dimension(:), public :: q_cons_temp
 
    type(scalar_field), allocatable, dimension(:), public :: q_prim_vf !<
    !! Primitive variables
 
    type(integer_field), allocatable, dimension(:, :), public :: bc_type !<
    !! Boundary condition identifiers
 
    type(scalar_field), public :: q_t_sf !<
    !! Temperature field
 
    ! type(scalar_field), public :: ib_markers !<
    type(integer_field), public :: ib_markers
 
    procedure(s_read_abstract_data_files), pointer :: s_read_data_files => null()
 
contains
 
    !> Helper subroutine to read grid data files for a given direction
    !!  @param t_step_dir Directory containing the time-step data
    !!  @param direction Direction name ('x', 'y', 'z')
    !!  @param cb_array Cell boundary array to populate
    !!  @param d_array Cell width array to populate
    !!  @param cc_array Cell center array to populate
    !!  @param size_dim Size of the dimension
    impure subroutine s_read_grid_data_direction(t_step_dir, direction, cb_array, d_array, cc_array, size_dim)
 
        character(len=*), intent(in) :: t_step_dir
        character(len=1), intent(in) :: direction
        real(wp), dimension(-1:), intent(out) :: cb_array
        real(wp), dimension(0:), intent(out) :: d_array
        real(wp), dimension(0:), intent(out) :: cc_array
        integer, intent(in) :: size_dim
 
        character(LEN=len_trim(t_step_dir) + 10) :: file_loc
        logical :: file_check
 
        ! Checking whether direction_cb.dat exists
        file_loc = trim(t_step_dir)//'/'//direction//'_cb.dat'
        inquire (file=trim(file_loc), exist=file_check)
 
        ! Reading direction_cb.dat if it exists, exiting otherwise
        if (file_check) then
            open (1, file=trim(file_loc), form='unformatted', &
                  status='old', action='read')
            read (1) cb_array(-1:size_dim)
            close (1)
        else
            call s_mpi_abort('File '//direction//'_cb.dat is missing in '// &
                             trim(t_step_dir)//'. Exiting.')
        end if
 
        ! Computing the cell-width distribution
        d_array(0:size_dim) = cb_array(0:size_dim) - cb_array(-1:size_dim - 1)
 
        ! Computing the cell-center locations
        cc_array(0:size_dim) = cb_array(-1:size_dim - 1) + d_array(0:size_dim)/2._wp
 
    end subroutine s_read_grid_data_direction
 
#ifdef MFC_MPI
    !> Helper subroutine to setup MPI data I/O parameters
    !!  @param data_size Local array size (output)
    !!  @param m_MOK, n_MOK, p_MOK MPI offset kinds for dimensions (output)
    !!  @param WP_MOK, MOK, str_MOK, NVARS_MOK Other MPI offset kinds (output)
    impure subroutine s_setup_mpi_io_params(data_size, m_MOK, n_MOK, p_MOK, WP_MOK, MOK, str_MOK, NVARS_MOK)
 
        integer, intent(out) :: data_size
        integer(KIND=MPI_OFFSET_KIND), intent(out) :: m_mok, n_mok, p_mok
        integer(KIND=MPI_OFFSET_KIND), intent(out) :: wp_mok, mok, str_mok, nvars_mok
 
        ! Initialize MPI data I/O
        if (ib) then
            call s_initialize_mpi_data(q_cons_vf, ib_markers)
        else
            call s_initialize_mpi_data(q_cons_vf)
        end if
 
        ! Size of local arrays
        data_size = (m + 1)*(n + 1)*(p + 1)
 
        ! Resize some integers so MPI can read even the biggest file
        m_mok = int(m_glb + 1, mpi_offset_kind)
        n_mok = int(n_glb + 1, mpi_offset_kind)
        p_mok = int(p_glb + 1, mpi_offset_kind)
        wp_mok = int(8._wp, mpi_offset_kind)
        mok = int(1._wp, mpi_offset_kind)
        str_mok = int(name_len, mpi_offset_kind)
        nvars_mok = int(sys_size, mpi_offset_kind)
 
    end subroutine s_setup_mpi_io_params
#endif
 
    !> Helper subroutine to read IB data files
    !!  @param file_loc_base Base file location for IB data
    !!  @param t_step Time step index
    impure subroutine s_read_ib_data_files(file_loc_base, t_step)
 
        character(len=*), intent(in) :: file_loc_base
        integer, intent(in), optional :: t_step
 
        character(LEN=len_trim(file_loc_base) + 20) :: file_loc
        logical :: file_exist
        integer :: ifile, ierr, data_size, var_mok
 
#ifdef MFC_MPI
        integer, dimension(MPI_STATUS_SIZE) :: status
        integer(KIND=MPI_OFFSET_KIND) :: disp
        integer :: m_mok, n_mok, p_mok, mok, wp_mok, save_index
 
#endif
        if (.not. ib) return
 
        if (parallel_io) then
            write (file_loc, '(A)') trim(file_loc_base)//'ib.dat'
        else
            write (file_loc, '(A)') trim(file_loc_base)//'/ib.dat'
        end if
        inquire (file=trim(file_loc), exist=file_exist)
 
        if (file_exist) then
            if (parallel_io) then
#ifdef MFC_MPI
                call mpi_file_open(mpi_comm_world, file_loc, mpi_mode_rdonly, mpi_info_int, ifile, ierr)
 
                m_mok = int(m_glb + 1, mpi_offset_kind)
                n_mok = int(n_glb + 1, mpi_offset_kind)
                p_mok = int(p_glb + 1, mpi_offset_kind)
                mok = int(1._wp, mpi_offset_kind)
                wp_mok = int(8._wp, mpi_offset_kind)
                save_index = t_step/t_step_save ! get the number of saves done to this point
 
                data_size = (m + 1)*(n + 1)*(p + 1)
                var_mok = int(sys_size + 1, mpi_offset_kind)
                if (t_step == 0) then
                    disp = 0
                else
                    disp = m_mok*max(mok, n_mok)*max(mok, p_mok)*wp_mok*(var_mok - 1 + int(save_index, mpi_offset_kind))
                end if
 
                call mpi_file_set_view(ifile, disp, mpi_integer, mpi_io_ib_data%view, &
                                       'native', mpi_info_int, ierr)
                call mpi_file_read(ifile, mpi_io_ib_data%var%sf, data_size, &
                                   mpi_integer, status, ierr)
 
                call mpi_file_close(ifile, ierr)
#endif
            else
                open (2, file=trim(file_loc), &
                      form='unformatted', &
                      action='read', &
                      status='old')
                read (2) ib_markers%sf(0:m, 0:n, 0:p)
                close (2)
            end if
        else
            call s_mpi_abort('File '//trim(file_loc)//' is missing. Exiting.')
        end if
 
    end subroutine s_read_ib_data_files
 
    !> Helper subroutine to allocate field arrays for given dimensionality
    !!  @param local_start_idx Starting index for allocation
    !!  @param end_x End index for x dimension
    !!  @param end_y End index for y dimension
    !!  @param end_z End index for z dimension
    impure subroutine s_allocate_field_arrays(local_start_idx, end_x, end_y, end_z)
 
        integer, intent(in) :: local_start_idx, end_x, end_y, end_z
        integer :: i
 
        do i = 1, sys_size
            allocate (q_cons_vf(i)%sf(local_start_idx:end_x, local_start_idx:end_y, local_start_idx:end_z))
            allocate (q_prim_vf(i)%sf(local_start_idx:end_x, local_start_idx:end_y, local_start_idx:end_z))
        end do
 
        if (ib) then
            allocate (ib_markers%sf(local_start_idx:end_x, local_start_idx:end_y, local_start_idx:end_z))
        end if
 
        if (chemistry) then
            allocate (q_t_sf%sf(local_start_idx:end_x, local_start_idx:end_y, local_start_idx:end_z))
        end if
 
    end subroutine s_allocate_field_arrays
 
    !>  This subroutine is called at each time-step that has to
        !!      be post-processed in order to read the raw data files
        !!      present in the corresponding time-step directory and to
        !!      populate the associated grid and conservative variables.
        !!  @param t_step Current time-step
    impure subroutine s_read_serial_data_files(t_step)
 
        integer, intent(in) :: t_step
 
        character(LEN=len_trim(case_dir) + 2*name_len) :: t_step_dir !<
            !! Location of the time-step directory associated with t_step
 
        character(LEN=len_trim(case_dir) + 3*name_len) :: file_loc !<
            !! Generic string used to store the location of a particular file
 
        character(len= &
                  int(floor(log10(real(sys_size, wp)))) + 1) :: file_num !<
            !! Used to store the variable position, in character form, of the
            !! currently manipulated conservative variable file
 
        character(LEN=len_trim(case_dir) + 2*name_len) :: t_step_ib_dir !<
        !! Location of the time-step directory associated with t_step
 
        logical :: dir_check !<
            !! Generic logical used to test the existence of a particular folder
 
        logical :: file_check  !<
            !! Generic logical used to test the existence of a particular file
 
        integer :: i !< Generic loop iterator
 
        ! Setting location of time-step folder based on current time-step
        write (t_step_dir, '(A,I0,A,I0)') '/p_all/p', proc_rank, '/', t_step
        t_step_dir = trim(case_dir)//trim(t_step_dir)
 
        ! Inquiring as to the existence of the time-step directory
        file_loc = trim(t_step_dir)//'/.'
        call my_inquire(file_loc, dir_check)
 
        ! If the time-step directory is missing, the post-process exits.
        if (dir_check .neqv. .true.) then
            call s_mpi_abort('Time-step folder '//trim(t_step_dir)// &
                             ' is missing. Exiting.')
        end if
 
        if (bc_io) then
            call s_read_serial_boundary_condition_files(t_step_dir, bc_type)
        else
            call s_assign_default_bc_type(bc_type)
        end if
 
        ! Reading the Grid Data Files using helper subroutine
        call s_read_grid_data_direction(t_step_dir, 'x', x_cb, dx, x_cc, m)
 
        if (n > 0) then
            call s_read_grid_data_direction(t_step_dir, 'y', y_cb, dy, y_cc, n)
 
            if (p > 0) then
                call s_read_grid_data_direction(t_step_dir, 'z', z_cb, dz, z_cc, p)
            end if
        end if
 
        ! Reading the Conservative Variables Data Files
        do i = 1, sys_size
 
            ! Checking whether the data file associated with the variable
            ! position of currently manipulated conservative variable exists
            write (file_num, '(I0)') i
            file_loc = trim(t_step_dir)//'/q_cons_vf'// &
                       trim(file_num)//'.dat'
            inquire (file=trim(file_loc), exist=file_check)
 
            ! Reading the data file if it exists, exiting otherwise
            if (file_check) then
                open (1, file=trim(file_loc), form='unformatted', &
                      status='old', action='read')
                read (1) q_cons_vf(i)%sf(0:m, 0:n, 0:p)
                close (1)
            else if (bubbles_lagrange .and. i == beta_idx) then
                ! beta (Lagrangian void fraction) is not written by pre_process
                ! for t_step_start; initialize to zero.
                q_cons_vf(i)%sf(0:m, 0:n, 0:p) = 0._wp
            else
                call s_mpi_abort('File q_cons_vf'//trim(file_num)// &
                                 '.dat is missing in '//trim(t_step_dir)// &
                                 '. Exiting.')
            end if
 
        end do
 
        ! Reading IB data using helper subroutine
        call s_read_ib_data_files(t_step_dir)
 
    end subroutine s_read_serial_data_files
 
    !>  This subroutine is called at each time-step that has to
        !!      be post-processed in order to parallel-read the raw data files
        !!      present in the corresponding time-step directory and to
        !!      populate the associated grid and conservative variables.
        !!  @param t_step Current time-step
    impure subroutine s_read_parallel_data_files(t_step)
 
        integer, intent(in) :: t_step
 
#ifdef MFC_MPI
 
        real(wp), allocatable, dimension(:) :: x_cb_glb, y_cb_glb, z_cb_glb
 
        integer :: ifile, ierr, data_size, filetype, stride
        integer, dimension(MPI_STATUS_SIZE) :: status
 
        integer(KIND=MPI_OFFSET_KIND) :: disp
        integer(KIND=MPI_OFFSET_KIND) :: m_mok, n_mok, p_mok
        integer(KIND=MPI_OFFSET_KIND) :: wp_mok, var_mok, str_mok
        integer(KIND=MPI_OFFSET_KIND) :: nvars_mok
        integer(KIND=MPI_OFFSET_KIND) :: mok
        integer(kind=MPI_OFFSET_KIND) :: offset
        real(wp) :: delx, dely, delz
 
        character(LEN=path_len + 2*name_len) :: file_loc
        logical :: file_exist
 
        character(len=10) :: t_step_string
 
        integer :: i
 
        allocate (x_cb_glb(-1:m_glb))
        allocate (y_cb_glb(-1:n_glb))
        allocate (z_cb_glb(-1:p_glb))
 
        if (down_sample) then
            stride = 3
        else
            stride = 1
        end if
 
        ! Read in cell boundary locations in x-direction
        file_loc = trim(case_dir)//'/restart_data'//trim(mpiiofs)//'x_cb.dat'
        inquire (file=trim(file_loc), exist=file_exist)
 
        if (file_exist) then
            data_size = m_glb + 2
            call mpi_file_open(mpi_comm_world, file_loc, mpi_mode_rdonly, mpi_info_int, ifile, ierr)
 
            call mpi_type_vector(data_size, 1, stride, mpi_p, filetype, ierr)
            call mpi_type_commit(filetype, ierr)
 
            offset = 0
            call mpi_file_set_view(ifile, offset, mpi_p, filetype, 'native', mpi_info_int, ierr)
 
            call mpi_file_read(ifile, x_cb_glb, data_size, mpi_p, status, ierr)
            call mpi_file_close(ifile, ierr)
        else
            call s_mpi_abort('File '//trim(file_loc)//' is missing. Exiting.')
        end if
 
        ! Assigning local cell boundary locations
        x_cb(-1:m) = x_cb_glb((start_idx(1) - 1):(start_idx(1) + m))
        ! Computing the cell width distribution
        dx(0:m) = x_cb(0:m) - x_cb(-1:m - 1)
        ! Computing the cell center location
        x_cc(0:m) = x_cb(-1:m - 1) + dx(0:m)/2._wp
 
        if (n > 0) then
            ! Read in cell boundary locations in y-direction
            file_loc = trim(case_dir)//'/restart_data'//trim(mpiiofs)//'y_cb.dat'
            inquire (file=trim(file_loc), exist=file_exist)
 
            if (file_exist) then
                data_size = n_glb + 2
                call mpi_file_open(mpi_comm_world, file_loc, mpi_mode_rdonly, mpi_info_int, ifile, ierr)
 
                call mpi_type_vector(data_size, 1, stride, mpi_p, filetype, ierr)
                call mpi_type_commit(filetype, ierr)
 
                offset = 0
                call mpi_file_set_view(ifile, offset, mpi_p, filetype, 'native', mpi_info_int, ierr)
 
                call mpi_file_read(ifile, y_cb_glb, data_size, mpi_p, status, ierr)
                call mpi_file_close(ifile, ierr)
            else
                call s_mpi_abort('File '//trim(file_loc)//' is missing. Exiting.')
            end if
 
            ! Assigning local cell boundary locations
            y_cb(-1:n) = y_cb_glb((start_idx(2) - 1):(start_idx(2) + n))
            ! Computing the cell width distribution
            dy(0:n) = y_cb(0:n) - y_cb(-1:n - 1)
            ! Computing the cell center location
            y_cc(0:n) = y_cb(-1:n - 1) + dy(0:n)/2._wp
 
            if (p > 0) then
                ! Read in cell boundary locations in z-direction
                file_loc = trim(case_dir)//'/restart_data'//trim(mpiiofs)//'z_cb.dat'
                inquire (file=trim(file_loc), exist=file_exist)
 
                if (file_exist) then
                    data_size = p_glb + 2
                    call mpi_file_open(mpi_comm_world, file_loc, mpi_mode_rdonly, mpi_info_int, ifile, ierr)
 
                    call mpi_type_vector(data_size, 1, stride, mpi_p, filetype, ierr)
                    call mpi_type_commit(filetype, ierr)
 
                    offset = 0
                    call mpi_file_set_view(ifile, offset, mpi_p, filetype, 'native', mpi_info_int, ierr)
 
                    call mpi_file_read(ifile, z_cb_glb, data_size, mpi_p, status, ierr)
                    call mpi_file_close(ifile, ierr)
                else
                    call s_mpi_abort('File '//trim(file_loc)//' is missing. Exiting.')
                end if
 
                ! Assigning local cell boundary locations
                z_cb(-1:p) = z_cb_glb((start_idx(3) - 1):(start_idx(3) + p))
                ! Computing the cell width distribution
                dz(0:p) = z_cb(0:p) - z_cb(-1:p - 1)
                ! Computing the cell center location
                z_cc(0:p) = z_cb(-1:p - 1) + dz(0:p)/2._wp
            end if
        end if
 
        call s_read_parallel_conservative_data(t_step, m_mok, n_mok, p_mok, wp_mok, mok, str_mok, nvars_mok)
 
        deallocate (x_cb_glb, y_cb_glb, z_cb_glb)
 
        if (bc_io) then
            call s_read_parallel_boundary_condition_files(bc_type)
        else
            call s_assign_default_bc_type(bc_type)
        end if
 
#endif
 
    end subroutine s_read_parallel_data_files
 
#ifdef MFC_MPI
    !> Helper subroutine to read parallel conservative variable data
    !!  @param t_step Current time-step
    !!  @param m_MOK, n_MOK, p_MOK MPI offset kinds for dimensions
    !!  @param WP_MOK, MOK, str_MOK, NVARS_MOK Other MPI offset kinds
    impure subroutine s_read_parallel_conservative_data(t_step, m_MOK, n_MOK, p_MOK, WP_MOK, MOK, str_MOK, NVARS_MOK)
 
        integer, intent(in) :: t_step
        integer(KIND=MPI_OFFSET_KIND), intent(inout) :: m_mok, n_mok, p_mok
        integer(KIND=MPI_OFFSET_KIND), intent(inout) :: wp_mok, mok, str_mok, nvars_mok
 
        integer :: ifile, ierr, data_size
        integer, dimension(MPI_STATUS_SIZE) :: status
        integer(KIND=MPI_OFFSET_KIND) :: disp, var_mok
        character(LEN=path_len + 2*name_len) :: file_loc
        logical :: file_exist
        character(len=10) :: t_step_string
        integer :: i
 
        if (file_per_process) then
            call s_int_to_str(t_step, t_step_string)
            ! Open the file to read conservative variables
            write (file_loc, '(I0,A1,I7.7,A)') t_step, '_', proc_rank, '.dat'
            file_loc = trim(case_dir)//'/restart_data/lustre_'//trim(t_step_string)//trim(mpiiofs)//trim(file_loc)
            inquire (file=trim(file_loc), exist=file_exist)
 
            if (file_exist) then
                call mpi_file_open(mpi_comm_self, file_loc, mpi_mode_rdonly, mpi_info_int, ifile, ierr)
 
                if (down_sample) then
                    call s_initialize_mpi_data_ds(q_cons_temp)
                else
                    ! Initialize MPI data I/O
                    if (ib) then
                        call s_initialize_mpi_data(q_cons_vf, ib_markers)
                    else
                        call s_initialize_mpi_data(q_cons_vf)
                    end if
                end if
 
                if (down_sample) then
                    ! Size of local arrays
                    data_size = (m + 3)*(n + 3)*(p + 3)
                else
                    ! Size of local arrays
                    data_size = (m + 1)*(n + 1)*(p + 1)
                end if
 
                ! Resize some integers so MPI can read even the biggest file
                m_mok = int(m_glb + 1, mpi_offset_kind)
                n_mok = int(n_glb + 1, mpi_offset_kind)
                p_mok = int(p_glb + 1, mpi_offset_kind)
                wp_mok = int(8._wp, mpi_offset_kind)
                mok = int(1._wp, mpi_offset_kind)
                str_mok = int(name_len, mpi_offset_kind)
                nvars_mok = int(sys_size, mpi_offset_kind)
 
                ! Read the data for each variable
                if (bubbles_euler .or. elasticity .or. mhd) then
                    do i = 1, sys_size
                        var_mok = int(i, mpi_offset_kind)
                        call mpi_file_read_all(ifile, mpi_io_data%var(i)%sf, data_size*mpi_io_type, &
                                               mpi_io_p, status, ierr)
                    end do
                else
                    do i = 1, sys_size
                        var_mok = int(i, mpi_offset_kind)
                        call mpi_file_read_all(ifile, mpi_io_data%var(i)%sf, data_size*mpi_io_type, &
                                               mpi_io_p, status, ierr)
                    end do
                end if
 
                call s_mpi_barrier()
                call mpi_file_close(ifile, ierr)
 
                if (down_sample) then
                    do i = 1, sys_size
                        q_cons_vf(i)%sf(0:m, 0:n, 0:p) = q_cons_temp(i)%sf(0:m, 0:n, 0:p)
                    end do
                end if
 
                call s_read_ib_data_files(trim(case_dir)//'/restart_data'//trim(mpiiofs), t_step)
            else
                call s_mpi_abort('File '//trim(file_loc)//' is missing. Exiting.')
            end if
        else
            ! Open the file to read conservative variables
            write (file_loc, '(I0,A)') t_step, '.dat'
            file_loc = trim(case_dir)//'/restart_data'//trim(mpiiofs)//trim(file_loc)
            inquire (file=trim(file_loc), exist=file_exist)
 
            if (file_exist) then
                call mpi_file_open(mpi_comm_world, file_loc, mpi_mode_rdonly, mpi_info_int, ifile, ierr)
 
                call s_setup_mpi_io_params(data_size, m_mok, n_mok, p_mok, wp_mok, mok, str_mok, nvars_mok)
 
                ! Read the data for each variable
                do i = 1, sys_size
                    var_mok = int(i, mpi_offset_kind)
 
                    ! Initial displacement to skip at beginning of file
                    disp = m_mok*max(mok, n_mok)*max(mok, p_mok)*wp_mok*(var_mok - 1)
 
                    call mpi_file_set_view(ifile, disp, mpi_p, mpi_io_data%view(i), &
                                           'native', mpi_info_int, ierr)
                    call mpi_file_read_all(ifile, mpi_io_data%var(i)%sf, data_size*mpi_io_type, &
                                           mpi_io_p, status, ierr)
                end do
 
                call s_mpi_barrier()
                call mpi_file_close(ifile, ierr)
 
                call s_read_ib_data_files(trim(case_dir)//'/restart_data'//trim(mpiiofs), t_step)
            else
                call s_mpi_abort('File '//trim(file_loc)//' is missing. Exiting.')
            end if
        end if
    end subroutine s_read_parallel_conservative_data
#endif
 
    !>  Computation of parameters, allocation procedures, and/or
        !!      any other tasks needed to properly setup the module
    impure subroutine s_initialize_data_input_module
 
        integer :: i !< Generic loop iterator
 
        ! Allocating the parts of the conservative and primitive variables
        ! that do not require the direct knowledge of the dimensionality of
        ! the simulation
        allocate (q_cons_vf(1:sys_size))
        allocate (q_prim_vf(1:sys_size))
        allocate (q_cons_temp(1:sys_size))
 
        ! Allocating the parts of the conservative and primitive variables
        ! that do require the direct knowledge of the dimensionality of
        ! the simulation using helper subroutine
 
        ! Simulation is at least 2D
        if (n > 0) then
            ! Simulation is 3D
            if (p > 0) then
                call s_allocate_field_arrays(-buff_size, m + buff_size, n + buff_size, p + buff_size)
                if (down_sample) then
                    do i = 1, sys_size
                        allocate (q_cons_temp(i)%sf(-1:m + 1, -1:n + 1, -1:p + 1))
                    end do
                end if
            else
                ! Simulation is 2D
                call s_allocate_field_arrays(-buff_size, m + buff_size, n + buff_size, 0)
            end if
        else
            ! Simulation is 1D
            call s_allocate_field_arrays(-buff_size, m + buff_size, 0, 0)
        end if
 
        ! Allocating arrays to store the bc types
        allocate (bc_type(1:num_dims, 1:2))
 
        allocate (bc_type(1, 1)%sf(0:0, 0:n, 0:p))
        allocate (bc_type(1, 2)%sf(0:0, 0:n, 0:p))
        if (n > 0) then
            allocate (bc_type(2, 1)%sf(-buff_size:m + buff_size, 0:0, 0:p))
            allocate (bc_type(2, 2)%sf(-buff_size:m + buff_size, 0:0, 0:p))
            if (p > 0) then
                allocate (bc_type(3, 1)%sf(-buff_size:m + buff_size, -buff_size:n + buff_size, 0:0))
                allocate (bc_type(3, 2)%sf(-buff_size:m + buff_size, -buff_size:n + buff_size, 0:0))
            end if
        end if
 
        if (parallel_io .neqv. .true.) then
            s_read_data_files => s_read_serial_data_files
        else
            s_read_data_files => s_read_parallel_data_files
        end if
 
    end subroutine s_initialize_data_input_module
 
    !> Deallocation procedures for the module
    impure subroutine s_finalize_data_input_module
 
        integer :: i !< Generic loop iterator
 
        ! Deallocating the conservative and primitive variables
        do i = 1, sys_size
            deallocate (q_cons_vf(i)%sf)
            deallocate (q_prim_vf(i)%sf)
            if (down_sample) then
                deallocate (q_cons_temp(i)%sf)
            end if
        end do
 
        deallocate (q_cons_vf)
        deallocate (q_prim_vf)
        deallocate (q_cons_temp)
 
        if (ib) then
            deallocate (ib_markers%sf)
        end if
 
        if (chemistry) then
            deallocate (q_t_sf%sf)
        end if
 
        deallocate (bc_type(1, 1)%sf, bc_type(1, 2)%sf)
        if (n > 0) then
            deallocate (bc_type(2, 1)%sf, bc_type(2, 2)%sf)
            if (p > 0) then
                deallocate (bc_type(3, 1)%sf, bc_type(3, 2)%sf)
            end if
        end if
 
        deallocate (bc_type)
 
        s_read_data_files => null()
 
    end subroutine s_finalize_data_input_module
 
end module m_data_input