Small world (Watts-Strogatz) initializer

[MartinuzziFrancesco]

Description

Given that the standard initialization of the reservoir follows Erdős–Rényi it makes sense to try different topologies. In Kawai2019 the authors provide a very detailed analysis of the small world topology, and it could be used as base for the implementation

Possible implementation

Of course an implementation is already available at Graphs.jl but that seems a big dependency to call for just this one function.

using Random

# watts_strogatz or small_world?
function small_world(rng::AbstractRNG, ::Type{T}, dims::Integer...;
        num_neighbors::Integer, rewiring_probability::Float64, radius::AbstractFloat=1.0) where T
    num_nodes = if length(dims) == 1
        dims[1]
    elseif length(dims) == 2
        if dims[1] != dims[2]
            error("The provided dimensions must form a square matrix.")
        else
            dims[1]
        end
    else
        error("dims must be one integer or two equal integers for a square matrix.")
    end
    if num_neighbors % 2 != 0 || num_neighbors >= num_nodes
        error("num_neighbors must be an even number and less than the number of nodes.")
    end
    reservoir_matrix = zeros(T, num_nodes, num_nodes)
    half_neighbors = div(num_neighbors, 2)
    for current_node in 1:num_nodes
        for offset in 1:half_neighbors
            neighbor_node = mod1(current_node + offset, num_nodes)
            reservoir_matrix[current_node, neighbor_node] = one(T)
            reservoir_matrix[neighbor_node, current_node] = one(T)
        end
    end
    for current_node in 1:num_nodes
        for offset in 1:half_neighbors
            neighbor_node = mod1(current_node + offset, num_nodes)
            if rand(rng) < rewiring_probability
                potential_new_neighbors = [candidate for candidate in 1:num_nodes if candidate != current_node && reservoir_matrix[current_node, candidate] == 0]
                if !isempty(potential_new_neighbors)
                    new_neighbor_node = rand(rng, potential_new_neighbors)
                    reservoir_matrix[current_node, neighbor_node] = zero(T)
                    reservoir_matrix[neighbor_node, current_node] = zero(T)
                    reservoir_matrix[current_node, new_neighbor_node] = one(T)
                    reservoir_matrix[new_neighbor_node, current_node] = one(T)
                end
            end
        end
    end
    scale_radius!(reservoir_matrix, radius)
    return reservoir_matrix
end