module box
    !This module contains information on the properties of the current box.
    use parameters
    use functions
    implicit none 

    real(kind=dp) :: box_bd(6) !Global box boundaries
    character(len=3) :: box_bc !Box boundary conditions (periodic or shrinkwrapped)
    logical :: bound_called 
    !The subbox variables contain values for each subbox, being the boxes read in through some
    !command. Currently only mode_merge will require sub_boxes, for mode_create it will always 
    !allocate to only 1 sub_box
    integer :: sub_box_num = 0
    real(kind=dp), allocatable :: sub_box_ori(:,:,:)!Orientations for each of the subboxes
    real(kind=dp), allocatable :: sub_box_bd(:,:)!Boundaries for each of the sub_boxes

    !Below are some simulation parameters which may be adjusted if reading in restart files
    integer :: timestep=0
    real(kind=dp) :: total_time=0.0_dp


    public 
    contains 

    subroutine box_init
        !Initialize some box functions
        box_bd(:) = 0.0_dp
        box_bc = 'ppp'
        bound_called=.false.
    end subroutine box_init

    subroutine alloc_sub_box(n)
        !Allocate the sub_box variables

        integer, intent(in) :: n

        integer :: i

        allocate(sub_box_ori(3,3,n), sub_box_bd(6,n))
        do i = 1, n
            sub_box_ori(:,:,i) = identity_mat(3)
            sub_box_bd(:,i) = 0.0_dp
        end do
    end subroutine alloc_sub_box

    subroutine grow_sub_box(n)
        !Grows sub box arrays, this is only called when a new file is read in
        integer, intent(in) :: n

        integer, allocatable :: temp_array_bd(:,:,:), temp_file(:)
        real(kind=dp), allocatable :: temp_ori(:,:,:), temp_bd(:,:)
        !Allocate temporary arrays
        allocate(temp_ori(3,3,sub_box_num+n),temp_bd(6,sub_box_num+n), &
                 temp_array_bd(2,2,sub_box_num+n), temp_file(sub_box_num+n))

        !Move allocation for all sub_box_arrays
        temp_ori(:,:,1:sub_box_num) = sub_box_ori
        temp_ori(:,:,sub_box_num+1:) = 0.0_dp
        call move_alloc(temp_ori, sub_box_ori)

        temp_bd(:, 1:sub_box_num) = sub_box_bd
        temp_bd(:, sub_box_num+1:) = 0.0_dp
        call move_alloc(temp_bd, sub_box_bd)

        return
    end subroutine grow_sub_box

    subroutine grow_box(temp_box_bd)
        !This function takes in a temporary box boundary and adjusts the overall box boundaries
        !to include  it

        real(kind=dp), dimension(6), intent(in) :: temp_box_bd

        integer :: i

        if(all(abs(box_bd) < lim_zero)) then 
            box_bd = temp_box_bd
        else
            do i = 1, 3
                if(temp_box_bd(2*i-1) < box_bd(2*i-1)) box_bd(2*i-1) = temp_box_bd(2*i-1)
                if(temp_box_bd(2*i) > box_bd(2*i)) box_bd(2*i) = temp_box_bd(2*i)
            end do
        end if

        return
    end subroutine grow_box


    subroutine in_sub_box(r, which_sub_box)
        !This returns which sub_box a point is in. It returns the first sub_box with boundaries which 
        !contain the point.
        real(kind=dp), dimension(3), intent(in) :: r
        integer, intent(out) :: which_sub_box

        integer :: i 

        do i = 1, sub_box_num
            if( in_block_bd(r, sub_box_bd(:,i))) then 
                which_sub_box = i
                exit
            end if
        end do 
    return 
    end subroutine in_sub_box

    subroutine reset_box
        !This subroutine just resets the box boundary and position
        box_bc = "ppp"
        box_bd(:) = 0
    end subroutine reset_box

    pure function box_volume()
        real(kind = dp) :: box_volume

        box_volume = (box_bd(2) - box_bd(1)) * (box_bd(4) - box_bd(3))*(box_bd(6) - box_bd(5))
        return
    end function
end module box