module opt_redef_box

    use box
    use elements
    use subroutines
    implicit none
    
    character(len=3) :: redef_dim, new_bc
    real(kind=dp) :: new_bd(6)
    public 
    contains

    subroutine redef_box(arg_pos)
        !This is the main calling function for opt_redef_box
        integer, intent(inout) :: arg_pos
        integer :: i, inod, ibasis, iatom, delete_list(atom_num), delete_num, type_interp(max_basisnum*max_esize**3)
        real(kind=dp):: r_interp(3, max_basisnum*max_esize**3)
        logical :: node_out(8)

        !First parse the argument
        call parse_redef_box(arg_pos)

        print *, '------------------------------------------------------------'
        print *, 'Option redef_box'
        print *, '------------------------------------------------------------'

        !Now first filter atoms that don't fit in the new box bounds and delete them
        delete_num = 0
        do i = 1, atom_num 
            if(.not.in_block_bd(r_atom(:,i),new_bd)) then 
                delete_num = delete_num + 1
                delete_list(delete_num) = i
            end if
        end do
        call delete_atoms(delete_num, delete_list(1:delete_num))

        !Now loop over all elements
        delete_num = 0
        delete_list(:) = 0
        do i = 1, ele_num
            !Determine if all nodes are within the new boundaries
            node_out(:) = .false.
            do inod = 1, ng_node(lat_ele(i))
                do ibasis = 1, basisnum(lat_ele(i))
                    if(.not.in_block_bd(r_node(:,ibasis,inod,i), new_bd)) then 
                        node_out(inod) = .true.
                    end if
                end do
            end do

            !If all nodes are out just add the element to the delete list
            if(all(node_out)) then 
                delete_num = delete_num +1
                delete_list(delete_num) = i

            !If any nodes are out we add the element to the delete list, but then loop over the interpoalted atoms to figure out
            !which ones fit inside the boundary to keep the box rectangular 
            else if (any(node_out)) then 
                delete_num = delete_num +1
                delete_list(delete_num) = i

                call interpolate_atoms(type_ele(i), size_ele(i), lat_ele(i), r_node(:,:,:,i), type_interp, r_interp)

                !loop over all interpolated atoms and add them to the system
                do iatom = 1, basisnum(lat_ele(i))*size_ele(i)**3
                    if(in_block_bd(r_interp(:,iatom), new_bd)) then 
                        call add_atom(0,type_interp(iatom), sbox_ele(i), r_interp(:,iatom))
                    end if
                end do 
            end if
        end do

        call delete_elements(delete_num, delete_list(1:delete_num))

        print *, "Old box_bd: ", box_bd, " is redefined to new box boundaries: ", new_bd
        box_bd=new_bd
        box_bc = new_bc

    end subroutine redef_box

        subroutine parse_redef_box(arg_pos)
        !Parse the command
        integer, intent(inout) :: arg_pos

        integer :: i, j, arglen
        character(len=100) textholder
        

        !First read in the dimensions that we are redefining
        redef_dim = ''
        arg_pos=arg_pos+1
        call get_command_argument(arg_pos, redef_dim, arglen)
        select case(trim(adjustl(redef_dim))) 
        case('x','y','z','xy','yx','xz','zx','yz','zy','xyz','yxz','xzy','zyx','zxy','yzx')
            continue
        case default
            print *, "Incorrect redef_dim ", redef_dim, "please select any permuation of x, y, z, xy, yz, xz, xyz"
            stop 3
        end select

        !Now read in the new dimensions
        new_bd = box_bd
        new_bc = box_bc
        do i = 1, 3
            select case(trim(adjustl(redef_dim(i:i))))
            case('x')
                j = 1
            case('y')
                j = 2
            case('z')
                j = 3
            case default
                exit
            end select 

            arg_pos=arg_pos +1
            call get_command_argument(arg_pos, textholder, arglen)
            if(arglen == 0) stop "Missing a box dimension in opt_redef_box"
            call parse_pos(j, textholder,new_bd(2*j-1))

            arg_pos=arg_pos +1
            call get_command_argument(arg_pos, textholder, arglen)
            if(arglen == 0) stop "Missing a box dimension in opt_redef_box"
            call parse_pos(j, textholder,new_bd(2*j))
            new_bc(j:j) = 's'
        end do 

        arg_pos = arg_pos + 1
    end subroutine parse_redef_box

end module opt_redef_box