499 lines
19 KiB
C++
499 lines
19 KiB
C++
#include "matrix.hpp"
|
|
#include <stdint.h>
|
|
#include <string.h>
|
|
#include <time.h>
|
|
#include <stdlib.h>
|
|
#include <cstdint>
|
|
#include <vector>
|
|
#include <algorithm>
|
|
|
|
std::vector<std::vector<uint64_t>> random_adjacency(const uint64_t vertex_count) {
|
|
std::vector<std::vector<uint64_t>> adjacency_matrix(vertex_count, std::vector<uint64_t>(vertex_count, 0));
|
|
srand(time(NULL));
|
|
|
|
for (uint64_t row_index = 0; row_index < vertex_count; row_index++) {
|
|
for (uint64_t column_index = 0; column_index < vertex_count; column_index++) {
|
|
if (column_index != row_index) {
|
|
adjacency_matrix[row_index][column_index] = rand() % 2;
|
|
}
|
|
}
|
|
}
|
|
|
|
return adjacency_matrix;
|
|
}
|
|
|
|
std::vector<std::vector<uint64_t>> calculate_distance_matrix(const std::vector<std::vector<uint64_t>>& adjacency_matrix) {
|
|
uint64_t vertex_count = adjacency_matrix.size();
|
|
std::vector<std::vector<uint64_t>> distance_matrix(vertex_count, std::vector<uint64_t>(vertex_count, 0));
|
|
std::vector<std::vector<uint64_t>> power_matrix(vertex_count, std::vector<uint64_t>(vertex_count, 0));
|
|
std::copy(adjacency_matrix.begin(), adjacency_matrix.end(), power_matrix.begin());
|
|
|
|
for (uint64_t row_index = 0; row_index < vertex_count; row_index++) {
|
|
for (uint64_t column_index = 0; column_index < vertex_count; column_index++) {
|
|
if (row_index == column_index) {
|
|
continue;
|
|
} else if (adjacency_matrix[row_index][column_index] == 1) {
|
|
distance_matrix[row_index][column_index] = 1;
|
|
} else {
|
|
distance_matrix[row_index][column_index] = UINT64_MAX;
|
|
}
|
|
}
|
|
}
|
|
|
|
for(uint64_t k = 2; k <= vertex_count; k++) {
|
|
power_matrix = gemm_basic(adjacency_matrix, power_matrix);
|
|
|
|
for (uint64_t row_index = 0; row_index < vertex_count; row_index++) {
|
|
for (uint64_t column_index = 0; column_index < vertex_count; column_index++) {
|
|
if (power_matrix[row_index][column_index] != 0 && distance_matrix[row_index][column_index] == UINT64_MAX) {
|
|
distance_matrix[row_index][column_index] = k;
|
|
}
|
|
}
|
|
}
|
|
}
|
|
|
|
return distance_matrix;
|
|
}
|
|
|
|
std::vector<uint64_t> get_eccentricities(const std::vector<std::vector<uint64_t>>& distance_matrix) {
|
|
uint64_t vertex_count = distance_matrix.size();
|
|
std::vector<uint64_t> eccentricities(vertex_count, UINT64_MAX);
|
|
uint64_t eccentricity;
|
|
|
|
for (uint64_t row_index = 0; row_index < vertex_count; row_index++) {
|
|
eccentricity = 0;
|
|
|
|
for (uint64_t column_index = 0; column_index < vertex_count; column_index++) {
|
|
if (distance_matrix[row_index][column_index] > eccentricity) {
|
|
eccentricity = distance_matrix[row_index][column_index];
|
|
}
|
|
}
|
|
|
|
if (eccentricity == 0) {
|
|
break;
|
|
}
|
|
eccentricities[row_index] = eccentricity;
|
|
}
|
|
|
|
return eccentricities;
|
|
}
|
|
|
|
uint64_t get_radius(const std::vector<uint64_t>& eccentricities) {
|
|
uint64_t radius = UINT64_MAX;
|
|
|
|
for (uint64_t index = 0; index < eccentricities.size(); index++) {
|
|
if (eccentricities[index] < radius) {
|
|
radius = eccentricities[index];
|
|
}
|
|
}
|
|
|
|
return radius;
|
|
}
|
|
|
|
uint64_t get_diameter(const std::vector<uint64_t>& eccentricities) {
|
|
uint64_t diamter = 0;
|
|
|
|
for (uint64_t index = 0; index < eccentricities.size(); index++) {
|
|
if (eccentricities[index] > diamter) {
|
|
diamter = eccentricities[index];
|
|
}
|
|
}
|
|
|
|
return diamter;
|
|
}
|
|
|
|
std::vector<uint64_t> get_centre(const std::vector<uint64_t>& eccentricities) {
|
|
std::vector<uint64_t> centre;
|
|
uint64_t radius = get_radius(eccentricities);
|
|
|
|
for (uint64_t index = 0; index < eccentricities.size(); index++) {
|
|
if (eccentricities[index] == radius) {
|
|
centre.push_back(index + 1);
|
|
}
|
|
}
|
|
|
|
centre.shrink_to_fit();
|
|
return centre;
|
|
}
|
|
|
|
std::vector<std::vector<uint64_t>> calculate_path_matrix(const std::vector<std::vector<uint64_t>>& adjacency_matrix) {
|
|
uint64_t vertex_count = adjacency_matrix.size();
|
|
std::vector<std::vector<uint64_t>> path_matrix(vertex_count, std::vector<uint64_t>(vertex_count, 0));
|
|
std::vector<std::vector<uint64_t>> power_matrix(vertex_count, std::vector<uint64_t>(vertex_count, 0));
|
|
std::copy(adjacency_matrix.begin(), adjacency_matrix.end(), power_matrix.begin());
|
|
|
|
for (uint64_t row_index = 0; row_index < vertex_count; row_index++) {
|
|
for (uint64_t column_index = 0; column_index < vertex_count; column_index++) {
|
|
if (row_index == column_index || adjacency_matrix[row_index][column_index] == 1) {
|
|
path_matrix[row_index][column_index] = 1;
|
|
}
|
|
}
|
|
}
|
|
|
|
for(uint64_t k = 2; k <= vertex_count; k++) {
|
|
power_matrix = gemm_basic(adjacency_matrix, power_matrix);
|
|
|
|
for (uint64_t row_index = 0; row_index < vertex_count; row_index++) {
|
|
for (uint64_t column_index = 0; column_index < vertex_count; column_index++) {
|
|
if (power_matrix[row_index][column_index] != 0) {
|
|
path_matrix[row_index][column_index] = 1;
|
|
}
|
|
}
|
|
}
|
|
}
|
|
|
|
return path_matrix;
|
|
}
|
|
|
|
std::vector<std::vector<uint64_t>> find_components_basic(const std::vector<std::vector<uint64_t>>& path_matrix) {
|
|
uint64_t vertex_count = path_matrix.size();
|
|
std::vector<std::vector<uint64_t>> components;
|
|
std::vector<uint64_t> component;
|
|
components.reserve(0);
|
|
component.reserve(vertex_count);
|
|
|
|
for (uint64_t row_index = 0; row_index < vertex_count; row_index++) {
|
|
component.clear();
|
|
|
|
for (uint64_t column_index = 0; column_index < vertex_count; column_index++) {
|
|
if (path_matrix[row_index][column_index] == 1) {
|
|
component.push_back(column_index + 1);
|
|
}
|
|
}
|
|
|
|
component.shrink_to_fit();
|
|
std::sort(component.begin(), component.end());
|
|
|
|
auto iterator = find(components.begin(), components.end(), component);
|
|
|
|
if (iterator == components.end()) {
|
|
components.push_back(component);
|
|
}
|
|
}
|
|
|
|
components.shrink_to_fit();
|
|
return components;
|
|
}
|
|
|
|
void dfs(const std::vector<std::vector<uint64_t>>& adjacency_matrix, const uint64_t vertex, std::vector<bool>& visited) {
|
|
visited[vertex] = true;
|
|
|
|
for (uint64_t neighbor_vertex = 0; neighbor_vertex < adjacency_matrix.size(); neighbor_vertex++) {
|
|
if (adjacency_matrix[vertex][neighbor_vertex] != 1) {
|
|
continue;
|
|
}
|
|
if (visited[neighbor_vertex]) {
|
|
continue;
|
|
}
|
|
dfs(adjacency_matrix, neighbor_vertex, visited);
|
|
}
|
|
}
|
|
|
|
std::vector<std::vector<uint64_t>> find_components_dfs(const std::vector<std::vector<uint64_t>>& adjacency_matrix) {
|
|
uint64_t vertex_count = adjacency_matrix.size();
|
|
std::vector<std::vector<uint64_t>> components;
|
|
std::vector<uint64_t> component;
|
|
std::vector<bool> visited(vertex_count, false);
|
|
components.reserve(0);
|
|
component.reserve(vertex_count);
|
|
|
|
for (uint64_t vertex = 0; vertex < vertex_count; vertex++) {
|
|
component.clear();
|
|
|
|
dfs(adjacency_matrix, vertex, visited);
|
|
|
|
for (uint64_t index = 0; index < vertex_count; index++) {
|
|
if (visited[index]) {
|
|
component.push_back(index + 1);
|
|
}
|
|
}
|
|
|
|
component.shrink_to_fit();
|
|
std::sort(component.begin(), component.end());
|
|
|
|
auto iterator = find(components.begin(), components.end(), component);
|
|
|
|
if (iterator == components.end()) {
|
|
components.push_back(component);
|
|
}
|
|
}
|
|
|
|
components.shrink_to_fit();
|
|
return components;
|
|
}
|
|
|
|
std::vector<std::vector<uint64_t>> find_bridges_basic(const std::vector<std::vector<uint64_t>>& adjacency_matrix) {
|
|
uint64_t vertex_count = adjacency_matrix.size();
|
|
std::vector<std::vector<uint64_t>> path_matrix(vertex_count, std::vector<uint64_t>(vertex_count, 0));
|
|
std::vector<std::vector<uint64_t>> temp_adjacency_matrix(vertex_count, std::vector<uint64_t>(vertex_count, 0));
|
|
std::vector<std::vector<uint64_t>> components = find_components_basic(path_matrix);
|
|
std::vector<std::vector<uint64_t>> temp_components;
|
|
std::vector<std::vector<uint64_t>> bridges;
|
|
std::vector<uint64_t> bridge(2, 0);
|
|
|
|
for (uint64_t row_index = 0; row_index < vertex_count; row_index++) {
|
|
for (uint64_t column_index = 0; column_index < vertex_count; column_index++) {
|
|
if (row_index == column_index) {
|
|
continue;
|
|
}
|
|
|
|
std::copy(adjacency_matrix.begin(), adjacency_matrix.end(), temp_adjacency_matrix.begin());
|
|
|
|
temp_adjacency_matrix[row_index][column_index] = 0;
|
|
temp_adjacency_matrix[column_index][row_index] = 0;
|
|
|
|
path_matrix = calculate_path_matrix(temp_adjacency_matrix);
|
|
temp_components = find_components_basic(path_matrix);
|
|
|
|
if (temp_components.size() <= components.size()) {
|
|
continue;
|
|
}
|
|
|
|
if (column_index < row_index) {
|
|
bridge[0] = column_index + 1;
|
|
bridge[1] = row_index + 1;
|
|
} else {
|
|
bridge[0] = row_index + 1;
|
|
bridge[1] = column_index + 1;
|
|
}
|
|
|
|
auto iterator = find(bridges.begin(), bridges.end(), bridge);
|
|
|
|
if (iterator == bridges.end()) {
|
|
bridges.push_back(bridge);
|
|
}
|
|
}
|
|
}
|
|
|
|
bridges.shrink_to_fit();
|
|
return bridges;
|
|
}
|
|
|
|
std::vector<std::vector<uint64_t>> find_bridges_dfs_v1(const std::vector<std::vector<uint64_t>>& adjacency_matrix) {
|
|
uint64_t vertex_count = adjacency_matrix.size();
|
|
std::vector<std::vector<uint64_t>> temp_adjacency_matrix(vertex_count, std::vector<uint64_t>(vertex_count, 0));
|
|
std::vector<std::vector<uint64_t>> components = find_components_dfs(adjacency_matrix);
|
|
std::vector<std::vector<uint64_t>> temp_components;
|
|
std::vector<std::vector<uint64_t>> bridges;
|
|
std::vector<uint64_t> bridge(2, 0);
|
|
|
|
for (uint64_t row_index = 0; row_index < vertex_count; row_index++) {
|
|
for (uint64_t column_index = 0; column_index < vertex_count; column_index++) {
|
|
if (row_index == column_index) {
|
|
continue;
|
|
}
|
|
|
|
std::copy(adjacency_matrix.begin(), adjacency_matrix.end(), temp_adjacency_matrix.begin());
|
|
|
|
temp_adjacency_matrix[row_index][column_index] = 0;
|
|
temp_adjacency_matrix[column_index][row_index] = 0;
|
|
|
|
temp_components = find_components_dfs(temp_adjacency_matrix);
|
|
|
|
if (temp_components.size() <= components.size()) {
|
|
continue;
|
|
}
|
|
|
|
if (column_index < row_index) {
|
|
bridge[0] = column_index + 1;
|
|
bridge[1] = row_index + 1;
|
|
} else {
|
|
bridge[0] = row_index + 1;
|
|
bridge[1] = column_index + 1;
|
|
}
|
|
|
|
auto iterator = find(bridges.begin(), bridges.end(), bridge);
|
|
|
|
if (iterator == bridges.end()) {
|
|
bridges.push_back(bridge);
|
|
}
|
|
}
|
|
}
|
|
|
|
bridges.shrink_to_fit();
|
|
return bridges;
|
|
}
|
|
|
|
void dfs_bridges(const std::vector<std::vector<uint64_t>>& adjacency_matrix, const uint64_t vertex, const uint64_t parent_vertex, std::vector<bool>& visited,
|
|
uint64_t current_time, std::vector<uint64_t>& discovery_time, std::vector<uint64_t>& lowest_time, std::vector<std::vector<uint64_t>>& bridges) {
|
|
current_time++;
|
|
visited[vertex] = true;
|
|
discovery_time[vertex] = current_time;
|
|
lowest_time[vertex] = current_time;
|
|
std::vector<uint64_t> bridge(2, 0);
|
|
|
|
for (uint64_t neighbor_vertex = 0; neighbor_vertex < adjacency_matrix.size(); neighbor_vertex++) {
|
|
if (adjacency_matrix[vertex][neighbor_vertex] != 1) {
|
|
continue;
|
|
}
|
|
if (parent_vertex != neighbor_vertex && visited[neighbor_vertex]) {
|
|
if (lowest_time[vertex] > discovery_time[neighbor_vertex]) {
|
|
lowest_time[vertex] = discovery_time[neighbor_vertex];
|
|
}
|
|
continue;
|
|
}
|
|
if (visited[neighbor_vertex]) {
|
|
continue;
|
|
}
|
|
|
|
dfs_bridges(adjacency_matrix, neighbor_vertex, vertex, visited, current_time, discovery_time, lowest_time, bridges);
|
|
|
|
if (lowest_time[vertex] > lowest_time[neighbor_vertex]) {
|
|
lowest_time[vertex] = lowest_time[neighbor_vertex];
|
|
}
|
|
if (discovery_time[vertex] >= lowest_time[neighbor_vertex]) {
|
|
continue;
|
|
}
|
|
bridge[0] = vertex + 1;
|
|
bridge[1] = neighbor_vertex + 1;
|
|
|
|
auto iterator = find(bridges.begin(), bridges.end(), bridge);
|
|
|
|
if (iterator == bridges.end()) {
|
|
bridges.push_back(bridge);
|
|
}
|
|
}
|
|
}
|
|
|
|
std::vector<std::vector<uint64_t>> find_bridges_dfs_v2(const std::vector<std::vector<uint64_t>>& adjacency_matrix) {
|
|
uint64_t vertex_count = adjacency_matrix.size();
|
|
std::vector<bool> visited(vertex_count, false);
|
|
uint64_t current_time = 0;
|
|
std::vector<uint64_t> discovery_time(vertex_count, 0);
|
|
std::vector<uint64_t> lowest_time(vertex_count, 0);
|
|
std::vector<std::vector<uint64_t>> bridges;
|
|
|
|
for (uint64_t vertex = 0; vertex < vertex_count; vertex++) {
|
|
if (!visited[vertex]) {
|
|
dfs_bridges(adjacency_matrix, vertex, UINT64_MAX, visited, current_time, discovery_time, lowest_time, bridges);
|
|
}
|
|
}
|
|
|
|
bridges.shrink_to_fit();
|
|
return bridges;
|
|
}
|
|
|
|
std::vector<uint64_t> find_articulations_basic(const std::vector<std::vector<uint64_t>>& adjacency_matrix) {
|
|
uint64_t vertex_count = adjacency_matrix.size();
|
|
std::vector<std::vector<uint64_t>> path_matrix(vertex_count, std::vector<uint64_t>(vertex_count, 0));
|
|
std::vector<std::vector<uint64_t>> temp_adjacency_matrix(vertex_count, std::vector<uint64_t>(vertex_count, 0));
|
|
std::vector<std::vector<uint64_t>> components = find_components_basic(path_matrix);
|
|
std::vector<std::vector<uint64_t>> temp_components;
|
|
std::vector<uint64_t> articulations;
|
|
|
|
for (uint64_t i = 0; i < vertex_count; i++) {
|
|
std::copy(adjacency_matrix.begin(), adjacency_matrix.end(), temp_adjacency_matrix.begin());
|
|
|
|
for (uint64_t row_index = 0; row_index < vertex_count; row_index++) {
|
|
for (uint64_t column_index = 0; column_index < vertex_count; column_index++) {
|
|
temp_adjacency_matrix[row_index][i] = 0;
|
|
temp_adjacency_matrix[i][column_index] = 0;
|
|
}
|
|
}
|
|
|
|
path_matrix = calculate_path_matrix(temp_adjacency_matrix);
|
|
temp_components = find_components_basic(path_matrix);
|
|
|
|
// the + 1 is needed because I am not removing the vertex, I am just removing all of its edges
|
|
// removing all of its edges, means it itself becomes a component, which needs to be accounted for
|
|
if (temp_components.size() > components.size() + 1) {
|
|
articulations.push_back(i + 1);
|
|
}
|
|
}
|
|
|
|
articulations.shrink_to_fit();
|
|
std::sort(articulations.begin(), articulations.end());
|
|
return articulations;
|
|
}
|
|
|
|
std::vector<uint64_t> find_articulations_dfs_v1(const std::vector<std::vector<uint64_t>>& adjacency_matrix) {
|
|
uint64_t vertex_count = adjacency_matrix.size();
|
|
std::vector<std::vector<uint64_t>> path_matrix(vertex_count, std::vector<uint64_t>(vertex_count, 0));
|
|
std::vector<std::vector<uint64_t>> temp_adjacency_matrix(vertex_count, std::vector<uint64_t>(vertex_count, 0));
|
|
std::vector<std::vector<uint64_t>> components = find_components_dfs(adjacency_matrix);
|
|
std::vector<std::vector<uint64_t>> temp_components;
|
|
std::vector<uint64_t> articulations;
|
|
|
|
for (uint64_t i = 0; i < vertex_count; i++) {
|
|
std::copy(adjacency_matrix.begin(), adjacency_matrix.end(), temp_adjacency_matrix.begin());
|
|
|
|
for (uint64_t row_index = 0; row_index < vertex_count; row_index++) {
|
|
for (uint64_t column_index = 0; column_index < vertex_count; column_index++) {
|
|
temp_adjacency_matrix[row_index][i] = 0;
|
|
temp_adjacency_matrix[i][column_index] = 0;
|
|
}
|
|
}
|
|
|
|
temp_components = find_components_dfs(temp_adjacency_matrix);
|
|
|
|
// the + 1 is needed because I am not removing the vertex, I am just removing all of its edges
|
|
// removing all of its edges, means it itself becomes a component, which needs to be accounted for
|
|
if (temp_components.size() > components.size() + 1) {
|
|
articulations.push_back(i + 1);
|
|
}
|
|
}
|
|
|
|
articulations.shrink_to_fit();
|
|
std::sort(articulations.begin(), articulations.end());
|
|
return articulations;
|
|
}
|
|
|
|
void dfs_articulations(const std::vector<std::vector<uint64_t>>& adjacency_matrix, const uint64_t vertex, const uint64_t parent_vertex, std::vector<bool>& visited,
|
|
uint64_t current_time, std::vector<uint64_t>& discovery_time, std::vector<uint64_t>& lowest_time, std::vector<uint64_t>& articulations) {
|
|
current_time++;
|
|
visited[vertex] = true;
|
|
discovery_time[vertex] = current_time;
|
|
lowest_time[vertex] = current_time;
|
|
uint64_t child_count = 0;
|
|
bool is_articulation = false;
|
|
|
|
for (uint64_t neighbor_vertex = 0; neighbor_vertex < adjacency_matrix.size(); neighbor_vertex++) {
|
|
if (adjacency_matrix[vertex][neighbor_vertex] != 1) {
|
|
continue;
|
|
}
|
|
if (visited[neighbor_vertex]) {
|
|
if (lowest_time[vertex] > discovery_time[neighbor_vertex]) {
|
|
lowest_time[vertex] = discovery_time[neighbor_vertex];
|
|
}
|
|
continue;
|
|
}
|
|
|
|
child_count++;
|
|
|
|
dfs_articulations(adjacency_matrix, neighbor_vertex, vertex, visited, current_time, discovery_time, lowest_time, articulations);
|
|
|
|
if (lowest_time[vertex] > lowest_time[neighbor_vertex]) {
|
|
lowest_time[vertex] = lowest_time[neighbor_vertex];
|
|
}
|
|
|
|
if (parent_vertex != UINT64_MAX && discovery_time[vertex] <= lowest_time[neighbor_vertex]) {
|
|
is_articulation = true;
|
|
}
|
|
}
|
|
if (parent_vertex == UINT64_MAX && child_count > 1) {
|
|
is_articulation = true;
|
|
}
|
|
if (is_articulation) {
|
|
articulations.push_back(vertex + 1);
|
|
}
|
|
}
|
|
|
|
std::vector<uint64_t> find_articulations_dfs_v2(const std::vector<std::vector<uint64_t>>& adjacency_matrix) {
|
|
uint64_t vertex_count = adjacency_matrix.size();
|
|
std::vector<bool> visited(vertex_count, false);
|
|
uint64_t current_time = 0;
|
|
std::vector<uint64_t> discovery_time(vertex_count, 0);
|
|
std::vector<uint64_t> lowest_time(vertex_count, 0);
|
|
std::vector<uint64_t> articulations;
|
|
|
|
for (uint64_t vertex = 0; vertex < vertex_count; vertex++) {
|
|
if (!visited[vertex]) {
|
|
dfs_articulations(adjacency_matrix, vertex, UINT64_MAX, visited, current_time, discovery_time, lowest_time, articulations);
|
|
}
|
|
}
|
|
|
|
articulations.shrink_to_fit();
|
|
std::sort(articulations.begin(), articulations.end());
|
|
return articulations;
|
|
}
|