Next Previous Table of Contents
Imagine you want to find the maximum flow from one corner to the oposite corner in a simple cubic lattice. You would first create a Graph. The constructor takes a Geometry and a size argument. To create a simple cubic lattice with L=25 use:
Graph G = Graph(Graph::SIMPLE_CUBIC, 25)
cout << G.maxFlow(G.getFirstVertex(), G.getLastVertex()) << endl;
An Ising model can be visualized as a set of little magnets that can only point north or south, in other words they have two distinct states. Often these magnets are referred to as spins and instead of north or south up and down are used to indicate the state of the spins.
Two magnets that are close usually interact with each other. In this model this interaction is represented by a bond. The strength of the bond indicates how strongly two spins interact. In the random bond model this strength can vary but is always positive. Positive means here that the spins want to be aligned.
If we put such a spin on each site in a simple cubic lattice and don't do anything else they would all align and form a ferromagnet. To make it more interesting we can fix the boundary conditions. All spins on one face of the cube are be fixed in an up position, all spins on the opposite face are fixed in a down position. We now want to calculate the interface between the up and the down domain.
Each spin is mapped to a vertex and each bond to an edge where the strength of the bond becomes the capacity of the edge. In addition two so called ghost sites s and t are used. These impose the boundary conditions on the system.
Each spin on the face where all the spins face upward is connected to s with a large, ideally infinite, capacity edge. All spins on the opposite face are connected to t again with a large capacity edge.
Now by finding the maximum flow from s to t we can find the interface between the two domains.
First we need to create a Graph object. SmaGL provides an inbuilt routine to create a graph representing a simple cubic lattice.
Graph G = Graph(Graph::SIMPLE_CUBIC, 25);
Vertex& s = G.addVertex(); Vertex& t = G.addVertex();
Next s and t need to be connected with the proper faces of the cube. The simple cube is created by adding edges row by row. Therefore the top face has the vertices 1 to L^2. We can get an iterator for the vector in which the vertices are stored. This allows us to just go through the first L^2 vertices and connect them to s.
/* Connect the top plane to s */
vector<Vertex*>::iterator p1 = G.vertices.begin();
for (int i = 0; i < (length * length); i++){
G.addEdge(s, *(*p1),500,1);
G.addEdge(*(*p1), s, 0, 1); // Adding the reverse edge is very important!
p1++; // Go to the next Vertex in the vector
}
The same can be done with the bottom face. Here we have to be careful! We added two more vertices (s and t) after creating the simple cube.
p1 = G.vertices.end(); // Points to the position after the last Vertex
p1 = p1 - 3; // p1 now points to last vertex before s
for (int i = 0; i < (length * length); i++){
G.addEdge(*(*p1), t, 500, 1);
G.addEdge(t, *(*p1), 0, 1);
p1--; // Going backwards through the vector of vertices.
}
Creating the right graph is usually the hardest part. So far all edges within the cube have capacity or strength one. We do want a random network though. We therefore need to randomize the capacity of the edges within the simple cube without changing the ones connecting to s and t. The following code does this:
/* Randomize the bonds */
int seed = time(0);
srand(seed); // Initialize random number generator
p1 = G.vertices.begin(); // Get a pointer to the first Vertex in the Graph.
/* All vertices but the last two belong to the cube.*/
while( p1 < G.vertices.end() - 2 ){
typedef list<Edge*>::const_iterator LI;
/* Each Vertex has the outgoing edges stored in a list. Here we go through
* those edges
*/
for(LI p2 = (*p1)->outEdges.begin(); p2 != (*p1)->outEdges.end(); ++p2){
double cap = (double(rand()) / RAND_MAX) * 50; // random double [0,49]
/* if the target of the Edge is either t or s don't do anything.
if ( (*p2)->getTarget() == t || (*p2)->getTarget() == s) continue;
/* Set new capcity. Reverse Edge should have same capacity. */
(*p2)->setCapacity((int) cap);
G.reverseEdge((*p2))->setCapacity((int) cap);
}
p1++; // Go to next Vertex
}
cout << G.maxFlow(s , t) << endl;
Next Previous Table of Contents