graphtreepartition.cc.svn-base

Probabilistic graphical models in matlab.

	// *********** CODE actually begins here ************ //
	
	//cout << "******************************************************" << endl;
	//cout << "We are starting the growing of a new tree." << endl;
	
	tie(v, v_end) = vertices(g);
	
	if (v == v_end)
	{
		return false;
	}
	
	//++teton;
	/*
	if (teton ==6)
	{
		exit (0);	
	}
	*/
	
	// Initialize ALL vertices to white, and they should be their OWN parents
	
	for (tie(v, v_end) = vertices(g); v != v_end; ++v)
	{
		put(color,*v, Color::white());
		put(parent,*v, *v);
	}
	
	/*
	 
	 I am not sure if we should do any sort of pruning for Factor Graphs
	 
	 prune_out_degree_1(g, color, parent);
	 
	 //cout << "end of pruning_1" << endl;
	 
	 prune_out_degree_2(g, color, parent);
	 
	 cout << "End of pruning_2" << endl;
	 
	 */
	
	for (tie(v, v_end) = vertices(g); v != v_end; ++v)
	{
		if (get(color,*v) == Color::green())
		{
			continue;
		}
		
		vertex_with_lowest_degree = *v;
		b = out_degree(vertex_with_lowest_degree, g);
		break;
	}
	
	for (tie(v, v_end) = vertices(g); v != v_end; ++v)
	{
		free_edges.push_back(out_degree(*v, g));
		
		if (get(color,*v) == Color::green())
		{
			continue;
		}
		
		a = out_degree(*v, g);
		
		if ( (a < b) && (!get_phantom_potential(g[*v])) )
		{
			b = a;
			vertex_with_lowest_degree = *v;
		}
		
	}
	
	// Obtain vertex with lowest degree
	
	gray_vertices.push(make_tuple (vertex_with_lowest_degree, compute_possible_choices(g, color, vertex_with_lowest_degree), out_degree(vertex_with_lowest_degree, g), 0, 0));
	
	put(color, vertex_with_lowest_degree, Color::gray());
	
	unsigned step = 0;
	
	VertexDescriptor candidate;
	
	// Follow subroutine for a given subgraph
	
	//cout << "Entering main loop" << endl;
	
	while(!gray_vertices.empty())
	{
		++step;
		
		bool non_trapped_vertex(false);
		
		// It is *only* a candidate because of the backtrack mechanism
		
		candidate = gray_vertices.top().get<0>();
		
		gray_vertices.pop();
		
		// Here we test if it is still a candidate, e.g. if it hasn't been pruned out already
		
		if (get(color, candidate) != Color::gray())
		{
			//cerr << "Candidate " << candidate+1 << " was refused " << endl;
			
			continue;
		}
		
		if (get_phantom_potential(g[candidate]))
		{
			//cerr << "Candidate was refused as a phantom potential." << endl;
			//exit(0);
			
			continue;
		}
		
		VariableMessageNode < DiscreteRandomVariable, DiscretePotential, VertexDescriptor> * ptr = dynamic_cast < VariableMessageNode < DiscreteRandomVariable, DiscretePotential, VertexDescriptor> * > (g [candidate]);
		
		if (ptr)
		{
			//cout << "We explore variable " << ptr->get_random_variable()->get_index() << endl;
			non_trapped_vertex = check_rv(g, candidate, backtrack_list_gray, color, parent);
		}
		
		else
		{
			//PotentialMessageNode < DiscreteRandomVariable, DiscretePotential, VertexDescriptor> * ptr = dynamic_cast < PotentialMessageNode < DiscreteRandomVariable, DiscretePotential, VertexDescriptor> * > (g [candidate]);
			
			//cout << "We explore potential " << g [candidate]->index << endl;
			non_trapped_vertex = check_potential (g, candidate, backtrack_list_gray, backtrack_list_black, color, parent, free_edges);
		}
		
		
		if (!non_trapped_vertex) // nothing happens... we don't choose this vertex
		{
			//cout << "Backtracking." << endl;
			backtrack_list_black.clear();
			backtrack_list_gray.clear();
			//put(color, candidate, Color::black());
		}
		
		else // incorporate backtrack stuff into main	
		{
			if (ptr)
			{
				//cout << "We added vertex " << ptr->get_random_variable()->get_index() << endl;
			}
			
			else
				
			{
				//cout << "We added potential " << g [candidate]->index << endl;
			}
			
			put(color, candidate, Color::red());
			
			for ( list<VertexDescriptor>::iterator it = backtrack_list_black.begin() ; it != backtrack_list_black.end(); ++it)
			{
				put(color, *it, Color::black());
				free_edges [*it] = free_edges [*it] -1;
			}
			
			for ( list<VertexDescriptor>::iterator it = backtrack_list_gray.begin() ; it != backtrack_list_gray.end(); ++it)
			{
				put(color, *it, Color::gray());
				
				gray_vertices.push(make_tuple (*it,  compute_possible_choices(g,  color, *it), out_degree(*it, g), step, 0));
				
				free_edges [*it] = free_edges [*it] -1;
			}
			
			backtrack_list_black.clear();
			backtrack_list_gray.clear();
		}
		
    }
	
	return true;
}

template <>
inline unsigned int check_phantom_potential_status(boost::graph_traits<SubDiscretePairwiseGraph>::vertex_descriptor , DiscretePairwiseGraphColorMap &, SubDiscretePairwiseGraph &)
{
	return 0;
}

template <>
unsigned int check_phantom_potential_status(boost::graph_traits<SubDiscreteFactorGraph>::vertex_descriptor u, DiscreteFactorGraphColorMap & color, SubDiscreteFactorGraph & g)
{
	typedef graph_traits <DiscreteFactorGraph>::vertex_descriptor VertexDescriptor;

	PotentialMessageNode < DiscreteRandomVariable, DiscretePotential, VertexDescriptor> * ptr = dynamic_cast < PotentialMessageNode < DiscreteRandomVariable, DiscretePotential, VertexDescriptor> * > (g [u]);
	
	if (!ptr)
	{
		return 0;
	}
	
	graph_traits <DiscreteFactorGraph>::adjacency_iterator v, v_end;
	typedef property_traits<DiscreteFactorGraphColorMap>::value_type ColorValue;
	typedef color_traits<ColorValue> Color;
	
	if ( (get(color,u) != Color::red() ) &&  (!ptr->phantom_potential ))
	{
		return 0;
	} 
	
	bool initial_status = ptr->phantom_potential;
	
	ptr->phantom_potential = false;
	
	unsigned int black_neighbors (0);
	
	//cout << "Examining potential " << ptr->index << " for phantom. ";
	
	for (tie (v, v_end) = adjacent_vertices(u, g); v != v_end; ++v)
	{
		//cout << "Neighbour " << g[*v]->get_random_variable()->get_index();
		
		if ( get(color,*v) == Color::red() )
		{
			//cout << ". Color is red." << endl;
			continue;
		} 
		

		if ( get(color,*v) == Color::black() )
		{
			//cout << ". Color is black." << endl;
		} 
		
		if ( get(color,*v) == Color::gray() )
		{
			//cout << ". Color is gray." << endl;
		} 
		
		if ( get(color,*v) == Color::white() )
		{
			//cout << ". Color is white." << endl;
		} 
		
		//cout << "Error: we encountered a strange color." << endl;
		
		++black_neighbors;
	}
	
	if ( black_neighbors > 1)
	{
		ptr->phantom_potential = true;
	}
	
	if (ptr->phantom_potential )
	{
		//cout << "Phantom." << endl;		
	}
	
	else
	{
		//cout << "Refused." << endl;	
	}
	
	if ((!ptr->phantom_potential) && (initial_status))
	{
		return 1;
	}
	
	if ((!ptr->phantom_potential) && (!initial_status))
	{
		cout << "Strange " << endl;
		return 0;
	}	
	
	if ((ptr->phantom_potential) && (!initial_status))
	{
		assert (get(color,u) == Color::red());
		return 2;
	}
		
	if ((ptr->phantom_potential) && (initial_status))
	{
		return 3;
	}
	
	cout << "FATAL ERROR" << endl;
		
	return 4;
}

template <class SubGraph, class ColorMap, class ParentMap>
void partition_graph_into_trees_implementation (SubGraph & g, ColorMap color, ParentMap parent)
{	
	typedef typename SubGraph::graph_type Graph;
	
	typedef typename property_traits<ColorMap>::value_type ColorValue;
	typedef color_traits<ColorValue> Color;
	
	typename graph_traits <Graph>::vertex_iterator u, u_end;
	typename graph_traits <Graph>::adjacency_iterator v, v_end;
	
	
	cout << "Starting the partition... " << flush;
	
	SubGraph * root_subgraph = new SubGraph();
	copy_graph(g,*root_subgraph);
	
	Graph * working_graph = new Graph();
	copy_graph(*root_subgraph,*working_graph);
	
	SubGraph & working_subgraph = root_subgraph->create_subgraph();
	
	working_subgraph = *root_subgraph; // assignment
	
	// The working graph is a GRAPH, not a subgraph. We make all the color operations on this graph.
	// We create it at each iteration.
	
	while ( get_tree(*working_graph, color, parent) )
	{
		
		// New tree will only be used to contain the correct stuff. It is directly linked to the parent graph as one of the subgraph.
		// We need not to bother with it.
		
		SubGraph & new_tree = g.create_subgraph();
		
		// The root subgraph will contain 
		
		SubGraph & new_working_subgraph = root_subgraph->create_subgraph();
		
		for (tie(u, u_end) = vertices(working_subgraph); u != u_end; ++u)
		{
			//cout << "nago" << endl;
			
			switch (check_phantom_potential_status(*u, color, working_subgraph))
			{
				case 0: break;
					
				case 1: // the potential was a phantom potential, but has just become 'normal'
				{
					//cout << "ARA" << endl;
					continue;
				}
					break;
					
				case 2: // the potential has just become a phantom potential
				{
					//cout << "ARA" << endl;
					
					add_vertex( working_subgraph.local_to_global(*u), new_tree);
					add_vertex( working_subgraph.local_to_global(*u), new_working_subgraph);
					continue;
				}
					break;
					
				case 3: // the potential was a phantom potential and continues to be a phantom potential
				{
					add_vertex( working_subgraph.local_to_global(*u), new_working_subgraph);
					continue;
				}	
					break;
					
				default: break;
			}
					
			//cout << "bobo" << endl;
			
			if (get (color,*u) == Color::red() )
			{
				//cout << "We encountered red vertex " << *u+1 << " ( " << working_subgraph.local_to_global(*u) + 1<< " )" << endl;
				add_vertex( working_subgraph.local_to_global(*u), new_tree);
				continue;
			}
			
			if (get (color,*u) == Color::green() )
			{
				//cout << "We encoutered green vertex " << *u+1 << " ( " << working_subgraph.local_to_global(*u) + 1<< " )" << endl;
				
				// Here we should test if the parent is Red Or Not
				
				if (get( color, get(parent, *u)) == Color::red() )
				{
					//cout << "Adding it because parent is Red" << endl;
					add_vertex( working_subgraph.local_to_global(*u), new_tree);
					continue;
				}
				
			}
			
			// Else...
			
			add_vertex( working_subgraph.local_to_global(*u), new_working_subgraph);
		}
		
		// Clearing the Graph
		
		delete working_graph;
		working_graph = new Graph();
		
		// Here the operations are very important, we create a new working subgraph, and make it identical to 
		
		copy_graph(new_working_subgraph, *working_graph);
		
		working_subgraph = new_working_subgraph; // trying assignment operator
	}
	
	delete working_graph;
	delete root_subgraph;
	
	//display_subgraph(g);
	
	typename subgraph<Graph>::children_iterator s, s_end;	
	
	unsigned int total_vertices = 0;
	
	for ( tie(s, s_end) = g.children(); s != s_end; ++s)
	{
		total_vertices += num_vertices(*s);
	}
	
	assert (total_vertices == num_vertices(g));
	
	cout << "We are finished with the partition." << endl;
}

template <class Graph>
void partition_graph_into_trees (Graph & g)
{
	vector <typename graph_traits <Graph>::vertex_descriptor> v_1 (num_vertices(g));
	vector <default_color_type> v_2 (num_vertices(g));
	
	partition_graph_into_trees_implementation(g, make_iterator_property_map(v_2.begin(), get(vertex_index,g)), make_iterator_property_map( v_1.begin(), get(vertex_index,g)));
}