GPU621/DPS921 T-eaSyPeasy
Revision as of 07:39, 14 April 2016 by Bernard Wouter de Vries Robles (talk | contribs) (→Source code for the assignment)
T-eaSyPeasy
Team Members
Progress
Introduction
For my project instead of picking a new method to learn to parallelise with, I decided to try to paralellise an algorithm that is fairly hard to write, to see if after following this course I could actually parallelise an algorithm. The thing is, in the workshops we have only been parallelising algorithms that have been specifically written to be parallelised. I wanted to change that and see if I could still do it. The algorithm I am parallelising is a Branch & Bound Algorithm by Little from 1966. A good explanation of the algorithm can be found here: https://www.youtube.com/watch?v=nN4K8xA8ShM, the explanation of the branch and bound I wrote starts at 34:30.
Assignment 2
Source
#include <stdint.h>
#include <vector>
#include <algorithm>
#include <fstream>
#include <iostream>
#include <cmath>
#include <time.h>
struct Selection {
uint16_t row;
uint16_t col;
uint16_t val;
Selection(uint16_t x, uint16_t y, uint16_t v);
bool operator==(const Selection& rhs)
{
return (col == rhs.col && row == rhs.row);
}
bool operator!=(const Selection& rhs)
{
return !(col == rhs.col && row == rhs.row);
}
Selection() : row(UINT16_MAX), col(UINT16_MAX), val(UINT16_MAX) {}
};
Selection::Selection(uint16_t r, uint16_t c, uint16_t v) {
row = r;
col = c;
val = v;
}
struct Zero {
uint16_t row;
uint16_t col;
uint16_t horPen;
uint16_t vertPen;
Zero(uint16_t row, uint16_t col);
Zero() {};
};
Zero::Zero(uint16_t r, uint16_t c) {
row = r;
col = c;
};
struct Node {
int lowerBound;
int matrixSize;
bool solutionFound = false;
bool invalidSolution = false;
uint16_t* originalMatrix;
uint16_t* workingCopy;
Node(int amountOfCities, uint16_t* matrix, std::vector<std::vector<Selection>> solution, uint16_t lb);
Node(int amountOfCities, uint16_t* matrix, std::vector<std::vector<Selection>> solution);
std::vector<std::vector<Selection>> solution;
Node* left;
Node* right;
void evaluate(){
std::vector<Selection> rowMinima;
std::vector<Selection> colMinima;
std::vector<Zero> zeroes;
std::vector<bool> zeroInColumn(matrixSize, false);
lowerBound = sumCurrentSolution();
createWorkingCopy();
findRowMinima(rowMinima, zeroInColumn);
lowerBound += sumMinima(rowMinima);
subtractMinima(rowMinima, true);
findColMinima(colMinima, zeroInColumn);
lowerBound += sumMinima(colMinima);
subtractMinima(colMinima, false);
gatherZeroes(rowMinima, colMinima, zeroes);
if (zeroes.empty()) {
invalidSolution = true;
}
else {
findPenalties(zeroes);
Zero assignment = assignLargestPenalty(zeroes);
createLeftNode(solution, assignment);
int subtourIndex = addToResult(solution, assignment);
createRightNode(solution, assignment, subtourIndex);
checkSolutionFound();
}
free(workingCopy);
}
int sumCurrentSolution() {
int result = 0;
std::for_each(solution.begin(), solution.end(), [&](std::vector<Selection> v) {
std::for_each(v.begin(), v.end(), [&](Selection s) {
result += s.val;
});
});
return result;
}
void checkSolutionFound() {
solutionFound = (solution.size() == 1 && solution[0].size() == matrixSize - 1);
}
bool isSolutionFound() {
return solutionFound;
}
std::vector<std::vector<Selection>> getSolution() {
return solution;
}
void createWorkingCopy() {
workingCopy = (uint16_t*)malloc(sizeof(uint16_t) * matrixSize * matrixSize);
memcpy(workingCopy, originalMatrix, (sizeof(uint16_t) * matrixSize * matrixSize));
}
void findRowMinima(std::vector<Selection>& rowMinima, std::vector<bool>& zeroInColumn) {
for (int row = 0; row < matrixSize; row++) {
int rowOffset = row * matrixSize;
uint16_t currentMinimum = UINT16_MAX;
Selection currentSelection;
for (int col = 0; col < matrixSize; col++) {
if (workingCopy[rowOffset + col] < currentMinimum) {
currentMinimum = workingCopy[rowOffset + col];
currentSelection = Selection(row, col, currentMinimum);
zeroInColumn[col] = true;
}
}
if (currentMinimum != UINT16_MAX) {
rowMinima.push_back(currentSelection);
}
}
//std::for_each(rowMinima.begin(), rowMinima.end(), [&](Selection s) {
// std::cout << " ( " << s.row << ", " << s.col << ") val = " << s.val << " \n";
//});
}
int sumMinima(const std::vector<Selection>& minima) {
int sum_of_elems = 0;
std::for_each(minima.begin(), minima.end(), [&](Selection s) {
sum_of_elems += s.val;
});
return sum_of_elems;
}
void subtractMinima(std::vector<Selection> minima, const bool row) {
if (row) {
std::for_each(minima.begin(), minima.end(), [&](Selection s) {
if (s.val != 0) {
for (int j = 0; j < matrixSize; j++) {
if (workingCopy[s.row*matrixSize + j] != UINT16_MAX) {
workingCopy[s.row*matrixSize + j] -= s.val;
}
}
}
});
}
else {
std::for_each(minima.begin(), minima.end(), [&](Selection s) {
if (s.val != 0) {
for (int row = 0; row < matrixSize; row++) {
if (workingCopy[row*matrixSize + s.col] != UINT16_MAX) {
workingCopy[row*matrixSize + s.col] -= s.val;
}
}
}
});
}
}
void findColMinima(std::vector<Selection>& colMinima, const std::vector<bool>& zeroInColumn) {
for (int col = 0; col < matrixSize; col++) {
if (!zeroInColumn[col]) {
uint16_t currentMinimum = UINT16_MAX;
Selection currentSelection;
for (int row = 0; row < matrixSize; row++) {
if (workingCopy[row*matrixSize + col] < currentMinimum) {
currentMinimum = workingCopy[row*matrixSize + col];
currentSelection = Selection(row, col, currentMinimum);
}
}
if (currentMinimum != UINT16_MAX) {
colMinima.push_back(currentSelection);
}
}
}
//std::for_each(colMinima.begin(), colMinima.end(), [&](Selection s) {
// std::cout << " ( " << s.row << ", " << s.col << ") val = " << s.val << " \n";
//});
}
void printMatrix() {
std::cout << "\n Matrix: ";
for (int i = 0; i < matrixSize; i++) {
std::cout << "\n";
for (int j = 0; j < matrixSize; j++) {
std::cout << workingCopy[i*matrixSize + j] << " ";
}
}
std::cout << "\n";
}
void printOriginalMatrix() {
std::cout << "\nOriginal Matrix: ";
for (int i = 0; i < matrixSize; i++) {
std::cout << "\n";
for (int j = 0; j < matrixSize; j++) {
std::cout << originalMatrix[i*matrixSize + j] << " ";
}
}
std::cout << "\n";
}
void printSolution() {
std::cout << "\n intermediate :\n";
std::for_each(solution.begin(), solution.end(), [&](std::vector<Selection> v) {
std::for_each(v.begin(), v.end(), [&](Selection s) {
std::cout << "(" << s.row << ", " << s.col << ") value = " << s.val << " ";
});
std::cout << "\n";
});
std::cout << "\n";
}
void gatherZeroes(std::vector<Selection>& rowMinima, std::vector<Selection>& colMinima, std::vector<Zero>& zeroes) {
for (int i = 0; i < matrixSize; i++) {
for (int j = 0; j < matrixSize; j++) {
if (workingCopy[i*matrixSize + j] == 0) {
zeroes.push_back(Zero(i, j));
}
}
}
//std::for_each(zeroes.begin(), zeroes.end(), [&](Zero z) {
// std::cout << " zero found : row : " << z.row << " col : " << z.col << "\n";
//});
}
void findPenalties(std::vector<Zero>& zeroes) {
std::for_each(zeroes.begin(), zeroes.end(), [&](Zero& z) {
uint16_t horizontalMin = UINT16_MAX;
for (int col = 0; col < matrixSize; col++) {
if (col != z.col && workingCopy[z.row*matrixSize + col] < horizontalMin) {
horizontalMin = workingCopy[z.row*matrixSize + col];
}
}
z.horPen = horizontalMin;
uint16_t verticalMin = UINT16_MAX;
for (int row = 0; row < matrixSize; row++) {
if (row != z.row && workingCopy[row * matrixSize + z.col] < verticalMin) {
verticalMin = workingCopy[row * matrixSize + z.col];
}
}
z.vertPen = verticalMin;
});
//std::for_each(zeroes.begin(), zeroes.end(), [&](Zero z) {
// std::cout << " (" << z.row << ", " << z.col << ") Horizontal penalty: " << z.horPen << " Vertical penalty: " << z.vertPen << "\n";
//});
}
Zero assignLargestPenalty(std::vector<Zero>& zeroes) {
Zero assignment = zeroes[0];
std::for_each(zeroes.begin(), zeroes.end(), [&](Zero z) {
if ((z.vertPen + z.horPen) > (assignment.vertPen + assignment.horPen)) {
assignment = z;
}
});
return assignment;
}
int addToResult(std::vector<std::vector<Selection>>& solution, Zero assignment) {
int startOfSubtourIndex = -1;
int endOfSubtourIndex = -1;
int index = 0;
int result = 0;
//dealing with subtours....
std::for_each(solution.begin(), solution.end(), [&](std::vector<Selection>& subtour) {
if (subtour[0].row == assignment.col) {
startOfSubtourIndex = index;
}
else if (subtour.back().col == assignment.row) {
endOfSubtourIndex = index;
}
index++;
});
if (startOfSubtourIndex > -1 && endOfSubtourIndex > -1) {
//this case is true if our solution connects 2 subtours, by being at the start of one and at the end of the other.
//we combine them by concatenating the subtour that it is the start of after the subtour that it is the end of.
solution[endOfSubtourIndex].push_back(Selection(assignment.row, assignment.col, originalMatrix[assignment.row* matrixSize + assignment.col]));
solution[endOfSubtourIndex].insert(solution[endOfSubtourIndex].end(), solution[startOfSubtourIndex].begin(), solution[startOfSubtourIndex].end());
solution.erase(solution.begin() + startOfSubtourIndex);
result = endOfSubtourIndex;
if (startOfSubtourIndex < endOfSubtourIndex) {
result--;
}
}
else if (startOfSubtourIndex > -1) {
solution[startOfSubtourIndex].insert(solution[startOfSubtourIndex].begin(), Selection(assignment.row, assignment.col, originalMatrix[assignment.row * matrixSize + assignment.col]));
//subtour elimination, we can no longer go FROM(column) subtour.back().row TO(row) assignment.col.
result = startOfSubtourIndex;
}
else if (endOfSubtourIndex > -1) {
solution[endOfSubtourIndex].push_back(Selection(assignment.row, assignment.col, originalMatrix[assignment.row* matrixSize + assignment.col]));
//subtour elimination, we can no longer go FROM(column) assignment.row TO(row) subtour[0].col.
result = endOfSubtourIndex;
}
else {
solution.push_back(std::vector<Selection>(1, Selection(assignment.row, assignment.col, originalMatrix[assignment.row* matrixSize + assignment.col])));
result = solution.size() - 1;
}
return result;
}
void createLeftNode(std::vector<std::vector<Selection>>& solution, Zero assignment) {
uint16_t* leftSideMatrix = (uint16_t*)malloc(sizeof(uint16_t)*matrixSize*matrixSize);
for (int row = 0; row < matrixSize; row++) {
for (int col = 0; col < matrixSize; col++) {
if (row == assignment.row && col == assignment.col) {
leftSideMatrix[row*matrixSize + col] = UINT16_MAX;
}
else {
leftSideMatrix[row*matrixSize + col] = originalMatrix[row*matrixSize + col];
}
}
}
left = new Node(matrixSize, leftSideMatrix, solution, lowerBound + assignment.horPen + assignment.vertPen);
}
void createRightNode(std::vector<std::vector<Selection>>& solution, Zero assignment, int subtourIndex) {
uint16_t* rightSideMatrix = (uint16_t*)malloc(sizeof(uint16_t)*matrixSize*matrixSize);
for (int row = 0; row < matrixSize; row++) {
for (int col = 0; col < matrixSize; col++) {
if (assignment.row == row || assignment.col == col) {
rightSideMatrix[row*matrixSize + col] = UINT16_MAX;
}
else {
rightSideMatrix[row*matrixSize + col] = originalMatrix[row*matrixSize + col];
}
}
}
rightSideMatrix[solution[subtourIndex].back().col * matrixSize + solution[subtourIndex][0].row] = UINT16_MAX;
right = new Node(matrixSize, rightSideMatrix, solution);
}
Node getRight() {
return *right;
}
Node getLeft() {
return *left;
}
};
Node::Node(int amountOfCities, uint16_t* matrix, std::vector<std::vector<Selection>> solution) : solution(solution) {
matrixSize = amountOfCities;
originalMatrix = matrix;
}
Node::Node(int amountOfCities,uint16_t* matrix, std::vector<std::vector<Selection>> solution, uint16_t lb) : solution(solution){
matrixSize = amountOfCities;
originalMatrix = matrix;
lowerBound = lb;
}
class Tree {
Node *start;
int lowerBound;
int upperBound;
std::vector<Node> nodes;
};
void getMatrixFromFile(std::string filename, int matrixSize, uint16_t* originalMatrix) {
std::vector<int> xValues(matrixSize);
std::vector<int> yValues(matrixSize);
int n, m, i = 0;
std::ifstream read(filename);
while (i < matrixSize && read >> n >> m) {
xValues[i] = n;
yValues[i] = m;
i++;
}
for (int i = 0; i < matrixSize; i++) {
for (int j = 0; j < matrixSize; j++) {
if (i == j) {
originalMatrix[i*matrixSize + j] = UINT16_MAX;
}
else {
originalMatrix[i*matrixSize + j] = (uint16_t)sqrt(pow(xValues[i]-xValues[j], 2) + pow(yValues[i] - yValues[j],2));
}
}
}
}
void pruneNodeList(std::vector<Node>& nodeList, int upperBound) {
int i = 0;
while (i < nodeList.size()) {
if (nodeList[i].lowerBound >= upperBound) {
nodeList.erase(nodeList.begin() + i);
}
else {
++i;
}
}
}
int calculateSolution(uint16_t* originalMatrix, std::vector<std::vector<Selection>>& solution, int matrixSize) {
int sum = 0;
std::for_each(solution.begin(), solution.end(), [&](std::vector<Selection> v) {
std::for_each(v.begin(), v.end(), [&](Selection s) {
sum += s.val;
});
});
sum += originalMatrix[solution[0].back().col * matrixSize + solution[0][0].row];
return sum;
}
int findNextEvaluation(std::vector<Node>& nodeList) {
Node eval = nodeList[0];
int index = 0;
int currentIndex = 0;
std::for_each(nodeList.begin(), nodeList.end(), [&](Node n) {
if (n.lowerBound < eval.lowerBound) {
currentIndex = index;
}
index++;
});
return currentIndex;
}
int main() {
int amountOfCities = 15;
uint16_t* originalMatrix = (uint16_t*)malloc(sizeof(uint16_t)*amountOfCities*amountOfCities);
getMatrixFromFile("cities.txt", amountOfCities, originalMatrix);
std::vector<std::vector<Selection>> solution;
Node node = Node(amountOfCities, originalMatrix, solution);
std::vector<Node> nodeList;
std::vector<std::vector<Selection>> currentSolution;
std::vector<std::vector<Selection>> tempSolution;
int upperBound = UINT16_MAX;
int count = 0;
clock_t t1, t2;
t1 = clock();
while (!node.isSolutionFound()) {
node.evaluate();
count++;
if (!node.isSolutionFound()) {
nodeList.push_back(node.getLeft());
node = node.getRight();
}
}
currentSolution = node.getSolution();
upperBound = calculateSolution(originalMatrix, currentSolution, amountOfCities);
pruneNodeList(nodeList, upperBound);
while (!nodeList.empty()) {
int nextNode = findNextEvaluation(nodeList);
node = nodeList[nextNode];
nodeList.erase(nodeList.begin() + nextNode);
while (!node.isSolutionFound() && !node.invalidSolution) {
node.evaluate();
count++;
if (!node.isSolutionFound() && !node.invalidSolution && node.lowerBound < upperBound) {
if (node.getLeft().lowerBound < upperBound) {
nodeList.push_back(node.getLeft());
}
node = node.getRight();
}
else if (node.lowerBound >= upperBound) {
break;
} else if(!node.invalidSolution){
tempSolution = node.getSolution();
int value = calculateSolution(originalMatrix, tempSolution, amountOfCities);
if (value < upperBound) {
currentSolution = tempSolution;
upperBound = value;
pruneNodeList(nodeList, upperBound);
}
}
}
}
t2 = clock();
float diff((float)t2 - (float)t1);
std::cout << "time spent: " << diff << "\n";
std::cout << "count = " << count << "\n";
std::cout << "\n Solution :\n";
std::for_each(currentSolution.begin(), currentSolution.end(), [&](std::vector<Selection> v) {
std::for_each(v.begin(), v.end(), [&](Selection s) {
std::cout << "(" << s.row << ", " << s.col << ") ";
});
std::cout << "\n";
});
std::cout << "\nResult = " << upperBound;
}