Finds bonded clusters of atoms in the structure with periodic boundary conditions.
    If there are atoms that are not bonded to anything, returns [0,1,0].(For faster computation time in FindDimension())

    Args:
        struct (Structure): Input structure
        connected_list: Must be made from the same structure with FindConnected() function.
            An array of shape (number of bonded pairs, 2); each row of is of the form [atomi, atomj].

    Returns:
        max_cluster: the size of the largest cluster in the crystal structure
        min_cluster: the size of the smallest cluster in the crystal structure
        clusters: list of bonded clusters found here, clusters are formatted as sets of indices of atoms

def find_clusters(struct, connected_list):
    """
    Finds bonded clusters of atoms in the structure with periodic boundary conditions.
    If there are atoms that are not bonded to anything, returns [0,1,0].(For faster computation time in FindDimension())

    Args:
        struct (Structure): Input structure
        connected_list: Must be made from the same structure with FindConnected() function.
            An array of shape (number of bonded pairs, 2); each row of is of the form [atomi, atomj].

    Returns:
        max_cluster: the size of the largest cluster in the crystal structure
        min_cluster: the size of the smallest cluster in the crystal structure
        clusters: list of bonded clusters found here, clusters are formatted as sets of indices of atoms
    """
    n_atoms = len(struct.species)
    if len(np.unique(connected_list)) != n_atoms:
        return [0, 1, 0]
    if n_atoms == 0:
        return [0, 0, 0]
    cluster_sizes = []
    clusters = []
    for atom in range(n_atoms):
        connected_inds = np.where(connected_list == atom)[0]
        atom_cluster = np.unique(connected_list[connected_inds])
        atom_cluster = set(atom_cluster)
        if len(clusters) == 0:
            new_clusters = [atom_cluster]
            new_cluster_sizes = [len(atom_cluster)]
        else:
            clusters_w_atom = [atom_cluster]
            clusters_noatom = []
            clusters_noatom_sizes = []
            for cluster in clusters:
                if len(cluster.intersection(atom_cluster)) > 0:
                    clusters_w_atom.append(cluster)
                else:
                    clusters_noatom.append(cluster)
                    clusters_noatom_sizes.append(len(cluster))
            if len(clusters_w_atom) > 1:
                clusters_w_atom = [set.union(*clusters_w_atom)]
            new_clusters = clusters_noatom + clusters_w_atom
            new_cluster_sizes = clusters_noatom_sizes + [len(clusters_w_atom[0])]
        clusters = list(new_clusters)
        cluster_sizes = list(new_cluster_sizes)
        if n_atoms in cluster_sizes:
            break
    max_cluster = max(cluster_sizes)
    min_cluster = min(cluster_sizes)
    return [max_cluster, min_cluster, clusters]