Skip to content

Commit

Permalink
Small world network implementation. WIP commit prior to FLAMEGPU_INIT…
Browse files Browse the repository at this point in the history
…_FUNC refactor
  • Loading branch information
ptheywood committed Dec 10, 2024
1 parent bfee518 commit 53b1fdb
Show file tree
Hide file tree
Showing 13 changed files with 722 additions and 86 deletions.
14 changes: 7 additions & 7 deletions data/inputs/sample.csv
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
rng_seed,param_id,duration,n_total,n_seed_infection,population_0_9,population_10_19,population_20_29,population_30_39,population_40_49,population_50_59,population_60_69,population_70_79,population_80,household_size_1,household_size_2,household_size_3,household_size_4,household_size_5,household_size_6,p_interaction_susceptible_to_exposed,mean_time_to_infected,sd_time_to_infected,mean_time_to_recovered,sd_time_to_recovered,mean_time_to_susceptible,sd_time_to_susceptible,relative_susceptibility_0_9,relative_susceptibility_10_19,relative_susceptibility_20_29,relative_susceptibility_30_39,relative_susceptibility_40_49,relative_susceptibility_50_59,relative_susceptibility_60_69,relative_susceptibility_70_79,relative_susceptibility_80,child_network_adults,elderly_network_adults,relative_transmission_household,relative_transmission_occupation,relative_transmission_random,daily_fraction_work,mean_random_interactions_0_19,sd_random_interactions_0_19,mean_random_interactions_20_69,sd_random_interactions_20_69,mean_random_interactions_70plus,sd_random_interactions_70plus
12,0,60,32,4,331,370,374,371,361,415,377,270,142,442,492,208,165,51,27,0.10,3.5,2.0,7.0,2.0,28,10.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,0.2,0.2,2.0,0.5,0.5,0.001,2,2,4,4,3,3
12,1,365,32,4,331,370,374,371,361,415,377,270,142,442,492,208,165,51,27,0.10,3.5,2.0,7.0,2.0,28,10.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,0.2,0.2,2.0,0.5,0.5,0.001,2,2,4,4,3,3
12,2,365,1024,8,331,370,374,371,361,415,377,270,142,442,492,208,165,51,27,0.10,3.5,2.0,7.0,2.0,28,10.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,0.2,0.2,2.0,0.5,0.5,0.001,2,2,4,4,3,3
12,3,365,4096,8,331,370,374,371,361,415,377,270,142,442,492,208,165,51,27,0.10,3.5,2.0,7.0,2.0,28,10.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,0.2,0.2,2.0,0.5,0.5,0.001,2,2,4,4,3,3
12,4,365,4096,8,331,370,374,371,361,415,377,270,142,442,492,208,165,51,27,0.10,3.5,2.0,7.0,2.0,28,10.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,0.2,0.2,2.0,0.5,0.5,0.001,2,2,4,4,3,3
12,5,365,65536,8,331,370,374,371,361,415,377,270,142,442,492,208,165,51,27,0.10,3.5,2.0,7.0,2.0,28,10.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,0.2,0.2,2.0,0.5,0.5,0.001,2,2,4,4,3,3
rng_seed,param_id,duration,n_total,n_seed_infection,population_0_9,population_10_19,population_20_29,population_30_39,population_40_49,population_50_59,population_60_69,population_70_79,population_80,household_size_1,household_size_2,household_size_3,household_size_4,household_size_5,household_size_6,p_interaction_susceptible_to_exposed,mean_time_to_infected,sd_time_to_infected,mean_time_to_recovered,sd_time_to_recovered,mean_time_to_susceptible,sd_time_to_susceptible,relative_susceptibility_0_9,relative_susceptibility_10_19,relative_susceptibility_20_29,relative_susceptibility_30_39,relative_susceptibility_40_49,relative_susceptibility_50_59,relative_susceptibility_60_69,relative_susceptibility_70_79,relative_susceptibility_80,child_network_adults,elderly_network_adults,relative_transmission_household,relative_transmission_occupation,relative_transmission_random,mean_work_interactions_child,mean_work_interactions_adult,mean_work_interactions_elderly,daily_fraction_work,work_network_rewire,mean_random_interactions_0_19,sd_random_interactions_0_19,mean_random_interactions_20_69,sd_random_interactions_20_69,mean_random_interactions_70plus,sd_random_interactions_70plus
12,0,60,32,4,331,370,374,371,361,415,377,270,142,442,492,208,165,51,27,0.10,3.5,2.0,7.0,2.0,28,10.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,0.2,0.2,2.0,0.5,0.5,10,7,3,0.5,0.1,2,2,4,4,3,3
12,1,365,32,4,331,370,374,371,361,415,377,270,142,442,492,208,165,51,27,0.10,3.5,2.0,7.0,2.0,28,10.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,0.2,0.2,2.0,0.5,0.5,10,7,3,0.5,0.1,2,2,4,4,3,3
12,2,365,1024,8,331,370,374,371,361,415,377,270,142,442,492,208,165,51,27,0.10,3.5,2.0,7.0,2.0,28,10.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,0.2,0.2,2.0,0.5,0.5,10,7,3,0.5,0.1,2,2,4,4,3,3
12,3,365,4096,8,331,370,374,371,361,415,377,270,142,442,492,208,165,51,27,0.10,3.5,2.0,7.0,2.0,28,10.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,0.2,0.2,2.0,0.5,0.5,10,7,3,0.5,0.1,2,2,4,4,3,3
12,4,365,4096,8,331,370,374,371,361,415,377,270,142,442,492,208,165,51,27,0.10,3.5,2.0,7.0,2.0,28,10.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,0.2,0.2,2.0,0.5,0.5,10,7,3,0.5,0.1,2,2,4,4,3,3
12,5,365,65536,8,331,370,374,371,361,415,377,270,142,442,492,208,165,51,27,0.10,3.5,2.0,7.0,2.0,28,10.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,0.2,0.2,2.0,0.5,0.5,10,7,3,0.5,0.1,2,2,4,4,3,3
2 changes: 2 additions & 0 deletions src/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -21,6 +21,7 @@ SET(LIBRARY_SRC
${CMAKE_CURRENT_SOURCE_DIR}/exateppabm/disease/SEIR.h
${CMAKE_CURRENT_SOURCE_DIR}/exateppabm/exatepp_abm.h
${CMAKE_CURRENT_SOURCE_DIR}/exateppabm/input.h
${CMAKE_CURRENT_SOURCE_DIR}/exateppabm/network.h
${CMAKE_CURRENT_SOURCE_DIR}/exateppabm/output.h
${CMAKE_CURRENT_SOURCE_DIR}/exateppabm/output/OutputFile.h
${CMAKE_CURRENT_SOURCE_DIR}/exateppabm/output/PerformanceFile.h
Expand All @@ -36,6 +37,7 @@ SET(LIBRARY_SRC
${CMAKE_CURRENT_SOURCE_DIR}/exateppabm/disease/SEIR.cu
${CMAKE_CURRENT_SOURCE_DIR}/exateppabm/exatepp_abm.cu
${CMAKE_CURRENT_SOURCE_DIR}/exateppabm/input.cu
${CMAKE_CURRENT_SOURCE_DIR}/exateppabm/network.cu
${CMAKE_CURRENT_SOURCE_DIR}/exateppabm/output.cu
${CMAKE_CURRENT_SOURCE_DIR}/exateppabm/output/PerformanceFile.cpp
${CMAKE_CURRENT_SOURCE_DIR}/exateppabm/output/TimeSeriesFile.cpp
Expand Down
1 change: 0 additions & 1 deletion src/exateppabm/demographics.h
Original file line number Diff line number Diff line change
Expand Up @@ -45,7 +45,6 @@ void define(flamegpu::ModelDescription& model, const exateppabm::input::config&
* Get an array containing one of each age demographic enum
*
* This is a workaround for the lack of reflection in c++17, used to simplify code elsewhere
* @todo constexpr? return refernce?
* @return std::array which is the inverse of the Age enum.
*/
std::array<demographics::Age, demographics::AGE_COUNT> getAllAgeDemographics();
Expand Down
22 changes: 22 additions & 0 deletions src/exateppabm/input.cu
Original file line number Diff line number Diff line change
Expand Up @@ -15,6 +15,8 @@ namespace exateppabm {
namespace input {
namespace {

constexpr char header[] = "rng_seed,param_id,duration,n_total,n_seed_infection,population_0_9,population_10_19,population_20_29,population_30_39,population_40_49,population_50_59,population_60_69,population_70_79,population_80,household_size_1,household_size_2,household_size_3,household_size_4,household_size_5,household_size_6,p_interaction_susceptible_to_exposed,mean_time_to_infected,sd_time_to_infected,mean_time_to_recovered,sd_time_to_recovered,mean_time_to_susceptible,sd_time_to_susceptible,relative_susceptibility_0_9,relative_susceptibility_10_19,relative_susceptibility_20_29,relative_susceptibility_30_39,relative_susceptibility_40_49,relative_susceptibility_50_59,relative_susceptibility_60_69,relative_susceptibility_70_79,relative_susceptibility_80,child_network_adults,elderly_network_adults,relative_transmission_household,relative_transmission_occupation,relative_transmission_random,mean_work_interactions_child,mean_work_interactions_adult,mean_work_interactions_elderly,daily_fraction_work,work_network_rewire,mean_random_interactions_0_19,sd_random_interactions_0_19,mean_random_interactions_20_69,sd_random_interactions_20_69,mean_random_interactions_70plus,sd_random_interactions_70plus"; //NOLINT

/**
* Get the next value from a csv row as a specified type
*
Expand Down Expand Up @@ -78,6 +80,10 @@ std::shared_ptr<exateppabm::input::config> read(const std::filesystem::path p, c
if (!std::getline(fs, line)) {
throw std::runtime_error("failed to read the header line @todo nicer error message");
}
// Warn if the header is not the expected value
if (strcmp(line.c_str(), header) != 0) {
fmt::print("Warning: {} header does not match the expected value. Errors may occur, or parameters may be used incorrectly.\nExpected header:\n{}\n", p.c_str(), header);
}

// Discard rows until the line number is the target line number
for (int currentLine = 1; currentLine < lineNumber; currentLine++) {
Expand Down Expand Up @@ -214,9 +220,21 @@ std::shared_ptr<exateppabm::input::config> read(const std::filesystem::path p, c
if (!valueFromCSVLine(line, c->relative_transmission_random)) {
throw std::runtime_error("bad value for relative_transmission_random during csv parsing @todo\n");
}
if (!valueFromCSVLine(line, c->mean_work_interactions_child)) {
throw std::runtime_error("bad value for mean_work_interactions_child during csv parsing @todo\n");
}
if (!valueFromCSVLine(line, c->mean_work_interactions_adult)) {
throw std::runtime_error("bad value for mean_work_interactions_adult during csv parsing @todo\n");
}
if (!valueFromCSVLine(line, c->mean_work_interactions_elderly)) {
throw std::runtime_error("bad value for mean_work_interactions_elderly during csv parsing @todo\n");
}
if (!valueFromCSVLine(line, c->daily_fraction_work)) {
throw std::runtime_error("bad value for daily_fraction_work during csv parsing @todo\n");
}
if (!valueFromCSVLine(line, c->work_network_rewire)) {
throw std::runtime_error("bad value for work_network_rewire during csv parsing @todo\n");
}
if (!valueFromCSVLine(line, c->mean_random_interactions_0_19)) {
throw std::runtime_error("bad value for mean_random_interactions_0_19 during csv parsing @todo\n");
}
Expand Down Expand Up @@ -292,7 +310,11 @@ void print(exateppabm::input::config config) {
fmt::print(" relative_transmission_household = {}\n", config.relative_transmission_household);
fmt::print(" relative_transmission_occupation = {}\n", config.relative_transmission_occupation);
fmt::print(" relative_transmission_random = {}\n", config.relative_transmission_occupation);
fmt::print(" mean_work_interactions_child = {}\n", config.mean_work_interactions_child);
fmt::print(" mean_work_interactions_adult = {}\n", config.mean_work_interactions_adult);
fmt::print(" mean_work_interactions_elderly = {}\n", config.mean_work_interactions_elderly);
fmt::print(" daily_fraction_work = {}\n", config.daily_fraction_work);
fmt::print(" work_network_rewire = {}\n", config.work_network_rewire);
fmt::print(" mean_random_interactions_0_19 = {}\n", config.mean_random_interactions_0_19);
fmt::print(" sd_random_interactions_0_19 = {}\n", config.sd_random_interactions_0_19);
fmt::print(" mean_random_interactions_20_69 = {}\n", config.mean_random_interactions_20_69);
Expand Down
46 changes: 35 additions & 11 deletions src/exateppabm/input.h
Original file line number Diff line number Diff line change
Expand Up @@ -100,7 +100,7 @@ struct config {
/**
* The mean time in days from exposed to infected state
* Default value is arbitrary
* @todo - might be equivalent to mean_time_to_syptoms, depending on how asym works in the reference model.
* @todo - might be equivalent to mean_time_to_symptoms, depending on how the more complex disease model works in the reference model.
*/
float mean_time_to_infected = 4;
/**
Expand Down Expand Up @@ -197,39 +197,63 @@ struct config {
*/
float relative_transmission_random = 1.0f;
/**
* Fraction of people in work network interacted with per day (by rng sampling)
* Mean daily interactions at school (0-19)
* This parameter is used during small world network generation (k)
* Arbitrary default value
*/
double mean_work_interactions_child = 10;
/**
* Mean daily interactions at adult workplaces (20-69)
* This parameter is used during small world network generation (k)
* Arbitrary default value
*/
double mean_work_interactions_adult = 7.5;
/**
* Mean daily interactions at equivalent to workplace for older people (70+)
* This parameter is used during small world network generation (k)
* Arbitrary default value
*/
double mean_work_interactions_elderly = 3;
/**
* Fraction of people in work network interacted with per day
*
* @todo - this probably needs using differently with more realistic networks
* This modifies the small world generation initial mean to create larger clusters, allowing for random sampling by this factor on each day
*
* Arbitrary default value
*/
float daily_fraction_work = 0.5f;

double daily_fraction_work = 0.5;
/**
* Small world network rewire parameter (probability) for workplaces, in range [0, 1]
* Arbitrary default
*/
double work_network_rewire = 0.1;
/**
* Mean number of random interactions for 0-19 year olds
* Mean number of random interactions for 0-19 age individuals
* Arbitrary default value
*/
double mean_random_interactions_0_19 = 2u;
/**
* Standard deviation for the number of random interactions per day for 0-19 year olds
* Standard deviation for the number of random interactions per day for 0-19 age individuals
* Default value is arbitrary
*/
double sd_random_interactions_0_19 = 2u;
/**
* Mean number of random interactions for 20-69 year olds
* Mean number of random interactions for 20-69 age individuals
* Arbitrary default value
*/
double mean_random_interactions_20_69 = 4u;
/**
* Standard deviation for the number of random interactions per day for 20-69 year olds
* Standard deviation for the number of random interactions per day for 20-69 age individuals
* Default value is arbitrary
*/
double sd_random_interactions_20_69 = 4u;
/**
* Mean number of random interactions for 70+ year olds
* Mean number of random interactions for 70+ age individuals
* Arbitrary default value
*/
double mean_random_interactions_70plus = 2u;
/**
* Standard deviation for the number of random interactions per day for 70+ year olds
* Standard deviation for the number of random interactions per day for 70+ age individuals
* Default value is arbitrary
*/
double sd_random_interactions_70plus = 2u;
Expand Down
122 changes: 122 additions & 0 deletions src/exateppabm/network.cu
Original file line number Diff line number Diff line change
@@ -0,0 +1,122 @@
#include "exateppabm/network.h"

#include <exception>
#include <iostream>
#include <string>
#include <utility>
#include <vector>
#include <unordered_set>
#include <random>

#include "flamegpu/flamegpu.h"
#include "fmt/core.h"

namespace exateppabm {
namespace network {


UndirectedNetwork generateFullyConnectedUndirectedNetwork(std::vector<flamegpu::id_t> nodes) {
UndirectedNetwork network = {};
network.nodes = nodes;
if (nodes.size() <= 1) {
return network;
}
std::uint32_t N = static_cast<std::int32_t>(nodes.size());
std::uint32_t E = (N * (N - 1)) / 2;
network.edges.reserve(E);
// Generate edges in ascending order for the undirected graph.
for (std::uint32_t i = 0; i < N; ++i) {
for (std::uint32_t j = i + 1; j < N; ++j) {
network.edges.push_back({i, j});
}
}
return network;
}

// @todo - should the initial lattice be randomised rather than sequential?
// @todo - should the graph be directed or undirected? Data structure is directed, so do pairs need rewiring together?
UndirectedNetwork generateSmallWorldUndirectedNetwork(std::vector<flamegpu::id_t> nodes, std::uint32_t K, double p_rewire, std::mt19937_64 rng) {
// fmt::print("@todo - check if small world network should be directed or undirected\n");
UndirectedNetwork network = {};
network.nodes = nodes;

// Return early with an empty network if 0 or 1 nodes, so no room for non-self edges
if (network.nodes.size() <= 1) {
return network;
}

// If K is 0 (or 1), no edges so do nothing
if (K <= 1) {
// throw std::runtime_error("@todo values - small world network K (" + std::to_string(K) + ") must be > 1");
return network;
} else if (K == network.nodes.size()) {
// If K == Node count, the graph is fully connected, so return a fully connected graph
return generateFullyConnectedUndirectedNetwork(nodes);
} else if (K >= network.nodes.size()) {
// Raise an exception if K is too large
throw std::runtime_error(std::string("@todo values - small world network K(") + std::to_string(K) + ") must be less than |N| (" + std::to_string(network.nodes.size()) + ")");
}
// If K is odd, use K-1
if (K % 2 == 1) {
K = K - 1;
}

// p_rewire must be between 0 and 1.
if (p_rewire < 0 || p_rewire > 1) {
throw std::runtime_error("generateSmallWorldNetwork p must be in [0, 1]. @todo include value");
}

// Initialise the edges as a lattice network, going K edges in each direction.
// Use signed integers to make modulo possible
std::int32_t N = static_cast<std::int32_t>(nodes.size());
std::int32_t E = (N * K) / 2;
network.edges.reserve(E);
for (std::int32_t i = 0; i < N; ++i) {
for (std::int32_t j = 1; j <= static_cast<std::int32_t>(K / 2); ++j) {
// only add the positive undirected edges
std::uint32_t s = static_cast<std::uint32_t>(i);
std::uint32_t d = static_cast<std::uint32_t>((i + j) % N);
// Ensure that the edge source node is lower than the destination, to simplify avoiding duplicate edges
if (s > d) {
std::swap(s, d);
}
network.edges.push_back({s, d});
// network.edges.push_back({static_cast<std::uint32_t>(i), static_cast<std::uint32_t>((i + j) % N)});
// network.edges.push_back({static_cast<std::uint32_t>(i), static_cast<std::uint32_t>((i - j + N) % N)});
}
}

// If the network is fully connected (enough edges for each node to be connected to every other node), rewiring is not needed, so adjust the probability to prevent infinite loops.
std::int32_t MAX_E = (N * (N - 1)) / 2;
if (E >= MAX_E) {
p_rewire = 0.0;
}

// If p_rewire is 0, no need to loop over the edges
if (p_rewire >= 0.0) {
// Get a uniform dist [0, 1)
std::uniform_real_distribution<double> p_dist(0.0, 1.0);
// Get a uniform integer distribution from [0, N) for generating new edge indices
std::uniform_int_distribution<std::uint32_t> dest_dist(0, N-1);
// Randomly rewire edges
for (auto& edge : network.edges) {
if (p_dist(rng) < p_rewire) {
std::uint32_t newDest = dest_dist(rng);
// Repeat the generation until the edge is not a self-edge, nor duplicate
// Not using do while so the search for the edge is ok.
// @note - this will be a challenge to parallelise efficiently without atomics (i.e. stable matching)
// @note - this will be expensive for large networks, an edge matrix would be cheaper (at the cost of N*N memory)
// @todo - need to avoid directed edge duplicates.
while (newDest == edge.source || std::find(network.edges.begin(), network.edges.end(), UndirectedNetwork::Edge{edge.source, newDest}) != network.edges.end()) {
newDest = dest_dist(rng);
}
edge.dest = newDest;
}
}
}

return network;
}

} // namespace network
} // namespace exateppabm
Loading

0 comments on commit 53b1fdb

Please sign in to comment.