I have a graph with N vertices and M edges (N is between 1 and 15 and M is between 1 and N^2). The graph is directed and weighted (with a probability for that excact edge). You are given a start vertex and a number of edges. The program is then going to calculate the probability for each vertex being the end vertex.
Examle input:
3 3 //Number of vertices and number of edges
1 2 0.4 //Edge nr.1 from vertex 1 to 2 with a probability of 0.4
1 3 0.5 //Edge nr.2 from vertex 1 to 3 with a probability of 0.5
2 1 0.8 //Edge nr.3...
3 //Number of questions
2 1 //Start vertex, number of edges to visit
1 1
1 2
Output:
0.8 0.2 0.0 //The probability for vertex 1 beign the last vertex is 0.8 for vertex 2 it is 0.2 and for vertex 3 it is 0.0
0.1 0.4 0.5
0.33 0.12 0.55
I have used a DFS in my solution, but when number of edges to visit can be up to 1 billion, this is way too slow... I have been looking at DP but I am not sure about how to implement it for this particular problem (if it is even the right way to solve it). So I was hoping that some of you could suggest an alternative to DFS and/or perhaps a way of using/implementing DP.
(I know it might be a bit messy, I have only been programming in C++ for a month)
#include <iostream>
#include <vector>
#include <stack>
using namespace std;
struct bird {
int colour;
float probability;
};
struct path {
int from;
int to;
};
vector <vector <bird>> birdChanges;
vector <int> layer;
vector <double> savedAnswers;
stack <path> nextBirds;
int fromBird;
//Self loop
void selfLoop(){
float totalOut = 0;
for (int i = 0; i < birdChanges.size(); i++) {
for (int j = 0; j < birdChanges[i].size(); j++) {
totalOut += birdChanges[i][j].probability;
}
if (totalOut < 1) {
bird a;
a.colour = i;
a.probability = 1 - totalOut;
birdChanges[i].push_back(a);
}
totalOut = 0;
}
}
double fillingUp(double momentarilyProbability, long long int numberOfBerries){
int layernumber=0;
while (layer[numberOfBerries - (1+layernumber)] == 0) {
layernumber++;
if (numberOfBerries == layernumber) {
break;
}
}
layernumber = layer.size() - layernumber;
path direction;
int b;
if (layernumber != 0) {
b= birdChanges[nextBirds.top().from][nextBirds.top().to].colour;//Usikker
}
else {
b = fromBird;
}
while (layer[numberOfBerries - 1] == 0) {
//int a = birdChanges[nextBirds.top().from][nextBirds.top().to].colour;
if (layernumber != 0) {
momentarilyProbability *= birdChanges[nextBirds.top().from][nextBirds.top().to].probability;
//b = birdChanges[nextBirds.top().from][nextBirds.top().to].colour;
}
for (int i = 0; i < birdChanges[b].size(); i++) {
direction.from = b;
direction.to = i;
//cout << endl << "Stacking " << b << " " << birdChanges[b][i].colour;
nextBirds.push(direction);
layer[layernumber]++;
}
layernumber++;
b = birdChanges[nextBirds.top().from][nextBirds.top().to].colour;
}
//cout << "Returning" << endl;
return momentarilyProbability *= birdChanges[nextBirds.top().from][nextBirds.top().to].probability;;
}
//DFS
void depthFirstSearch(int fromBird, long long int numberOfBerries) {
//Stack for next birds (stack)
path a;
double momentarilyProbability = 1;//Momentarily probability (float)
momentarilyProbability=fillingUp(1, numberOfBerries);
//cout << "Back " << momentarilyProbability << endl;
//Previous probabilities (stack)
while (layer[0] != 0) {
//cout << "Entering" << endl;
while (layer[numberOfBerries - 1] != 0) {
savedAnswers[birdChanges[nextBirds.top().from][nextBirds.top().to].colour] += momentarilyProbability;
//cout << "Probability for " << birdChanges[nextBirds.top().from][nextBirds.top().to].colour << " is " << momentarilyProbability << endl;
momentarilyProbability = momentarilyProbability / birdChanges[nextBirds.top().from][nextBirds.top().to].probability;
nextBirds.pop();
layer[numberOfBerries - 1]--;
if (layer[numberOfBerries - 1] != 0) {
momentarilyProbability *= birdChanges[nextBirds.top().from][nextBirds.top().to].probability;
}
}
if (layer[0] != 0) {
int k = 1;
while (layer[layer.size() - k]==0&&k+1<=layer.size()) {
//cout << "start" << endl;
momentarilyProbability = momentarilyProbability / birdChanges[nextBirds.top().from][nextBirds.top().to].probability;
//cout << "Popping " << nextBirds.top().from << birdChanges[nextBirds.top().from][nextBirds.top().to].colour << endl;
nextBirds.pop();
//cout << "k " << k << endl;
layer[numberOfBerries - 1 - k]--;
k++;
//cout << "end" << endl;
}
}
if (layer[0] != 0) {
//cout << 1 << endl;
//cout << "Filling up from " << nextBirds.top().from << birdChanges[nextBirds.top().from][nextBirds.top().to].colour << endl;
momentarilyProbability = fillingUp(momentarilyProbability, numberOfBerries);
}
}
//Printing out
for (int i = 1; i < savedAnswers.size(); i++) {
cout << savedAnswers[i] << " ";
}
cout << endl;
}
int main() {
int numberOfColours;
int possibleColourchanges;
cin >> numberOfColours >> possibleColourchanges;
birdChanges.resize(numberOfColours+1);
int from, to;
float probability;
for (int i = 0; i < possibleColourchanges; i++) {
cin >> from >> to >> probability;
bird a;
a.colour = to;
a.probability = probability;
birdChanges[from].push_back(a);
}
selfLoop();
int numberOfQuestions;
cin >> numberOfQuestions;
long long int numberOfBerries;
for (int i = 0; i < numberOfQuestions; i++) {
cin >> fromBird >> numberOfBerries;
savedAnswers.assign(numberOfColours + 1, 0);
layer.resize(numberOfBerries, 0);
//DFS
depthFirstSearch(fromBird, numberOfBerries);
}
system("pause");
}
Fast explanation of how to do this with the concept of a Markov Chain:
Basic algorithm:
Input: starting configuration vector b of probabilities of
being in a vertex after 0 steps,
Matrix A that stores the probability weights,
in the scheme of an adjacency matrix
precision threshold epsilon
Output:
an ending configuration b_inf of probabilities after infinite steps
Pseudocode:
b_old = b
b_new = A*b
while(difference(b_old, b_new) > epsilon){
b_old = b_new
b_new = A*b_old
}
return b_new
In this algorithm, we essentially compute potencies of the probability matrix and look for when those become stable.
b are the probabilities to be at a vertex after no steps where taken (so, in your case, every entry being zero except for the start vertex, which is one)
A*b are those after one step was taken
A^2 * b are those after two steps were taken, A^n * b after n steps.
When A^n * b is nearly the same as A^n-1 * b, we assume that nothing big will happen to it any more, that it is basically the same as A^infinity * b
One can mock this algorithm with some examples, like an edge that leads in a subgraph with a very small probability that will result one being in the subgraph with probability 1 after infinite steps, but for example from reality, it will work.
For the difference, the euclidean distance should work well, but essentially any norm does, you could also go with maximum or manhattan.
Note that I present a pragmatic point of view, a mathematician would go far more into detail about under which properties of A it will converge how fast for which values of epsilon.
You might want to use a good library for matrices for that, like Eigen.
EDIT:
Reading the comment of Jarod42, I realize that your amount of steps are given. In that case, simply go with A^steps * b for the exact solution. Use a good library for a fast computation of the potency.