remove_vertices Module Subroutine

  module subroutine remove_vertices(this, indices)
    !! Remove vertices from the graph.
    implicit none

    ! Arguments
    class(graph_type), intent(inout) :: this
    !! Parent. Instance of the graph structure.
    integer, dimension(:), intent(in) :: indices
    !! Indices of the vertices to be removed.

    ! Local variables
    integer :: i, j, k
    !! Loop indices.
    integer, dimension(size(indices, dim=1)) :: vertex_indices
    !! Indices of the vertices to be removed.
    integer, dimension(:), allocatable :: edge_indices
    !! Indices of the edges to be removed.

    ! Find all edges connected to vertices being removed
    allocate(edge_indices(0))
    do i = 1, size(indices)
       do j = 1, this%num_edges
          if(any(this%edge(j)%index .eq. indices(i)) .or. &
              any(this%edge(j)%index .eq. -indices(i)))then
             edge_indices = [edge_indices, j]
          end if
       end do
    end do
    
    ! Remove connected edges first
    if(size(edge_indices) .gt. 0)then
       call this%remove_edges(edge_indices, update_adjacency=.false.)
    end if

    ! Remove vertices in descending order to maintain index validity
    vertex_indices = indices
    do i = 1, size(indices)
       k = maxval(vertex_indices)
       this%vertex = [this%vertex(1:k-1), this%vertex(k+1:this%num_vertices)]
       this%num_vertices = this%num_vertices - 1
       vertex_indices(maxloc(vertex_indices, 1)) = 0
       
       ! Update edge indices to reflect removed vertex
       do j = 1, this%num_edges
          if(this%edge(j)%index(1) .gt. k)then
             this%edge(j)%index(1) = this%edge(j)%index(1) - 1
          elseif(this%edge(j)%index(1) .lt. -k)then
             this%edge(j)%index(1) = this%edge(j)%index(1) + 1
          end if
          if(this%edge(j)%index(2) .gt. k)then
             this%edge(j)%index(2) = this%edge(j)%index(2) - 1
          elseif(this%edge(j)%index(2) .lt. -k)then
             this%edge(j)%index(2) = this%edge(j)%index(2) + 1
          end if
       end do
    end do
    call this%generate_adjacency()
  end subroutine remove_vertices